From d97cf22b50fa14f3a2d916c8a4531b77b3e46551 Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Tue, 21 May 2019 00:28:27 -0700 Subject: [PATCH 01/13] Add a new way to calculate the solarposition --- pvlib/solarposition.py | 136 ++++++++++++++- pvlib/test/test_solarposition.py | 290 ++++++++++++++++++++++--------- 2 files changed, 344 insertions(+), 82 deletions(-) diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index e2dc536936..f342abc70f 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -7,6 +7,7 @@ # Will Holmgren (@wholmgren), University of Arizona, 2014 # Tony Lorenzo (@alorenzo175), University of Arizona, 2015 # Cliff hansen (@cwhanse), Sandia National Laboratories, 2018 +# Daniel Lassahn (@meteoDaniel), meteocontrol Energy and Weather Services, 2019 from __future__ import division import os @@ -18,7 +19,6 @@ from imp import reload except ImportError: pass - import numpy as np import pandas as pd import warnings @@ -27,8 +27,12 @@ from pvlib.tools import datetime_to_djd, djd_to_datetime from pvlib._deprecation import deprecated - NS_PER_HR = 1.e9 * 3600. # nanoseconds per hour +JULIAN_2000 = 2451544.5 +DT_2000 = dt.datetime(2000, 1, 1) +TIMESTAMP_2000 = DT_2000.timestamp() +DAY_SECONDS = 60 * 60 * 24 +JULIAN_YEARS = 365.2425 def get_solarposition(time, latitude, longitude, @@ -67,6 +71,9 @@ def get_solarposition(time, latitude, longitude, 'nrel_c' uses the NREL SPA C code [3]: :py:func:`spa_c` + 'spencer_mc' uses the Spencer formula [4] :py:func:`spencer_mc` + + temperature : float, default 12 Degrees C. @@ -81,6 +88,9 @@ def get_solarposition(time, latitude, longitude, solar radiation applications. Solar Energy, vol. 81, no. 6, p. 838, 2007. [3] NREL SPA code: http://rredc.nrel.gov/solar/codesandalgorithms/spa/ + + [4] Spencer (1972) can be found in "Solar energy fundamentals and + modeling techniques" from Zekai Sen """ if altitude is None and pressure is None: @@ -114,6 +124,9 @@ def get_solarposition(time, latitude, longitude, elif method == 'ephemeris': ephem_df = ephemeris(time, latitude, longitude, pressure, temperature, **kwargs) + elif method == 'spencer_mc': + ephem_df = spencer_mc(time, latitude, longitude) + else: raise ValueError('Invalid solar position method') @@ -1440,3 +1453,122 @@ def sun_rise_set_transit_geometric(times, latitude, longitude, declination, sunset = _local_times_from_hours_since_midnight(times, sunset_hour) transit = _local_times_from_hours_since_midnight(times, transit_hour) return sunrise, sunset, transit + + +def spencer_mc(times, latitude, longitude): + """ + Calculate the solar position using a python implementation of the + Spencer (1972) formulation + + Parameters + ---------- + times : pandas.DatetimeIndex + Corresponding timestamps, must be localized to the timezone for the + ``latitude`` and ``longitude``. + latitude : float + Latitude in degrees, positive north of equator, negative to south + longitude : float + Longitude in degrees, positive east of prime meridian, negative to west + + Returns + ------- + DataFrame + The DataFrame will have the following columns: + zenith (degrees), + elevation (degrees), + azimuth (degrees), + equation_of_time (seconds), + eccentricity, + declination (degrees). + + References + ---------- + [4] Spencer (1972) can be found in + Sen, Zekai. Solar energy fundamentals and modeling techniques: + atmosphere, environment, climate change and renewable energy. + Springer Science & Business Media, 2008. + """ + + julians = datetime2julian(times) + julians_2000 = np.asarray(julians, dtype=np.float) - JULIAN_2000 + + lat, lat_deg = np.radians(latitude), latitude + lon, lon_deg = np.radians(longitude), longitude + + # Compute fractional year (gamma) in radians + gamma = 2 * np.pi * (julians_2000 % JULIAN_YEARS) / JULIAN_YEARS + cos_gamma = np.cos(gamma), np.cos(gamma * 2), np.cos(gamma * 3) + sin_gamma = np.sin(gamma), np.sin(gamma * 2), np.sin(gamma * 3) + day_time = (julians_2000 % 1) * 24 + + # Eccentricity: correction factor of the earth's orbit. + eccentricity = (1.00011 + 0.034221 * cos_gamma[0] + 0.001280 * + sin_gamma[0] + 0.000719 * cos_gamma[1] + 0.000077 + * sin_gamma[1]) + + # declination. + declination = (0.006918 - 0.399912 * cos_gamma[0] + + 0.070257 * sin_gamma[0] - + 0.006758 * cos_gamma[1] + 0.000907 * sin_gamma[1] - + 0.002697 * cos_gamma[2] + 0.001480 * sin_gamma[2]) + + # Equation of time (difference between standard time and solar time). + eot = (0.000075 + 0.001868 * cos_gamma[0] - 0.032077 * sin_gamma[0] - + 0.014615 * cos_gamma[1] - 0.040849 * sin_gamma[1]) * 229.18 + + # True local time + tlt = (day_time + lon_deg / 15 + eot / 60) % 24 - 12 + + # Solar hour angle + ha = np.radians(tlt * 15) + + # Calculate sun elevation. + sin_sun_elevation = ( + np.sin(declination) * np.sin(lat) + np.cos(declination) * + np.cos(lat) * np.cos(ha) + ) + + # Compute the sun's elevation and zenith angle. + elevation = np.arcsin(sin_sun_elevation) + zenith = np.pi / 2 - elevation + + # Compute the sun's azimuth angle. + y = -(np.sin(lat) * np.sin(elevation) - np.sin(declination)) / \ + (np.cos(lat) * np.cos(elevation)) + azimuth = np.arccos(y) + + # Convert azimuth angle from 0-pi to 0-2pi. + tlt_filter = 0 <= tlt + azimuth[tlt_filter] = 2 * np.pi - azimuth[tlt_filter] + + result = pd.DataFrame({'zenith': np.degrees(zenith), + 'elevation': np.degrees(elevation), + 'azimuth': np.degrees(azimuth), + 'eccentricity': eccentricity, + 'declination': declination, + 'equation_of_time': eot}, + index=times) + return result + + +def datetime2julian(times): + """ + Transforms pd.DateTimeIndex to Julian days + + Parameters + ---------- + times : pandas.DatetimeIndex + Corresponding timestamps, must be localized to the timezone for the + ``latitude`` and ``longitude``. + + Returns + ------- + Float64Index + The float index contains julian dates + """ + + delta = times - DT_2000 + return ( + JULIAN_2000 + + delta.days + (delta.seconds + delta.microseconds / 1e6) / DAY_SECONDS + ) diff --git a/pvlib/test/test_solarposition.py b/pvlib/test/test_solarposition.py index 2bfadb2c9d..f33273dfe0 100644 --- a/pvlib/test/test_solarposition.py +++ b/pvlib/test/test_solarposition.py @@ -16,12 +16,11 @@ from conftest import (fail_on_pvlib_version, requires_ephem, needs_pandas_0_17, requires_spa_c, requires_numba) - # setup times and locations to be tested. -times = pd.date_range(start=datetime.datetime(2014,6,24), - end=datetime.datetime(2014,6,26), freq='15Min') +times = pd.date_range(start=datetime.datetime(2014, 6, 24), + end=datetime.datetime(2014, 6, 26), freq='15Min') -tus = Location(32.2, -111, 'US/Arizona', 700) # no DST issues possible +tus = Location(32.2, -111, 'US/Arizona', 700) # no DST issues possible times_localized = times.tz_localize(tus.tz) tol = 5 @@ -58,7 +57,8 @@ def expected_solpos_multi(): 'apparent_zenith': [50.111622, 50.478260], 'azimuth': [194.340241, 194.311132], 'apparent_elevation': [39.888378, 39.521740]}, - index=[['2003-10-17T12:30:30Z', '2003-10-18T12:30:30Z']]) + index=[ + ['2003-10-17T12:30:30Z', '2003-10-18T12:30:30Z']]) @pytest.fixture() @@ -126,7 +126,7 @@ def test_deprecated_07(): @requires_spa_c def test_spa_c_physical(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003,10,17,12,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), periods=1, freq='D', tz=golden_mst.tz) ephem_data = solarposition.spa_c(times, golden_mst.latitude, golden_mst.longitude, @@ -138,7 +138,7 @@ def test_spa_c_physical(expected_solpos, golden_mst): @requires_spa_c def test_spa_c_physical_dst(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.spa_c(times, golden.latitude, golden.longitude, @@ -149,7 +149,7 @@ def test_spa_c_physical_dst(expected_solpos, golden): def test_spa_python_numpy_physical(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003,10,17,12,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), periods=1, freq='D', tz=golden_mst.tz) ephem_data = solarposition.spa_python(times, golden_mst.latitude, golden_mst.longitude, @@ -162,7 +162,7 @@ def test_spa_python_numpy_physical(expected_solpos, golden_mst): def test_spa_python_numpy_physical_dst(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.spa_python(times, golden.latitude, golden.longitude, @@ -176,7 +176,7 @@ def test_spa_python_numpy_physical_dst(expected_solpos, golden): @requires_numba def test_spa_python_numba_physical(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003,10,17,12,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), periods=1, freq='D', tz=golden_mst.tz) ephem_data = solarposition.spa_python(times, golden_mst.latitude, golden_mst.longitude, @@ -190,7 +190,7 @@ def test_spa_python_numba_physical(expected_solpos, golden_mst): @requires_numba def test_spa_python_numba_physical_dst(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.spa_python(times, golden.latitude, golden.longitude, pressure=82000, @@ -269,20 +269,20 @@ def test_sun_rise_set_transit_ephem(expected_rise_set_ephem, golden): columns=['sunrise', 'sunset'], dtype='datetime64[ns]') expected['sunrise'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise']]) expected['sunset'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunset']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunset']]) expected['transit'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit']]) result = solarposition.sun_rise_set_transit_ephem(times, golden.latitude, @@ -308,24 +308,29 @@ def test_sun_rise_set_transit_ephem(expected_rise_set_ephem, golden): columns=['sunrise', 'sunset'], dtype='datetime64[ns]') expected['sunrise'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise']]) expected['sunset'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset']]) expected['transit'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit']]) result = solarposition.sun_rise_set_transit_ephem(times, - golden.latitude, golden.longitude, next_or_previous='previous', - altitude=golden.altitude, pressure=0, temperature=11, horizon='-0:34') + golden.latitude, + golden.longitude, + next_or_previous='previous', + altitude=golden.altitude, + pressure=0, + temperature=11, + horizon='-0:34') # round to nearest minute result_rounded = pd.DataFrame(index=result.index) for col, data in result.iteritems(): @@ -338,8 +343,13 @@ def test_sun_rise_set_transit_ephem(expected_rise_set_ephem, golden): for col, data in expected.iteritems(): expected[col] = data.dt.tz_convert('UTC') result = solarposition.sun_rise_set_transit_ephem(times, - golden.latitude, golden.longitude, next_or_previous='previous', - altitude=golden.altitude, pressure=0, temperature=11, horizon='-0:34') + golden.latitude, + golden.longitude, + next_or_previous='previous', + altitude=golden.altitude, + pressure=0, + temperature=11, + horizon='-0:34') # round to nearest minute result_rounded = pd.DataFrame(index=result.index) for col, data in result.iteritems(): @@ -368,13 +378,16 @@ def test_sun_rise_set_transit_ephem_horizon(golden): ]).tz_localize('MST') # center of sun disk center = solarposition.sun_rise_set_transit_ephem(times, - latitude=golden.latitude, longitude=golden.longitude) + latitude=golden.latitude, + longitude=golden.longitude) edge = solarposition.sun_rise_set_transit_ephem(times, - latitude=golden.latitude, longitude=golden.longitude, horizon='-0:34') + latitude=golden.latitude, + longitude=golden.longitude, + horizon='-0:34') result_rounded = (edge['sunrise'] - center['sunrise']).dt.round('min') sunrise_delta = datetime.datetime(2016, 1, 3, 7, 17, 11) - \ - datetime.datetime(2016, 1, 3, 7, 21, 33) + datetime.datetime(2016, 1, 3, 7, 21, 33) expected = pd.Series(index=times, data=sunrise_delta, name='sunrise').dt.round('min') @@ -383,7 +396,7 @@ def test_sun_rise_set_transit_ephem_horizon(golden): @requires_ephem def test_pyephem_physical(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003,10,17,12,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), periods=1, freq='D', tz=golden_mst.tz) ephem_data = solarposition.pyephem(times, golden_mst.latitude, golden_mst.longitude, pressure=82000, @@ -395,7 +408,8 @@ def test_pyephem_physical(expected_solpos, golden_mst): @requires_ephem def test_pyephem_physical_dst(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.pyephem(times, golden.latitude, golden.longitude, pressure=82000, @@ -411,7 +425,7 @@ def test_calc_time(): import math # validation from USNO solar position calculator online - epoch = datetime.datetime(1970,1,1) + epoch = datetime.datetime(1970, 1, 1) epoch_dt = pytz.utc.localize(epoch) loc = tus @@ -427,21 +441,21 @@ def test_calc_time(): actual_timestamp = (actual_time - epoch_dt).total_seconds() assert_allclose((alt.replace(second=0, microsecond=0) - - epoch_dt).total_seconds(), actual_timestamp) + epoch_dt).total_seconds(), actual_timestamp) assert_allclose((az.replace(second=0, microsecond=0) - - epoch_dt).total_seconds(), actual_timestamp) + epoch_dt).total_seconds(), actual_timestamp) @requires_ephem def test_earthsun_distance(): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), periods=1, freq='D') distance = solarposition.pyephem_earthsun_distance(times).values[0] assert_allclose(1, distance, atol=0.1) def test_ephemeris_physical(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003,10,17,12,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), periods=1, freq='D', tz=golden_mst.tz) ephem_data = solarposition.ephemeris(times, golden_mst.latitude, golden_mst.longitude, @@ -454,7 +468,7 @@ def test_ephemeris_physical(expected_solpos, golden_mst): def test_ephemeris_physical_dst(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.ephemeris(times, golden.latitude, golden.longitude, pressure=82000, @@ -466,7 +480,7 @@ def test_ephemeris_physical_dst(expected_solpos, golden): def test_ephemeris_physical_no_tz(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003,10,17,19,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 19, 30, 30), periods=1, freq='D') ephem_data = solarposition.ephemeris(times, golden_mst.latitude, golden_mst.longitude, @@ -479,7 +493,7 @@ def test_ephemeris_physical_no_tz(expected_solpos, golden_mst): def test_get_solarposition_error(golden): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), periods=1, freq='D', tz=golden.tz) with pytest.raises(ValueError): ephem_data = solarposition.get_solarposition(times, golden.latitude, @@ -492,14 +506,14 @@ def test_get_solarposition_error(golden): @pytest.mark.parametrize("pressure, expected", [ (82000, _expected_solpos_df()), (90000, pd.DataFrame( - np.array([[ 39.88997, 50.11003, 194.34024, 39.87205, 14.64151, - 50.12795]]), + np.array([[39.88997, 50.11003, 194.34024, 39.87205, 14.64151, + 50.12795]]), columns=['apparent_elevation', 'apparent_zenith', 'azimuth', 'elevation', 'equation_of_time', 'zenith'], index=['2003-10-17T12:30:30Z'])) - ]) +]) def test_get_solarposition_pressure(pressure, expected, golden): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.get_solarposition(times, golden.latitude, golden.longitude, @@ -515,14 +529,14 @@ def test_get_solarposition_pressure(pressure, expected, golden): @pytest.mark.parametrize("altitude, expected", [ (1830.14, _expected_solpos_df()), (2000, pd.DataFrame( - np.array([[ 39.88788, 50.11212, 194.34024, 39.87205, 14.64151, - 50.12795]]), + np.array([[39.88788, 50.11212, 194.34024, 39.87205, 14.64151, + 50.12795]]), columns=['apparent_elevation', 'apparent_zenith', 'azimuth', 'elevation', 'equation_of_time', 'zenith'], index=['2003-10-17T12:30:30Z'])) - ]) +]) def test_get_solarposition_altitude(altitude, expected, golden): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.get_solarposition(times, golden.latitude, golden.longitude, @@ -543,10 +557,10 @@ def test_get_solarposition_altitude(altitude, expected, golden): marks=[pytest.mark.xfail( reason='spa.calculate_deltat not implemented for numba yet')]), (67.0, 'nrel_numba') - ]) +]) def test_get_solarposition_deltat(delta_t, method, expected_solpos_multi, golden): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), periods=2, freq='D', tz=golden.tz) ephem_data = solarposition.get_solarposition(times, golden.latitude, golden.longitude, @@ -562,7 +576,7 @@ def test_get_solarposition_deltat(delta_t, method, expected_solpos_multi, def test_get_solarposition_no_kwargs(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), + times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.get_solarposition(times, golden.latitude, golden.longitude) @@ -587,7 +601,7 @@ def test_get_solarposition_method_pyephem(expected_solpos, golden): def test_nrel_earthsun_distance(): times = pd.DatetimeIndex([datetime.datetime(2015, 1, 2), - datetime.datetime(2015, 8, 2),] + datetime.datetime(2015, 8, 2), ] ).tz_localize('MST') result = solarposition.nrel_earthsun_distance(times, delta_t=64.0) expected = pd.Series(np.array([0.983289204601, 1.01486146446]), @@ -676,27 +690,27 @@ def test_analytical_azimuth(): azimuth_2 = solarposition.solar_azimuth_analytical(lat_rad, hour_angle, decl, zenith) - idx = np.where(solar_zenith < np.pi/2) + idx = np.where(solar_zenith < np.pi / 2) assert np.allclose(azimuth_1[idx], solar_azimuth.values[idx], atol=0.01) assert np.allclose(azimuth_2[idx], solar_azimuth.values[idx], atol=0.017) # test for NaN values at boundary conditions (PR #431) test_angles = np.radians(np.array( - [[ 0., -180., -20.], - [ 0., 0., -5.], - [ 0., 0., 0.], - [ 0., 0., 15.], - [ 0., 180., 20.], - [ 30., 0., -20.], - [ 30., 0., -5.], - [ 30., 0., 0.], - [ 30., 180., 5.], - [ 30., 0., 10.], - [ -30., 0., -20.], - [ -30., 0., -15.], - [ -30., 0., 0.], - [ -30., -180., 5.], - [ -30., 180., 10.]])) + [[0., -180., -20.], + [0., 0., -5.], + [0., 0., 0.], + [0., 0., 15.], + [0., 180., 20.], + [30., 0., -20.], + [30., 0., -5.], + [30., 0., 0.], + [30., 180., 5.], + [30., 0., 10.], + [-30., 0., -20.], + [-30., 0., -15.], + [-30., 0., 0.], + [-30., -180., 5.], + [-30., 180., 10.]])) zeniths = solarposition.solar_zenith_analytical(*test_angles.T) azimuths = solarposition.solar_azimuth_analytical(*test_angles.T, @@ -772,3 +786,119 @@ def test_sun_rise_set_transit_geometric(expected_rise_set_spa, golden_mst): atol=np.abs(expected_sunset_error).max()) assert np.allclose(test_transit, expected_transit, atol=np.abs(expected_transit_error).max()) + + +def test_spencer_mc(): + """ test for the calculation based on spencer 1972 """ + times = pd.date_range('2018-01-01', periods=24, freq='H') + latitude = 48.367073 + longitude = 10.868378 + assert_frame_equal(solarposition.spencer_mc(times, latitude, longitude), + pd.DataFrame([[153.40940267, -63.40940267, 21.07562885, + 1.03506345, + -0.40158632, -3.18821393], + [147.78245317, -57.78245317, 47.01029634, + 1.03506418, + -0.40152783, -3.20679537], + [139.47329285, -49.47329285, 65.74488069, + 1.03506489, + -0.40146911, -3.22536965], + [129.96316706, -39.96316706, 79.89448997, + 1.03506558, + -0.40141017, -3.24393676], + [120.04134593, -30.04134593, 91.64276999, + 1.03506625, + -0.40135099, -3.26249664], + [110.16576717, -20.16576717, 102.30874903, + 1.0350669, + -0.40129158, -3.28104929], + [100.67550276, -10.67550276, 112.71473915, + 1.03506753, + -0.40123194, -3.29959466], + [91.89177037, -1.89177037, 123.43188389, + 1.03506814, + -0.40117207, -3.31813273], + [84.16861912, 5.83138088, 134.88256108, + 1.03506874, + -0.40111197, -3.33666346], + [77.91130514, 12.08869486, 147.34079438, + 1.03506931, + -0.40105163, -3.35518683], + [73.55352662, 16.44647338, 160.85682789, + 1.03506987, + -0.40099107, -3.37370281], + [71.47824875, 18.52175125, 175.16512108, + 1.03507041, + -0.40093028, -3.39221136], + [71.89802747, 18.10197253, 189.6987874, + 1.03507093, + -0.40086925, -3.41071246], + [74.76804069, 15.23195931, 203.79458854, + 1.03507142, + -0.400808, -3.42920608], + [79.80562508, 10.19437492, 216.97266043, + 1.0350719, + -0.40074651, -3.44769219], + [86.59854728, 3.40145272, 229.07943767, + 1.03507237, + -0.40068479, -3.46617075], + [94.71660269, -4.71660269, 240.25140046, + 1.03507281, + -0.40062285, -3.48464174], + [103.77123543, -13.77123543, 250.81773352, + 1.03507323, + -0.40056067, -3.50310512], + [113.42256954, -23.42256954, 261.24657665, + 1.03507363, + -0.40049826, -3.52156088], + [123.35015695, -33.35015695, 272.17648044, + 1.03507402, + -0.40043563, -3.54000897], + [133.18907815, -43.18907815, 284.56071242, + 1.03507438, + -0.40037276, -3.55844937], + [142.4005988, -52.4005988, 299.96695581, + 1.03507473, + -0.40030966, -3.57688204], + [150.00459208, -60.00459208, 320.87720374, + 1.03507506, + -0.40024633, -3.59530697], + [154.24856743, -64.24856743, 349.27789739, + 1.03507536, + -0.40018277, -3.61372411]], + index=times, + columns=['zenith', + 'elevation', + 'azimuth', + 'eccentricity', + 'declination', + 'equation_of_time']) + ) + + +def test_datetime2julians(): + """ test transformation from datetime to julians """ + times = pd.date_range('2018-01-01', periods=24, freq='H') + julians = solarposition.datetime2julian(pd.to_datetime(times)) + np.testing.assert_array_almost_equal(julians, + np.array([2458119.5, 2458119.54166667, + 2458119.58333333, + 2458119.625, + 2458119.66666667, + 2458119.70833333, + 2458119.75, 2458119.79166667, + 2458119.83333333, + 2458119.875, + 2458119.91666667, + 2458119.95833333, + 2458120., 2458120.04166667, + 2458120.08333333, + 2458120.125, + 2458120.16666667, + 2458120.20833333, + 2458120.25, 2458120.29166667, + 2458120.33333333, + 2458120.375, + 2458120.41666667, + 2458120.45833333]) + ) From 07c9ea38a00d8ba81244c8a49ab92201e777065b Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Tue, 21 May 2019 00:47:51 -0700 Subject: [PATCH 02/13] Tidy --- pvlib/solarposition.py | 31 +++++++++++++++---------------- 1 file changed, 15 insertions(+), 16 deletions(-) diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index f342abc70f..e451be59c9 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -1492,8 +1492,7 @@ def spencer_mc(times, latitude, longitude): julians = datetime2julian(times) julians_2000 = np.asarray(julians, dtype=np.float) - JULIAN_2000 - lat, lat_deg = np.radians(latitude), latitude - lon, lon_deg = np.radians(longitude), longitude + lat = np.radians(latitude) # Compute fractional year (gamma) in radians gamma = 2 * np.pi * (julians_2000 % JULIAN_YEARS) / JULIAN_YEARS @@ -1502,30 +1501,30 @@ def spencer_mc(times, latitude, longitude): day_time = (julians_2000 % 1) * 24 # Eccentricity: correction factor of the earth's orbit. - eccentricity = (1.00011 + 0.034221 * cos_gamma[0] + 0.001280 * - sin_gamma[0] + 0.000719 * cos_gamma[1] + 0.000077 + eccentricity = (1.00011 + 0.034221 * cos_gamma[0] + 0.001280 + * sin_gamma[0] + 0.000719 * cos_gamma[1] + 0.000077 * sin_gamma[1]) # declination. - declination = (0.006918 - 0.399912 * cos_gamma[0] + - 0.070257 * sin_gamma[0] - - 0.006758 * cos_gamma[1] + 0.000907 * sin_gamma[1] - - 0.002697 * cos_gamma[2] + 0.001480 * sin_gamma[2]) + declination = (0.006918 - 0.399912 * cos_gamma[0] + + 0.070257 * sin_gamma[0] - 0.006758 + * cos_gamma[1] + 0.000907 * sin_gamma[1] + - 0.002697 * cos_gamma[2] + 0.001480 * sin_gamma[2]) # Equation of time (difference between standard time and solar time). - eot = (0.000075 + 0.001868 * cos_gamma[0] - 0.032077 * sin_gamma[0] - - 0.014615 * cos_gamma[1] - 0.040849 * sin_gamma[1]) * 229.18 + eot = (0.000075 + 0.001868 * cos_gamma[0] - 0.032077 * sin_gamma[0] + - 0.014615 * cos_gamma[1] - 0.040849 * sin_gamma[1]) * 229.18 # True local time - tlt = (day_time + lon_deg / 15 + eot / 60) % 24 - 12 + tlt = (day_time + longitude / 15 + eot / 60) % 24 - 12 # Solar hour angle ha = np.radians(tlt * 15) # Calculate sun elevation. sin_sun_elevation = ( - np.sin(declination) * np.sin(lat) + np.cos(declination) * - np.cos(lat) * np.cos(ha) + np.sin(declination) * np.sin(lat) + np.cos(declination) + * np.cos(lat) * np.cos(ha) ) # Compute the sun's elevation and zenith angle. @@ -1533,8 +1532,8 @@ def spencer_mc(times, latitude, longitude): zenith = np.pi / 2 - elevation # Compute the sun's azimuth angle. - y = -(np.sin(lat) * np.sin(elevation) - np.sin(declination)) / \ - (np.cos(lat) * np.cos(elevation)) + y = -(np.sin(lat) * np.sin(elevation) - np.sin(declination)) \ + / (np.cos(lat) * np.cos(elevation)) azimuth = np.arccos(y) # Convert azimuth angle from 0-pi to 0-2pi. @@ -1571,4 +1570,4 @@ def datetime2julian(times): return ( JULIAN_2000 + delta.days + (delta.seconds + delta.microseconds / 1e6) / DAY_SECONDS - ) + ) \ No newline at end of file From 98f65393807c5dda1966d109aaa97bc327dc2c21 Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Tue, 21 May 2019 01:10:41 -0700 Subject: [PATCH 03/13] Tidy --- pvlib/solarposition.py | 2 +- pvlib/test/test_solarposition.py | 281 ++++++++++++------------------- 2 files changed, 113 insertions(+), 170 deletions(-) diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index e451be59c9..717c6c38b0 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -1570,4 +1570,4 @@ def datetime2julian(times): return ( JULIAN_2000 + delta.days + (delta.seconds + delta.microseconds / 1e6) / DAY_SECONDS - ) \ No newline at end of file + ) diff --git a/pvlib/test/test_solarposition.py b/pvlib/test/test_solarposition.py index f33273dfe0..5a0a7604ad 100644 --- a/pvlib/test/test_solarposition.py +++ b/pvlib/test/test_solarposition.py @@ -16,11 +16,12 @@ from conftest import (fail_on_pvlib_version, requires_ephem, needs_pandas_0_17, requires_spa_c, requires_numba) + # setup times and locations to be tested. -times = pd.date_range(start=datetime.datetime(2014, 6, 24), - end=datetime.datetime(2014, 6, 26), freq='15Min') +times = pd.date_range(start=datetime.datetime(2014,6,24), + end=datetime.datetime(2014,6,26), freq='15Min') -tus = Location(32.2, -111, 'US/Arizona', 700) # no DST issues possible +tus = Location(32.2, -111, 'US/Arizona', 700) # no DST issues possible times_localized = times.tz_localize(tus.tz) tol = 5 @@ -57,8 +58,7 @@ def expected_solpos_multi(): 'apparent_zenith': [50.111622, 50.478260], 'azimuth': [194.340241, 194.311132], 'apparent_elevation': [39.888378, 39.521740]}, - index=[ - ['2003-10-17T12:30:30Z', '2003-10-18T12:30:30Z']]) + index=[['2003-10-17T12:30:30Z', '2003-10-18T12:30:30Z']]) @pytest.fixture() @@ -126,7 +126,7 @@ def test_deprecated_07(): @requires_spa_c def test_spa_c_physical(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,12,30,30), periods=1, freq='D', tz=golden_mst.tz) ephem_data = solarposition.spa_c(times, golden_mst.latitude, golden_mst.longitude, @@ -138,7 +138,7 @@ def test_spa_c_physical(expected_solpos, golden_mst): @requires_spa_c def test_spa_c_physical_dst(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.spa_c(times, golden.latitude, golden.longitude, @@ -149,7 +149,7 @@ def test_spa_c_physical_dst(expected_solpos, golden): def test_spa_python_numpy_physical(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,12,30,30), periods=1, freq='D', tz=golden_mst.tz) ephem_data = solarposition.spa_python(times, golden_mst.latitude, golden_mst.longitude, @@ -162,7 +162,7 @@ def test_spa_python_numpy_physical(expected_solpos, golden_mst): def test_spa_python_numpy_physical_dst(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.spa_python(times, golden.latitude, golden.longitude, @@ -176,7 +176,7 @@ def test_spa_python_numpy_physical_dst(expected_solpos, golden): @requires_numba def test_spa_python_numba_physical(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,12,30,30), periods=1, freq='D', tz=golden_mst.tz) ephem_data = solarposition.spa_python(times, golden_mst.latitude, golden_mst.longitude, @@ -190,7 +190,7 @@ def test_spa_python_numba_physical(expected_solpos, golden_mst): @requires_numba def test_spa_python_numba_physical_dst(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.spa_python(times, golden.latitude, golden.longitude, pressure=82000, @@ -269,20 +269,20 @@ def test_sun_rise_set_transit_ephem(expected_rise_set_ephem, golden): columns=['sunrise', 'sunset'], dtype='datetime64[ns]') expected['sunrise'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise']]) expected['sunset'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunset']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunset']]) expected['transit'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit']]) result = solarposition.sun_rise_set_transit_ephem(times, golden.latitude, @@ -308,29 +308,24 @@ def test_sun_rise_set_transit_ephem(expected_rise_set_ephem, golden): columns=['sunrise', 'sunset'], dtype='datetime64[ns]') expected['sunrise'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunrise'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'sunrise']]) expected['sunset'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'sunset']]) expected['transit'] = pd.Series(index=times, data= - [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], - expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit']]) + [expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 1), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 2), 'transit'], + expected_rise_set_ephem.loc[datetime.datetime(2015, 1, 3), 'transit']]) result = solarposition.sun_rise_set_transit_ephem(times, - golden.latitude, - golden.longitude, - next_or_previous='previous', - altitude=golden.altitude, - pressure=0, - temperature=11, - horizon='-0:34') + golden.latitude, golden.longitude, next_or_previous='previous', + altitude=golden.altitude, pressure=0, temperature=11, horizon='-0:34') # round to nearest minute result_rounded = pd.DataFrame(index=result.index) for col, data in result.iteritems(): @@ -343,13 +338,8 @@ def test_sun_rise_set_transit_ephem(expected_rise_set_ephem, golden): for col, data in expected.iteritems(): expected[col] = data.dt.tz_convert('UTC') result = solarposition.sun_rise_set_transit_ephem(times, - golden.latitude, - golden.longitude, - next_or_previous='previous', - altitude=golden.altitude, - pressure=0, - temperature=11, - horizon='-0:34') + golden.latitude, golden.longitude, next_or_previous='previous', + altitude=golden.altitude, pressure=0, temperature=11, horizon='-0:34') # round to nearest minute result_rounded = pd.DataFrame(index=result.index) for col, data in result.iteritems(): @@ -378,16 +368,13 @@ def test_sun_rise_set_transit_ephem_horizon(golden): ]).tz_localize('MST') # center of sun disk center = solarposition.sun_rise_set_transit_ephem(times, - latitude=golden.latitude, - longitude=golden.longitude) + latitude=golden.latitude, longitude=golden.longitude) edge = solarposition.sun_rise_set_transit_ephem(times, - latitude=golden.latitude, - longitude=golden.longitude, - horizon='-0:34') + latitude=golden.latitude, longitude=golden.longitude, horizon='-0:34') result_rounded = (edge['sunrise'] - center['sunrise']).dt.round('min') sunrise_delta = datetime.datetime(2016, 1, 3, 7, 17, 11) - \ - datetime.datetime(2016, 1, 3, 7, 21, 33) + datetime.datetime(2016, 1, 3, 7, 21, 33) expected = pd.Series(index=times, data=sunrise_delta, name='sunrise').dt.round('min') @@ -396,7 +383,7 @@ def test_sun_rise_set_transit_ephem_horizon(golden): @requires_ephem def test_pyephem_physical(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,12,30,30), periods=1, freq='D', tz=golden_mst.tz) ephem_data = solarposition.pyephem(times, golden_mst.latitude, golden_mst.longitude, pressure=82000, @@ -408,8 +395,7 @@ def test_pyephem_physical(expected_solpos, golden_mst): @requires_ephem def test_pyephem_physical_dst(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), - periods=1, + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.pyephem(times, golden.latitude, golden.longitude, pressure=82000, @@ -425,7 +411,7 @@ def test_calc_time(): import math # validation from USNO solar position calculator online - epoch = datetime.datetime(1970, 1, 1) + epoch = datetime.datetime(1970,1,1) epoch_dt = pytz.utc.localize(epoch) loc = tus @@ -441,21 +427,21 @@ def test_calc_time(): actual_timestamp = (actual_time - epoch_dt).total_seconds() assert_allclose((alt.replace(second=0, microsecond=0) - - epoch_dt).total_seconds(), actual_timestamp) + epoch_dt).total_seconds(), actual_timestamp) assert_allclose((az.replace(second=0, microsecond=0) - - epoch_dt).total_seconds(), actual_timestamp) + epoch_dt).total_seconds(), actual_timestamp) @requires_ephem def test_earthsun_distance(): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, freq='D') distance = solarposition.pyephem_earthsun_distance(times).values[0] assert_allclose(1, distance, atol=0.1) def test_ephemeris_physical(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,12,30,30), periods=1, freq='D', tz=golden_mst.tz) ephem_data = solarposition.ephemeris(times, golden_mst.latitude, golden_mst.longitude, @@ -468,7 +454,7 @@ def test_ephemeris_physical(expected_solpos, golden_mst): def test_ephemeris_physical_dst(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.ephemeris(times, golden.latitude, golden.longitude, pressure=82000, @@ -480,7 +466,7 @@ def test_ephemeris_physical_dst(expected_solpos, golden): def test_ephemeris_physical_no_tz(expected_solpos, golden_mst): - times = pd.date_range(datetime.datetime(2003, 10, 17, 19, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,19,30,30), periods=1, freq='D') ephem_data = solarposition.ephemeris(times, golden_mst.latitude, golden_mst.longitude, @@ -493,7 +479,7 @@ def test_ephemeris_physical_no_tz(expected_solpos, golden_mst): def test_get_solarposition_error(golden): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, freq='D', tz=golden.tz) with pytest.raises(ValueError): ephem_data = solarposition.get_solarposition(times, golden.latitude, @@ -506,14 +492,14 @@ def test_get_solarposition_error(golden): @pytest.mark.parametrize("pressure, expected", [ (82000, _expected_solpos_df()), (90000, pd.DataFrame( - np.array([[39.88997, 50.11003, 194.34024, 39.87205, 14.64151, - 50.12795]]), + np.array([[ 39.88997, 50.11003, 194.34024, 39.87205, 14.64151, + 50.12795]]), columns=['apparent_elevation', 'apparent_zenith', 'azimuth', 'elevation', 'equation_of_time', 'zenith'], index=['2003-10-17T12:30:30Z'])) -]) + ]) def test_get_solarposition_pressure(pressure, expected, golden): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.get_solarposition(times, golden.latitude, golden.longitude, @@ -529,14 +515,14 @@ def test_get_solarposition_pressure(pressure, expected, golden): @pytest.mark.parametrize("altitude, expected", [ (1830.14, _expected_solpos_df()), (2000, pd.DataFrame( - np.array([[39.88788, 50.11212, 194.34024, 39.87205, 14.64151, - 50.12795]]), + np.array([[ 39.88788, 50.11212, 194.34024, 39.87205, 14.64151, + 50.12795]]), columns=['apparent_elevation', 'apparent_zenith', 'azimuth', 'elevation', 'equation_of_time', 'zenith'], index=['2003-10-17T12:30:30Z'])) -]) + ]) def test_get_solarposition_altitude(altitude, expected, golden): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.get_solarposition(times, golden.latitude, golden.longitude, @@ -557,10 +543,10 @@ def test_get_solarposition_altitude(altitude, expected, golden): marks=[pytest.mark.xfail( reason='spa.calculate_deltat not implemented for numba yet')]), (67.0, 'nrel_numba') -]) + ]) def test_get_solarposition_deltat(delta_t, method, expected_solpos_multi, golden): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=2, freq='D', tz=golden.tz) ephem_data = solarposition.get_solarposition(times, golden.latitude, golden.longitude, @@ -576,7 +562,7 @@ def test_get_solarposition_deltat(delta_t, method, expected_solpos_multi, def test_get_solarposition_no_kwargs(expected_solpos, golden): - times = pd.date_range(datetime.datetime(2003, 10, 17, 13, 30, 30), + times = pd.date_range(datetime.datetime(2003,10,17,13,30,30), periods=1, freq='D', tz=golden.tz) ephem_data = solarposition.get_solarposition(times, golden.latitude, golden.longitude) @@ -601,7 +587,7 @@ def test_get_solarposition_method_pyephem(expected_solpos, golden): def test_nrel_earthsun_distance(): times = pd.DatetimeIndex([datetime.datetime(2015, 1, 2), - datetime.datetime(2015, 8, 2), ] + datetime.datetime(2015, 8, 2),] ).tz_localize('MST') result = solarposition.nrel_earthsun_distance(times, delta_t=64.0) expected = pd.Series(np.array([0.983289204601, 1.01486146446]), @@ -690,27 +676,27 @@ def test_analytical_azimuth(): azimuth_2 = solarposition.solar_azimuth_analytical(lat_rad, hour_angle, decl, zenith) - idx = np.where(solar_zenith < np.pi / 2) + idx = np.where(solar_zenith < np.pi/2) assert np.allclose(azimuth_1[idx], solar_azimuth.values[idx], atol=0.01) assert np.allclose(azimuth_2[idx], solar_azimuth.values[idx], atol=0.017) # test for NaN values at boundary conditions (PR #431) test_angles = np.radians(np.array( - [[0., -180., -20.], - [0., 0., -5.], - [0., 0., 0.], - [0., 0., 15.], - [0., 180., 20.], - [30., 0., -20.], - [30., 0., -5.], - [30., 0., 0.], - [30., 180., 5.], - [30., 0., 10.], - [-30., 0., -20.], - [-30., 0., -15.], - [-30., 0., 0.], - [-30., -180., 5.], - [-30., 180., 10.]])) + [[ 0., -180., -20.], + [ 0., 0., -5.], + [ 0., 0., 0.], + [ 0., 0., 15.], + [ 0., 180., 20.], + [ 30., 0., -20.], + [ 30., 0., -5.], + [ 30., 0., 0.], + [ 30., 180., 5.], + [ 30., 0., 10.], + [ -30., 0., -20.], + [ -30., 0., -15.], + [ -30., 0., 0.], + [ -30., -180., 5.], + [ -30., 180., 10.]])) zeniths = solarposition.solar_zenith_analytical(*test_angles.T) azimuths = solarposition.solar_azimuth_analytical(*test_angles.T, @@ -790,82 +776,40 @@ def test_sun_rise_set_transit_geometric(expected_rise_set_spa, golden_mst): def test_spencer_mc(): """ test for the calculation based on spencer 1972 """ - times = pd.date_range('2018-01-01', periods=24, freq='H') latitude = 48.367073 longitude = 10.868378 - assert_frame_equal(solarposition.spencer_mc(times, latitude, longitude), - pd.DataFrame([[153.40940267, -63.40940267, 21.07562885, - 1.03506345, - -0.40158632, -3.18821393], - [147.78245317, -57.78245317, 47.01029634, - 1.03506418, - -0.40152783, -3.20679537], - [139.47329285, -49.47329285, 65.74488069, - 1.03506489, - -0.40146911, -3.22536965], - [129.96316706, -39.96316706, 79.89448997, - 1.03506558, - -0.40141017, -3.24393676], - [120.04134593, -30.04134593, 91.64276999, - 1.03506625, - -0.40135099, -3.26249664], - [110.16576717, -20.16576717, 102.30874903, - 1.0350669, - -0.40129158, -3.28104929], - [100.67550276, -10.67550276, 112.71473915, - 1.03506753, - -0.40123194, -3.29959466], - [91.89177037, -1.89177037, 123.43188389, - 1.03506814, - -0.40117207, -3.31813273], - [84.16861912, 5.83138088, 134.88256108, - 1.03506874, - -0.40111197, -3.33666346], - [77.91130514, 12.08869486, 147.34079438, - 1.03506931, - -0.40105163, -3.35518683], - [73.55352662, 16.44647338, 160.85682789, - 1.03506987, - -0.40099107, -3.37370281], - [71.47824875, 18.52175125, 175.16512108, - 1.03507041, - -0.40093028, -3.39221136], - [71.89802747, 18.10197253, 189.6987874, - 1.03507093, - -0.40086925, -3.41071246], - [74.76804069, 15.23195931, 203.79458854, - 1.03507142, - -0.400808, -3.42920608], - [79.80562508, 10.19437492, 216.97266043, - 1.0350719, - -0.40074651, -3.44769219], - [86.59854728, 3.40145272, 229.07943767, - 1.03507237, - -0.40068479, -3.46617075], - [94.71660269, -4.71660269, 240.25140046, - 1.03507281, - -0.40062285, -3.48464174], - [103.77123543, -13.77123543, 250.81773352, - 1.03507323, - -0.40056067, -3.50310512], - [113.42256954, -23.42256954, 261.24657665, - 1.03507363, - -0.40049826, -3.52156088], - [123.35015695, -33.35015695, 272.17648044, - 1.03507402, - -0.40043563, -3.54000897], - [133.18907815, -43.18907815, 284.56071242, - 1.03507438, - -0.40037276, -3.55844937], - [142.4005988, -52.4005988, 299.96695581, - 1.03507473, - -0.40030966, -3.57688204], - [150.00459208, -60.00459208, 320.87720374, - 1.03507506, - -0.40024633, -3.59530697], - [154.24856743, -64.24856743, 349.27789739, - 1.03507536, - -0.40018277, -3.61372411]], + assert_frame_equal(solarposition.spencer_mc(times[:10], + latitude, + longitude), + pd.DataFrame([[107.60190741, -17.60190741, 9.95454033, + 0.96706069, 0.40900622, -2.0908958], + [107.09533954, -17.09533954, 13.51934494, + 0.96705971, + 0.4090031, -2.09315672], + [106.43889684, -16.43889684, 17.04622292, + 0.96705874, + 0.40899997, -2.0954175], + [105.63682427, -15.63682427, 20.52704485, + 0.96705777, + 0.40899683, -2.09767815], + [104.69411613, -14.69411613, 23.95490053, + 0.9670568, + 0.40899367, -2.09993867], + [103.61638656, -13.61638656, 27.32420466, + 0.96705583, + 0.4089905, -2.10219905], + [102.40973635, -12.40973635, 30.63074363, + 0.96705486, + 0.40898732, -2.1044593], + [101.08062211, -11.08062211, 33.87167089, + 0.9670539, + 0.40898412, -2.10671942], + [99.63573414, -9.63573414, 37.04545803, + 0.96705293, + 0.40898091, -2.10897939], + [98.0818867, -8.0818867, 40.15181379, + 0.96705197, + 0.40897769, -2.11123924]], index=times, columns=['zenith', 'elevation', @@ -878,7 +822,6 @@ def test_spencer_mc(): def test_datetime2julians(): """ test transformation from datetime to julians """ - times = pd.date_range('2018-01-01', periods=24, freq='H') julians = solarposition.datetime2julian(pd.to_datetime(times)) np.testing.assert_array_almost_equal(julians, np.array([2458119.5, 2458119.54166667, From 5d464726f1f3fd12014421b116b94487c3aea93d Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Tue, 21 May 2019 10:01:29 -0700 Subject: [PATCH 04/13] Use existent declination function --- pvlib/solarposition.py | 18 ++++++++++-------- 1 file changed, 10 insertions(+), 8 deletions(-) diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index 8a4a609255..3c895d65a0 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -1463,7 +1463,7 @@ def sun_rise_set_transit_geometric(times, latitude, longitude, declination, def spencer_mc(times, latitude, longitude): """ Calculate the solar position using a python implementation of the - Spencer (1972) formulation + Spencer (1972) formulation provided by meteocontrol Parameters ---------- @@ -1506,15 +1506,17 @@ def spencer_mc(times, latitude, longitude): day_time = (julians_2000 % 1) * 24 # Eccentricity: correction factor of the earth's orbit. - eccentricity = (1.00011 + 0.034221 * cos_gamma[0] + 0.001280 - * sin_gamma[0] + 0.000719 * cos_gamma[1] + 0.000077 - * sin_gamma[1]) + eccentricity = ( + 1.00011 + + 0.034221 + * cos_gamma[0] + 0.001280 + * sin_gamma[0] + 0.000719 + * cos_gamma[1] + 0.000077 + * sin_gamma[1] + ) # declination. - declination = (0.006918 - 0.399912 * cos_gamma[0] - + 0.070257 * sin_gamma[0] - 0.006758 - * cos_gamma[1] + 0.000907 * sin_gamma[1] - - 0.002697 * cos_gamma[2] + 0.001480 * sin_gamma[2]) + declination = np.array(declination_spencer71(times.dayofyear)) # Equation of time (difference between standard time and solar time). eot = (0.000075 + 0.001868 * cos_gamma[0] - 0.032077 * sin_gamma[0] From 6ae7dae22c7d9bea5951f3f66c87883ced30ae6a Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Tue, 21 May 2019 10:07:04 -0700 Subject: [PATCH 05/13] Documentation stuff --- docs/sphinx/source/api.rst | 2 ++ docs/sphinx/source/whatsnew/v0.7.1.rst | 11 +++++++++++ 2 files changed, 13 insertions(+) create mode 100644 docs/sphinx/source/whatsnew/v0.7.1.rst diff --git a/docs/sphinx/source/api.rst b/docs/sphinx/source/api.rst index 45ef9ec50a..640159383c 100644 --- a/docs/sphinx/source/api.rst +++ b/docs/sphinx/source/api.rst @@ -45,6 +45,7 @@ unless you know that you need a different function. solarposition.ephemeris solarposition.pyephem solarposition.spa_c + solarposition.spencer_mc Additional functions for quantities closely related to solar position. @@ -55,6 +56,7 @@ Additional functions for quantities closely related to solar position. solarposition.calc_time solarposition.pyephem_earthsun_distance solarposition.nrel_earthsun_distance + solarposition.datetime2julian spa.calculate_deltat diff --git a/docs/sphinx/source/whatsnew/v0.7.1.rst b/docs/sphinx/source/whatsnew/v0.7.1.rst new file mode 100644 index 0000000000..fe46d5da89 --- /dev/null +++ b/docs/sphinx/source/whatsnew/v0.7.1.rst @@ -0,0 +1,11 @@ +.. _whatsnew_0701: + +v0.7.1 (MONTH DAY, YEAR) +--------------------- + +- Add a faster way to calculate the solar position (`spencer_mc`) + + +Contributors +~~~~~~~~~~~~ +* Daniel Lassahn (:ghuser:`meteoDaniel`) \ No newline at end of file From a43cf78f255a42402a76b2a93d9690c8af8e8e93 Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Tue, 21 May 2019 10:40:04 -0700 Subject: [PATCH 06/13] Tidy --- pvlib/solarposition.py | 34 ++++++++++++++-------------------- pvlib/spa.py | 2 +- 2 files changed, 15 insertions(+), 21 deletions(-) diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index 3c895d65a0..c35de620b2 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -122,8 +122,7 @@ def get_solarposition(time, latitude, longitude, pressure=pressure, temperature=temperature, **kwargs) elif method == 'ephemeris': - ephem_df = ephemeris(time, latitude, longitude, pressure, temperature, - **kwargs) + ephem_df = ephemeris(time, latitude, longitude, pressure, temperature) elif method == 'spencer_mc': ephem_df = spencer_mc(time, latitude, longitude) @@ -508,7 +507,7 @@ def sun_rise_set_transit_ephem(times, latitude, longitude, Parameters ---------- - time : pandas.DatetimeIndex + times : :class:`pandas.DatetimeIndex` Must be localized latitude : float Latitude in degrees, positive north of equator, negative to south @@ -1467,7 +1466,7 @@ def spencer_mc(times, latitude, longitude): Parameters ---------- - times : pandas.DatetimeIndex + times : :class:`pandas.DatetimeIndex` Corresponding timestamps, must be localized to the timezone for the ``latitude`` and ``longitude``. latitude : float @@ -1488,7 +1487,7 @@ def spencer_mc(times, latitude, longitude): References ---------- - [4] Spencer (1972) can be found in + [1] Spencer (1972) can be found in Sen, Zekai. Solar energy fundamentals and modeling techniques: atmosphere, environment, climate change and renewable energy. Springer Science & Business Media, 2008. @@ -1508,14 +1507,12 @@ def spencer_mc(times, latitude, longitude): # Eccentricity: correction factor of the earth's orbit. eccentricity = ( 1.00011 + - 0.034221 - * cos_gamma[0] + 0.001280 - * sin_gamma[0] + 0.000719 - * cos_gamma[1] + 0.000077 - * sin_gamma[1] + 0.034221 * cos_gamma[0] + + 0.001280 * sin_gamma[0] + + 0.000719 * cos_gamma[1] + + 0.000077 * sin_gamma[1] ) - # declination. declination = np.array(declination_spencer71(times.dayofyear)) # Equation of time (difference between standard time and solar time). @@ -1529,10 +1526,8 @@ def spencer_mc(times, latitude, longitude): ha = np.radians(tlt * 15) # Calculate sun elevation. - sin_sun_elevation = ( - np.sin(declination) * np.sin(lat) + np.cos(declination) - * np.cos(lat) * np.cos(ha) - ) + sin_sun_elevation = (np.sin(declination) * np.sin(lat) + + np.cos(declination) * np.cos(lat) * np.cos(ha)) # Compute the sun's elevation and zenith angle. elevation = np.arcsin(sin_sun_elevation) @@ -1563,10 +1558,9 @@ def datetime2julian(times): Parameters ---------- - times : pandas.DatetimeIndex + times : :class:`pandas.DatetimeIndex` Corresponding timestamps, must be localized to the timezone for the - ``latitude`` and ``longitude``. - + ``latitude`` and ``longitude`` Returns ------- Float64Index @@ -1575,6 +1569,6 @@ def datetime2julian(times): delta = times - DT_2000 return ( - JULIAN_2000 + - delta.days + (delta.seconds + delta.microseconds / 1e6) / DAY_SECONDS + JULIAN_2000 + delta.days + + (delta.seconds + delta.microseconds / 1e6) / DAY_SECONDS ) diff --git a/pvlib/spa.py b/pvlib/spa.py index 3a44a2efed..93e111e34d 100644 --- a/pvlib/spa.py +++ b/pvlib/spa.py @@ -1161,7 +1161,7 @@ def transit_sunrise_sunset(dates, lat, lon, delta_t, numthreads): Parameters ---------- - dates : array + dates : np.array Numpy array of ints/floats corresponding to the Unix time for the dates of interest, must be midnight UTC (00:00+00:00) on the day of interest. From 1c1b2c6b1e0ee0d5a6cdf415bdb400b2917ac0a3 Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Tue, 21 May 2019 14:02:54 -0700 Subject: [PATCH 07/13] Adapt tests and place datetime_to_julian in tools --- pvlib/solarposition.py | 29 +---- pvlib/test/test_solarposition.py | 184 +++++++++++++++++++------------ pvlib/test/test_tools.py | 22 ++++ pvlib/tools.py | 26 +++++ 4 files changed, 165 insertions(+), 96 deletions(-) diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index c35de620b2..4f17d49b53 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -24,14 +24,11 @@ import warnings from pvlib import atmosphere -from pvlib.tools import datetime_to_djd, djd_to_datetime +from pvlib.tools import datetime_to_djd, djd_to_datetime, datetime_to_julian from pvlib._deprecation import deprecated NS_PER_HR = 1.e9 * 3600. # nanoseconds per hour JULIAN_2000 = 2451544.5 -DT_2000 = dt.datetime(2000, 1, 1) -TIMESTAMP_2000 = DT_2000.timestamp() -DAY_SECONDS = 60 * 60 * 24 JULIAN_YEARS = 365.2425 @@ -1493,7 +1490,7 @@ def spencer_mc(times, latitude, longitude): Springer Science & Business Media, 2008. """ - julians = datetime2julian(times) + julians = datetime_to_julian(times) julians_2000 = np.asarray(julians, dtype=np.float) - JULIAN_2000 lat = np.radians(latitude) @@ -1550,25 +1547,3 @@ def spencer_mc(times, latitude, longitude): 'equation_of_time': eot}, index=times) return result - - -def datetime2julian(times): - """ - Transforms pd.DateTimeIndex to Julian days - - Parameters - ---------- - times : :class:`pandas.DatetimeIndex` - Corresponding timestamps, must be localized to the timezone for the - ``latitude`` and ``longitude`` - Returns - ------- - Float64Index - The float index contains julian dates - """ - - delta = times - DT_2000 - return ( - JULIAN_2000 + delta.days + - (delta.seconds + delta.microseconds / 1e6) / DAY_SECONDS - ) diff --git a/pvlib/test/test_solarposition.py b/pvlib/test/test_solarposition.py index 1a1bdbf34a..95521355bd 100644 --- a/pvlib/test/test_solarposition.py +++ b/pvlib/test/test_solarposition.py @@ -753,80 +753,126 @@ def test_sun_rise_set_transit_geometric(expected_rise_set_spa, golden_mst): assert np.allclose(test_transit, expected_transit, atol=np.abs(expected_transit_error).max()) - + def test_spencer_mc(): """ test for the calculation based on spencer 1972 """ latitude = 48.367073 longitude = 10.868378 - assert_frame_equal(solarposition.spencer_mc(times[:10], - latitude, - longitude), - pd.DataFrame([[107.60190741, -17.60190741, 9.95454033, - 0.96706069, 0.40900622, -2.0908958], - [107.09533954, -17.09533954, 13.51934494, - 0.96705971, - 0.4090031, -2.09315672], - [106.43889684, -16.43889684, 17.04622292, - 0.96705874, - 0.40899997, -2.0954175], - [105.63682427, -15.63682427, 20.52704485, - 0.96705777, - 0.40899683, -2.09767815], - [104.69411613, -14.69411613, 23.95490053, - 0.9670568, - 0.40899367, -2.09993867], - [103.61638656, -13.61638656, 27.32420466, - 0.96705583, - 0.4089905, -2.10219905], - [102.40973635, -12.40973635, 30.63074363, - 0.96705486, - 0.40898732, -2.1044593], - [101.08062211, -11.08062211, 33.87167089, - 0.9670539, - 0.40898412, -2.10671942], - [99.63573414, -9.63573414, 37.04545803, - 0.96705293, - 0.40898091, -2.10897939], - [98.0818867, -8.0818867, 40.15181379, - 0.96705197, - 0.40897769, -2.11123924]], - index=times, - columns=['zenith', - 'elevation', - 'azimuth', - 'eccentricity', - 'declination', - 'equation_of_time']) - ) - -def test_datetime2julians(): - """ test transformation from datetime to julians """ - julians = solarposition.datetime2julian(pd.to_datetime(times)) - np.testing.assert_array_almost_equal(julians, - np.array([2458119.5, 2458119.54166667, - 2458119.58333333, - 2458119.625, - 2458119.66666667, - 2458119.70833333, - 2458119.75, 2458119.79166667, - 2458119.83333333, - 2458119.875, - 2458119.91666667, - 2458119.95833333, - 2458120., 2458120.04166667, - 2458120.08333333, - 2458120.125, - 2458120.16666667, - 2458120.20833333, - 2458120.25, 2458120.29166667, - 2458120.33333333, - 2458120.375, - 2458120.41666667, - 2458120.45833333]) - ) - -# put numba tests at end of file to minimize reloading + to_test = solarposition.spencer_mc(times, + latitude, + longitude) + np.testing.assert_array_almost_equal( + to_test.zenith.values, + np.array( + [107.5944161, 107.08772172, 106.43116829, 105.62900042, + 104.68621177, 103.60841558, 102.40171152, 101.07255494, + 99.62763476, 98.07376382, 96.41778336, 94.66648432, + 92.82654495, 90.90448339, 88.90662551, 86.83908587, + 84.70775942, 82.51832403, 80.27625186, 77.98682746, + 75.65517374, 73.28628434, 70.88506165, 68.45636234, + 66.0050506, 63.53605921, 61.05446218, 58.56556086, + 56.07498558, 53.58881946, 51.11374884, 48.65724635, + 46.22779685, 43.83517458, 41.49078009, 39.20804631, + 37.00291337, 34.89435623, 32.90492124, 31.06117096, + 29.39385375, 27.93751303, 26.72916357, 25.80569641, + 25.19996401, 24.93609994, 25.02533181, 25.46381181, + 26.23337616, 27.30490594, 28.64295418, 30.21015627, + 31.97050468, 33.89127116, 35.94381137, 38.10363121, + 40.35004899, 42.66568743, 45.03593362, 47.44843417, + 49.89265446, 52.35950849, 54.8410528, 57.33023629, + 59.8206975, 62.30659851, 64.78248863, 67.24319292, + 69.68371835, 72.09917489, 74.4847094, 76.83544826, + 79.14644817, 81.41265486, 83.62886755, 85.78970975, + 87.88960763, 89.92277472, 91.88320492, 93.76467562, + 95.5607608, 97.26485655, 98.8702211, 100.37002926, + 101.7574431, 103.02569948, 104.16821258, 105.17869036, + 106.05126175, 106.78060921, 107.36210092, 107.79191584, + 108.06715341, 108.18592116, 108.14739387, 107.9518398, + 107.61767954, 107.11306117, 106.45849746, 105.65822415, + 104.71722821, 103.64111837, 102.43599187, 101.10830346, + 99.66474308, 98.11212592, 96.45729665, 94.7070504, + 92.86807018, 90.9468792, 88.94980851, 86.88297787, + 84.75228725, 82.56341928, 80.32185056, 78.03286967, + 75.70160305, 73.33304741, 70.93210757, 68.50364202, + 66.05251605, 63.58366268, 61.10215525, 58.61329329, + 56.12270393, 53.63646542, 51.16125711, 48.70454205, + 46.27479211, 43.88176429, 41.53683649, 39.25341219, + 37.04739373, 34.93770849, 32.94684449, 31.10129575, + 29.43173548, 27.97263389, 26.7609499, 25.83355775, + 25.22335515, 24.95459651, 25.03869819, 25.47203708, + 26.23666504, 27.30363156, 28.63759291, 30.20122571, + 31.95851686, 33.87670365, 35.92709343, 38.08514007, + 40.33011209, 42.64458753, 45.01391492, 47.42570839, + 49.86940645, 52.33590108, 54.81723096, 57.30633064, + 59.79682723, 62.2828738, 64.75901269, 67.22006367, + 69.66102988, 72.0770186, 74.46317504, 76.81462477, + 79.12642433, 81.39351992, 83.61071161, 85.77262406, + 87.87368475, 89.90810859, 91.86989066, 93.75280924, + 95.55043872, 97.25617485, 98.86327462, 100.36491047, + 101.75424074, 103.024497, 104.16908649, 105.1817084, + 106.05648115, 106.78807492, 107.3718441, 107.80395262, + 108.08148404, 108.20252959, 108.1662479, 107.97289173, + 107.64770976] + ) + ) + + np.testing.assert_array_almost_equal( + to_test.azimuth.values, + np.array( + [9.95354884, 13.51797616, 17.04447274, 20.52491225, + 23.95238719, 27.32131453, 30.62748235, 33.86804529, + 37.04147561, 40.14748222, 43.18691015, 46.1616283, + 49.07441566, 51.92885441, 54.72923301, 57.48046422, + 60.18802201, 62.85789652, 65.49656941, 68.11101133, + 70.70870005, 73.29766135, 75.88653487, 78.48466488, + 81.10221958, 83.75034334, 86.44134434, 89.18892358, + 92.00845226, 94.91730112, 97.93522718, 101.08482003, + 104.39199808, 107.88653257, 111.60254901, 115.57890087, + 119.8592317, 124.49141579, 129.52588659, 135.01216661, + 140.99278628, 147.49394635, 154.51314061, 162.00583842, + 169.87594269, 177.97643578, 186.12466927, 194.13016535, + 201.82539566, 209.08806319, 215.84877252, 222.08574604, + 227.81263221, 233.06524098, 237.89058053, 242.3392101, + 246.46057573, 250.30054442, 253.90036735, 257.29648187, + 260.52076061, 263.60096725, 266.56127763, 269.42279422, + 272.20401931, 274.92127106, 277.58904086, 280.22029634, + 282.82673312, 285.41898225, 288.00677995, 290.59910223, + 293.20426943, 295.83002468, 298.48358643, 301.17167776, + 303.9005345, 306.67589094, 309.50294481, 312.38630353, + 315.32991074, 318.33695647, 321.40977511, 324.54973374, + 327.75711805, 331.0310254, 334.36927221, 337.76832766, + 341.22328647, 344.72788848, 348.27459447, 351.85472406, + 355.45865247, 359.07606104, 2.69623101, 6.30836085, + 9.90413618, 13.46983714, 16.99770579, 20.47959425, + 23.90857277, 27.27903624, 30.58675167, 33.82885422, + 37.00379863, 40.1112785, 43.15212578, 46.12819847, + 49.04226661, 51.89790514, 54.69939668, 57.45164925, + 60.1601328, 62.83083386, 65.47023059, 68.08528982, + 70.68348498, 73.27283649, 75.86197737, 78.46024348, + 81.07779237, 83.72575502, 86.41642271, 89.16347527, + 91.9822575, 94.8901073, 97.90674111, 101.05469855, + 104.35983703, 107.85185424, 111.56478847, 115.53739258, + 119.81319937, 124.43997055, 129.46804248, 134.94688671, + 140.91907623, 147.41101064, 154.42060141, 161.90398226, + 169.76591074, 177.86024433, 186.00498616, 194.00988407, + 201.70716439, 208.97394261, 215.74010434, 221.98320309, + 227.71636618, 232.97505067, 237.80605645, 242.2598406, + 246.38581686, 250.22986166, 253.83325807, 257.23248545, + 260.45946126, 263.54199277, 266.50429617, 269.36751012, + 272.15016893, 274.86861866, 277.53737498, 280.1694266, + 282.77648747, 285.3692047, 287.95732871, 290.54984827, + 293.15509539, 295.78082413, 298.43426346, 301.12214686, + 303.85072078, 306.62573055, 309.45238569, 312.3353063, + 315.27844985, 318.28502138, 321.35737159, 324.49688503, + 327.70386592, 330.97743082, 334.31541557, 337.71430836, + 341.1692218, 344.67391164, 348.22085175, 351.80137124, + 355.40585053, 359.02397146, 2.64501087, 6.2581581, + 9.85597076] + ) + ) + +# put numba tests at end of file to minimize reloading @requires_numba def test_spa_python_numba_physical(expected_solpos, golden_mst): times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), diff --git a/pvlib/test/test_tools.py b/pvlib/test/test_tools.py index 42eaedca58..b88372fa32 100644 --- a/pvlib/test/test_tools.py +++ b/pvlib/test/test_tools.py @@ -1,6 +1,10 @@ +import numpy as np +import pandas as pd import pytest +import pvlib.tools from pvlib import tools +from pvlib.test.test_solarposition import times @pytest.mark.parametrize('keys, input_dict, expected', [ @@ -12,3 +16,21 @@ def test_build_kwargs(keys, input_dict, expected): kwargs = tools._build_kwargs(keys, input_dict) assert kwargs == expected + + +def test_datetime_to_julian(): + """ test transformation from datetime to julians """ + julians = pvlib.tools.datetime_to_julian(pd.to_datetime(times)) + np.testing.assert_array_almost_equal(np.array(julians[:10]), + np.array([ + 2456832.5, + 2456832.5104166665, + 2456832.5208333335, + 2456832.53125, + 2456832.5416666665, + 2456832.5520833335, + 2456832.5625, + 2456832.5729166665, + 2456832.5833333335, + 2456832.59375]) + ) diff --git a/pvlib/tools.py b/pvlib/tools.py index cbd59a8541..1ce9186fea 100644 --- a/pvlib/tools.py +++ b/pvlib/tools.py @@ -9,6 +9,10 @@ import pandas as pd import pytz +JULIAN_2000 = 2451544.5 +DT_2000 = dt.datetime(2000, 1, 1) +DAY_SECONDS = 60 * 60 * 24 + def cosd(angle): """ @@ -425,3 +429,25 @@ def _golden_sect_DataFrame(params, VL, VH, func): raise Exception("EXCEPTION:iterations exceeded maximum (50)") return func(df, 'V1'), df['V1'] + + +def datetime_to_julian(times): + """ + Transforms pd.DateTimeIndex to Julian days + + Parameters + ---------- + times : :class:`pandas.DatetimeIndex` + Corresponding timestamps, must be localized to the timezone for the + ``latitude`` and ``longitude`` + Returns + ------- + Float64Index + The float index contains julian dates + """ + + delta = times - DT_2000 + return ( + JULIAN_2000 + delta.days + + (delta.seconds + delta.microseconds / 1e6) / DAY_SECONDS + ) From 0749c55e7da0b9ce35b26707185f0874c5f34c87 Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Tue, 21 May 2019 15:24:57 -0700 Subject: [PATCH 08/13] Linter optimisation --- pvlib/solarposition.py | 14 +++++--------- 1 file changed, 5 insertions(+), 9 deletions(-) diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index 4f17d49b53..5d5c9ad340 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -1502,19 +1502,15 @@ def spencer_mc(times, latitude, longitude): day_time = (julians_2000 % 1) * 24 # Eccentricity: correction factor of the earth's orbit. - eccentricity = ( - 1.00011 + - 0.034221 * cos_gamma[0] + - 0.001280 * sin_gamma[0] + - 0.000719 * cos_gamma[1] + - 0.000077 * sin_gamma[1] - ) + eccentricity = (1.00011 + 0.034221 * cos_gamma[0] + 0.001280 * + sin_gamma[0] + 0.000719 * cos_gamma[1] + 0.000077 * + sin_gamma[1]) declination = np.array(declination_spencer71(times.dayofyear)) # Equation of time (difference between standard time and solar time). - eot = (0.000075 + 0.001868 * cos_gamma[0] - 0.032077 * sin_gamma[0] - - 0.014615 * cos_gamma[1] - 0.040849 * sin_gamma[1]) * 229.18 + eot = (0.000075 + 0.001868 * cos_gamma[0] - 0.032077 * sin_gamma[0] - + 0.014615 * cos_gamma[1] - 0.040849 * sin_gamma[1]) * 229.18 # True local time tlt = (day_time + longitude / 15 + eot / 60) % 24 - 12 From bf88e3c03c52d5b3e954935acaa37eb504a699dd Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Tue, 21 May 2019 20:36:35 -0700 Subject: [PATCH 09/13] Linter optimisation --- pvlib/solarposition.py | 14 +++++++++----- pvlib/test/test_solarposition.py | 4 ++-- 2 files changed, 11 insertions(+), 7 deletions(-) diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index 5d5c9ad340..3e7593f651 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -1502,15 +1502,19 @@ def spencer_mc(times, latitude, longitude): day_time = (julians_2000 % 1) * 24 # Eccentricity: correction factor of the earth's orbit. - eccentricity = (1.00011 + 0.034221 * cos_gamma[0] + 0.001280 * - sin_gamma[0] + 0.000719 * cos_gamma[1] + 0.000077 * - sin_gamma[1]) + eccentricity = 1.00011 + 0.034221 * cos_gamma[0] + eccentricity += 0.001280 * sin_gamma[0] + eccentricity += 0.000719 * cos_gamma[1] + eccentricity += 0.000077 * sin_gamma[1] declination = np.array(declination_spencer71(times.dayofyear)) # Equation of time (difference between standard time and solar time). - eot = (0.000075 + 0.001868 * cos_gamma[0] - 0.032077 * sin_gamma[0] - - 0.014615 * cos_gamma[1] - 0.040849 * sin_gamma[1]) * 229.18 + eot = 0.000075 + 0.001868 * cos_gamma[0] + eot -= 0.032077 * sin_gamma[0] + eot -= 0.014615 * cos_gamma[1] + eot -= 0.040849 * sin_gamma[1] + eot *= 229.18 # True local time tlt = (day_time + longitude / 15 + eot / 60) % 24 - 12 diff --git a/pvlib/test/test_solarposition.py b/pvlib/test/test_solarposition.py index 95521355bd..dacdd4e1b8 100644 --- a/pvlib/test/test_solarposition.py +++ b/pvlib/test/test_solarposition.py @@ -759,8 +759,8 @@ def test_spencer_mc(): latitude = 48.367073 longitude = 10.868378 to_test = solarposition.spencer_mc(times, - latitude, - longitude) + latitude, + longitude) np.testing.assert_array_almost_equal( to_test.zenith.values, np.array( From b20ca5319fd2ea01b0261401d112b70f49146604 Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Tue, 21 May 2019 20:49:24 -0700 Subject: [PATCH 10/13] Linter optimisation --- pvlib/solarposition.py | 4 ++-- pvlib/tools.py | 4 ++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index 3e7593f651..522764dc3c 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -1523,8 +1523,8 @@ def spencer_mc(times, latitude, longitude): ha = np.radians(tlt * 15) # Calculate sun elevation. - sin_sun_elevation = (np.sin(declination) * np.sin(lat) + - np.cos(declination) * np.cos(lat) * np.cos(ha)) + sin_sun_elevation = np.sin(declination) * np.sin(lat) + sin_sun_elevation += np.cos(declination) * np.cos(lat) * np.cos(ha) # Compute the sun's elevation and zenith angle. elevation = np.arcsin(sin_sun_elevation) diff --git a/pvlib/tools.py b/pvlib/tools.py index 1ce9186fea..998055c36d 100644 --- a/pvlib/tools.py +++ b/pvlib/tools.py @@ -447,7 +447,7 @@ def datetime_to_julian(times): """ delta = times - DT_2000 + delta_julians = (delta.seconds + delta.microseconds / 1e6) return ( - JULIAN_2000 + delta.days + - (delta.seconds + delta.microseconds / 1e6) / DAY_SECONDS + JULIAN_2000 + delta.days + delta_julians / DAY_SECONDS ) From fe5bf63cdca8c84e58a6c1f091343acd74f1b2ee Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Wed, 22 May 2019 18:34:56 -0700 Subject: [PATCH 11/13] Rename function, use existend equation of time function, whatsnew in 7.0.0. --- docs/sphinx/source/api.rst | 2 +- docs/sphinx/source/whatsnew/v0.7.0.rst | 3 + docs/sphinx/source/whatsnew/v0.7.1.rst | 11 -- pvlib/solarposition.py | 44 ++---- pvlib/test/test_solarposition.py | 208 ++++++++++++------------- pvlib/test/test_tools.py | 6 +- 6 files changed, 127 insertions(+), 147 deletions(-) delete mode 100644 docs/sphinx/source/whatsnew/v0.7.1.rst diff --git a/docs/sphinx/source/api.rst b/docs/sphinx/source/api.rst index 640159383c..22e84d10a0 100644 --- a/docs/sphinx/source/api.rst +++ b/docs/sphinx/source/api.rst @@ -45,7 +45,7 @@ unless you know that you need a different function. solarposition.ephemeris solarposition.pyephem solarposition.spa_c - solarposition.spencer_mc + solarposition.spencer Additional functions for quantities closely related to solar position. diff --git a/docs/sphinx/source/whatsnew/v0.7.0.rst b/docs/sphinx/source/whatsnew/v0.7.0.rst index cf2adfc9e6..a5c6f52e48 100644 --- a/docs/sphinx/source/whatsnew/v0.7.0.rst +++ b/docs/sphinx/source/whatsnew/v0.7.0.rst @@ -8,8 +8,11 @@ recommend all users of v0.6.3 upgrade to this release. **Python 2.7 support ended on June 1, 2019**. (:issue:`501`) +- Add a faster way to calculate the solar position (`spencer`) + Contributors ~~~~~~~~~~~~ * Mark Campanellli (:ghuser:`markcampanelli`) * Will Holmgren (:ghuser:`wholmgren`) +* Daniel Lassahn (:ghuser:`meteoDaniel`) \ No newline at end of file diff --git a/docs/sphinx/source/whatsnew/v0.7.1.rst b/docs/sphinx/source/whatsnew/v0.7.1.rst deleted file mode 100644 index fe46d5da89..0000000000 --- a/docs/sphinx/source/whatsnew/v0.7.1.rst +++ /dev/null @@ -1,11 +0,0 @@ -.. _whatsnew_0701: - -v0.7.1 (MONTH DAY, YEAR) ---------------------- - -- Add a faster way to calculate the solar position (`spencer_mc`) - - -Contributors -~~~~~~~~~~~~ -* Daniel Lassahn (:ghuser:`meteoDaniel`) \ No newline at end of file diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index 522764dc3c..413fc7a212 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -68,7 +68,7 @@ def get_solarposition(time, latitude, longitude, 'nrel_c' uses the NREL SPA C code [3]: :py:func:`spa_c` - 'spencer_mc' uses the Spencer formula [4] :py:func:`spencer_mc` + 'spencer' uses the Spencer formula [4] :py:func:`spencer` temperature : float, default 12 @@ -120,8 +120,8 @@ def get_solarposition(time, latitude, longitude, temperature=temperature, **kwargs) elif method == 'ephemeris': ephem_df = ephemeris(time, latitude, longitude, pressure, temperature) - elif method == 'spencer_mc': - ephem_df = spencer_mc(time, latitude, longitude) + elif method == 'spencer': + ephem_df = spencer(time, latitude, longitude) else: raise ValueError('Invalid solar position method') @@ -1456,14 +1456,13 @@ def sun_rise_set_transit_geometric(times, latitude, longitude, declination, return sunrise, sunset, transit -def spencer_mc(times, latitude, longitude): +def spencer(times, latitude, longitude): """ - Calculate the solar position using a python implementation of the - Spencer (1972) formulation provided by meteocontrol + Calculate the solar position using a formulation by Spencer 1971/1972. Parameters ---------- - times : :class:`pandas.DatetimeIndex` + times : pandas.DatetimeIndex` Corresponding timestamps, must be localized to the timezone for the ``latitude`` and ``longitude``. latitude : float @@ -1479,7 +1478,6 @@ def spencer_mc(times, latitude, longitude): elevation (degrees), azimuth (degrees), equation_of_time (seconds), - eccentricity, declination (degrees). References @@ -1493,28 +1491,13 @@ def spencer_mc(times, latitude, longitude): julians = datetime_to_julian(times) julians_2000 = np.asarray(julians, dtype=np.float) - JULIAN_2000 - lat = np.radians(latitude) - - # Compute fractional year (gamma) in radians - gamma = 2 * np.pi * (julians_2000 % JULIAN_YEARS) / JULIAN_YEARS - cos_gamma = np.cos(gamma), np.cos(gamma * 2), np.cos(gamma * 3) - sin_gamma = np.sin(gamma), np.sin(gamma * 2), np.sin(gamma * 3) + latitude_radians = np.radians(latitude) day_time = (julians_2000 % 1) * 24 - # Eccentricity: correction factor of the earth's orbit. - eccentricity = 1.00011 + 0.034221 * cos_gamma[0] - eccentricity += 0.001280 * sin_gamma[0] - eccentricity += 0.000719 * cos_gamma[1] - eccentricity += 0.000077 * sin_gamma[1] - declination = np.array(declination_spencer71(times.dayofyear)) # Equation of time (difference between standard time and solar time). - eot = 0.000075 + 0.001868 * cos_gamma[0] - eot -= 0.032077 * sin_gamma[0] - eot -= 0.014615 * cos_gamma[1] - eot -= 0.040849 * sin_gamma[1] - eot *= 229.18 + eot = np.array(equation_of_time_spencer71(times.dayofyear)) # True local time tlt = (day_time + longitude / 15 + eot / 60) % 24 - 12 @@ -1523,16 +1506,18 @@ def spencer_mc(times, latitude, longitude): ha = np.radians(tlt * 15) # Calculate sun elevation. - sin_sun_elevation = np.sin(declination) * np.sin(lat) - sin_sun_elevation += np.cos(declination) * np.cos(lat) * np.cos(ha) + sin_sun_elevation = ( + np.sin(declination) * np.sin(latitude_radians) + + np.cos(declination) * np.cos(latitude_radians) * np.cos(ha) + ) # Compute the sun's elevation and zenith angle. elevation = np.arcsin(sin_sun_elevation) zenith = np.pi / 2 - elevation # Compute the sun's azimuth angle. - y = -(np.sin(lat) * np.sin(elevation) - np.sin(declination)) \ - / (np.cos(lat) * np.cos(elevation)) + y = -(np.sin(latitude_radians) * np.sin(elevation) - np.sin(declination)) \ + / (np.cos(latitude_radians) * np.cos(elevation)) azimuth = np.arccos(y) # Convert azimuth angle from 0-pi to 0-2pi. @@ -1542,7 +1527,6 @@ def spencer_mc(times, latitude, longitude): result = pd.DataFrame({'zenith': np.degrees(zenith), 'elevation': np.degrees(elevation), 'azimuth': np.degrees(azimuth), - 'eccentricity': eccentricity, 'declination': declination, 'equation_of_time': eot}, index=times) diff --git a/pvlib/test/test_solarposition.py b/pvlib/test/test_solarposition.py index dacdd4e1b8..d93df2046f 100644 --- a/pvlib/test/test_solarposition.py +++ b/pvlib/test/test_solarposition.py @@ -754,120 +754,120 @@ def test_sun_rise_set_transit_geometric(expected_rise_set_spa, golden_mst): atol=np.abs(expected_transit_error).max()) -def test_spencer_mc(): +def test_spencer(): """ test for the calculation based on spencer 1972 """ - latitude = 48.367073 - longitude = 10.868378 - to_test = solarposition.spencer_mc(times, - latitude, - longitude) + latitude = 32.2 + longitude = -111 + to_test = solarposition.spencer(times, + latitude, + longitude) np.testing.assert_array_almost_equal( to_test.zenith.values, np.array( - [107.5944161, 107.08772172, 106.43116829, 105.62900042, - 104.68621177, 103.60841558, 102.40171152, 101.07255494, - 99.62763476, 98.07376382, 96.41778336, 94.66648432, - 92.82654495, 90.90448339, 88.90662551, 86.83908587, - 84.70775942, 82.51832403, 80.27625186, 77.98682746, - 75.65517374, 73.28628434, 70.88506165, 68.45636234, - 66.0050506, 63.53605921, 61.05446218, 58.56556086, - 56.07498558, 53.58881946, 51.11374884, 48.65724635, - 46.22779685, 43.83517458, 41.49078009, 39.20804631, - 37.00291337, 34.89435623, 32.90492124, 31.06117096, - 29.39385375, 27.93751303, 26.72916357, 25.80569641, - 25.19996401, 24.93609994, 25.02533181, 25.46381181, - 26.23337616, 27.30490594, 28.64295418, 30.21015627, - 31.97050468, 33.89127116, 35.94381137, 38.10363121, - 40.35004899, 42.66568743, 45.03593362, 47.44843417, - 49.89265446, 52.35950849, 54.8410528, 57.33023629, - 59.8206975, 62.30659851, 64.78248863, 67.24319292, - 69.68371835, 72.09917489, 74.4847094, 76.83544826, - 79.14644817, 81.41265486, 83.62886755, 85.78970975, - 87.88960763, 89.92277472, 91.88320492, 93.76467562, - 95.5607608, 97.26485655, 98.8702211, 100.37002926, - 101.7574431, 103.02569948, 104.16821258, 105.17869036, - 106.05126175, 106.78060921, 107.36210092, 107.79191584, - 108.06715341, 108.18592116, 108.14739387, 107.9518398, - 107.61767954, 107.11306117, 106.45849746, 105.65822415, - 104.71722821, 103.64111837, 102.43599187, 101.10830346, - 99.66474308, 98.11212592, 96.45729665, 94.7070504, - 92.86807018, 90.9468792, 88.94980851, 86.88297787, - 84.75228725, 82.56341928, 80.32185056, 78.03286967, - 75.70160305, 73.33304741, 70.93210757, 68.50364202, - 66.05251605, 63.58366268, 61.10215525, 58.61329329, - 56.12270393, 53.63646542, 51.16125711, 48.70454205, - 46.27479211, 43.88176429, 41.53683649, 39.25341219, - 37.04739373, 34.93770849, 32.94684449, 31.10129575, - 29.43173548, 27.97263389, 26.7609499, 25.83355775, - 25.22335515, 24.95459651, 25.03869819, 25.47203708, - 26.23666504, 27.30363156, 28.63759291, 30.20122571, - 31.95851686, 33.87670365, 35.92709343, 38.08514007, - 40.33011209, 42.64458753, 45.01391492, 47.42570839, - 49.86940645, 52.33590108, 54.81723096, 57.30633064, - 59.79682723, 62.2828738, 64.75901269, 67.22006367, - 69.66102988, 72.0770186, 74.46317504, 76.81462477, - 79.12642433, 81.39351992, 83.61071161, 85.77262406, - 87.87368475, 89.90810859, 91.86989066, 93.75280924, - 95.55043872, 97.25617485, 98.86327462, 100.36491047, - 101.75424074, 103.024497, 104.16908649, 105.1817084, - 106.05648115, 106.78807492, 107.3718441, 107.80395262, - 108.08148404, 108.20252959, 108.1662479, 107.97289173, - 107.64770976] + [60.22983081, 63.34160549, 66.43402974, 69.50442757, + 72.54998112, 75.56768635, 78.55430954, 81.50634428, + 84.41996557, 87.29098172, 90.11478465, 92.88629624, + 95.5999124, 98.24944658, 100.8280723, 103.32826814, + 105.74177, 108.05953355, 110.2717148, 112.36767828, + 114.33604136, 116.16476744, 117.84132062, 119.35289073, + 120.68669583, 121.83036181, 122.77236599, 123.50252116, + 124.01246252, 124.29608957, 124.34991302, 124.17326279, + 123.76832862, 123.14002805, 122.29571959, 121.24479895, + 119.99822708, 118.56803917, 116.96687825, 115.2075839, + 113.30285269, 111.26497738, 109.10566098, 106.83589588, + 104.46589822, 102.00508468, 99.46207977, 96.84474584, + 94.16022733, 91.4150024, 88.61493961, 85.76535568, + 82.87107126, 79.93646571, 76.96552925, 73.9619112, + 70.92896677, 67.86980202, 64.78731633, 61.68424611, + 58.56321009, 55.42675712, 52.27742195, 49.11779228, + 45.95059256, 42.77879745, 39.60579061, 36.43559544, + 33.27322843, 30.12525958, 27.00073787, 23.9127918, + 20.88152421, 17.93949588, 15.1425391, 12.59129497, + 10.47023949, 9.08837272, 8.80250072, 9.71016099, + 11.53222951, 13.91268049, 16.61145941, 19.49478556, + 22.48953441, 25.55411892, 28.66354138, 31.80189434, + 34.95848314, 38.12572775, 41.29797352, 44.47078752, + 47.64052254, 50.80403616, 53.95850349, 57.10128521, + 60.19189057, 63.30411711, 66.39705174, 69.46801919, + 72.51420324, 75.53260234, 78.51998605, 81.472852, + 84.38737999, 87.25938384, 90.08426172, 92.85694249, + 95.5718297, 98.22274508, 100.80287095, 103.30469507, + 105.71996264, 108.03963849, 110.2538871, 112.35208031, + 114.32284084, 116.15413465, 117.83342459, 119.34789444, + 120.68475046, 121.83160023, 122.77689598, 123.51041877, + 124.02376646, 124.31079722, 124.36797868, 124.19459854, + 123.79280774, 123.16749074, 122.32598014, 121.27765375, + 120.03346257, 118.60543913, 117.0062299, 115.2486827, + 113.34550579, 111.30900578, 109.15090068, 106.88219818, + 104.51312938, 102.0531253, 99.51082375, 96.8940993, + 94.21010748, 91.46533641, 88.6656635, 85.81641338, + 82.92241365, 79.98804982, 77.01731753, 74.01387089, + 70.98106934, 67.92202265, 64.8396335, 61.73664121, + 58.6156671, 55.47926225, 52.3299633, 49.17035952, + 46.00317641, 42.83138906, 39.65838052, 36.48817188, + 33.32577437, 30.17774741, 27.05311932, 23.9649783, + 20.93334666, 17.99061869, 15.19226962, 12.63816039, + 10.51116169, 9.11789559, 8.81495244, 9.70524228, + 11.51457879, 13.88698965, 16.5808083, 19.46100386, + 22.45371889, 25.51695489, 28.62547834, 31.76324145, + 34.91946275, 38.08650657, 41.25868161, 44.43153022, + 47.60138827, 50.76510173, 53.91983782, 57.06295199, + 60.15714089] ) ) np.testing.assert_array_almost_equal( to_test.azimuth.values, np.array( - [9.95354884, 13.51797616, 17.04447274, 20.52491225, - 23.95238719, 27.32131453, 30.62748235, 33.86804529, - 37.04147561, 40.14748222, 43.18691015, 46.1616283, - 49.07441566, 51.92885441, 54.72923301, 57.48046422, - 60.18802201, 62.85789652, 65.49656941, 68.11101133, - 70.70870005, 73.29766135, 75.88653487, 78.48466488, - 81.10221958, 83.75034334, 86.44134434, 89.18892358, - 92.00845226, 94.91730112, 97.93522718, 101.08482003, - 104.39199808, 107.88653257, 111.60254901, 115.57890087, - 119.8592317, 124.49141579, 129.52588659, 135.01216661, - 140.99278628, 147.49394635, 154.51314061, 162.00583842, - 169.87594269, 177.97643578, 186.12466927, 194.13016535, - 201.82539566, 209.08806319, 215.84877252, 222.08574604, - 227.81263221, 233.06524098, 237.89058053, 242.3392101, - 246.46057573, 250.30054442, 253.90036735, 257.29648187, - 260.52076061, 263.60096725, 266.56127763, 269.42279422, - 272.20401931, 274.92127106, 277.58904086, 280.22029634, - 282.82673312, 285.41898225, 288.00677995, 290.59910223, - 293.20426943, 295.83002468, 298.48358643, 301.17167776, - 303.9005345, 306.67589094, 309.50294481, 312.38630353, - 315.32991074, 318.33695647, 321.40977511, 324.54973374, - 327.75711805, 331.0310254, 334.36927221, 337.76832766, - 341.22328647, 344.72788848, 348.27459447, 351.85472406, - 355.45865247, 359.07606104, 2.69623101, 6.30836085, - 9.90413618, 13.46983714, 16.99770579, 20.47959425, - 23.90857277, 27.27903624, 30.58675167, 33.82885422, - 37.00379863, 40.1112785, 43.15212578, 46.12819847, - 49.04226661, 51.89790514, 54.69939668, 57.45164925, - 60.1601328, 62.83083386, 65.47023059, 68.08528982, - 70.68348498, 73.27283649, 75.86197737, 78.46024348, - 81.07779237, 83.72575502, 86.41642271, 89.16347527, - 91.9822575, 94.8901073, 97.90674111, 101.05469855, - 104.35983703, 107.85185424, 111.56478847, 115.53739258, - 119.81319937, 124.43997055, 129.46804248, 134.94688671, - 140.91907623, 147.41101064, 154.42060141, 161.90398226, - 169.76591074, 177.86024433, 186.00498616, 194.00988407, - 201.70716439, 208.97394261, 215.74010434, 221.98320309, - 227.71636618, 232.97505067, 237.80605645, 242.2598406, - 246.38581686, 250.22986166, 253.83325807, 257.23248545, - 260.45946126, 263.54199277, 266.50429617, 269.36751012, - 272.15016893, 274.86861866, 277.53737498, 280.1694266, - 282.77648747, 285.3692047, 287.95732871, 290.54984827, - 293.15509539, 295.78082413, 298.43426346, 301.12214686, - 303.85072078, 306.62573055, 309.45238569, 312.3353063, - 315.27844985, 318.28502138, 321.35737159, 324.49688503, - 327.70386592, 330.97743082, 334.31541557, 337.71430836, - 341.1692218, 344.67391164, 348.22085175, 351.80137124, - 355.40585053, 359.02397146, 2.64501087, 6.2581581, - 9.85597076] + [280.45072709, 282.11633942, 283.7817931, 285.45643927, + 287.14893238, 288.86746085, 290.61992811, 292.41409693, + 294.2577037, 296.15854802, 298.12456058, 300.16384793, + 302.28471289, 304.49564761, 306.80529163, 309.22234814, + 311.75544988, 314.41296211, 317.20271279, 320.13164164, + 323.20536119, 326.42763345, 329.79977946, 333.32005237, + 336.98302936, 340.77909977, 344.69413863, 348.7094615, + 352.80214024, 356.94571294, 1.11126554, 5.26879407, + 9.38869489, 13.44320688, 17.40764042, 21.26126792, + 24.98782151, 28.57560944, 32.01730949, 35.30953094, + 38.45224149, 41.44814098, 44.3020479, 47.02034252, + 49.61048749, 52.08063548, 54.43932309, 56.69524105, + 58.85707173, 60.93338359, 62.93257078, 64.86283035, + 66.73217087, 68.54844619, 70.31941227, 72.05280693, + 73.75645209, 75.438383, 77.10701127, 78.7713305, + 80.44118034, 82.12759231, 83.84324978, 85.60311369, + 87.42529319, 89.33228209, 91.3527583, 93.52426951, + 95.8973478, 98.54199782, 101.55824322, 105.09379965, + 109.37452354, 114.75768924, 121.82299993, 131.5046235, + 145.14346316, 163.81906282, 186.0060095, 206.77205045, + 222.71272671, 234.05797936, 242.20138813, 248.2792225, + 253.019093, 256.86789112, 260.10425777, 262.90749472, + 265.39707784, 267.65536873, 269.74098333, 271.69690061, + 273.5555329, 275.34199145, 277.07624914, 278.77460915, + 280.40857909, 282.07488513, 283.74081962, 285.41575589, + 287.10836546, 288.8268498, 290.57912232, 292.37295337, + 294.21608515, 296.11632179, 298.08159771, 300.12002294, + 302.239904, 304.44973753, 306.75816911, 309.1739102, + 311.70560495, 314.36163417, 317.1498465, 320.07720857, + 323.14936689, 326.37012512, 329.7408537, 333.25986209, + 336.9217888, 340.71708646, 344.63169085, 348.64697055, + 352.74003696, 356.88444935, 1.05129299, 5.21054047, + 9.33254359, 13.38947986, 17.35658738, 21.21306193, + 24.94256075, 28.53332377, 31.9779703, 35.2730625, + 38.41853254, 41.41705548, 44.27343437, 46.99404174, + 49.58633845, 52.05847988, 54.4190084, 56.6766226, + 58.84001384, 60.91776009, 62.91826507, 64.84973506, + 66.72018726, 68.5374834, 70.30938637, 72.04363982, + 73.74807028, 75.43071613, 77.09999039, 78.7648859, + 80.43523879, 82.12207338, 83.83806078, 85.59814255, + 87.420398, 89.32727554, 91.34738372, 93.51816308, + 95.88997704, 98.53255847, 101.54548029, 105.07568588, + 109.34766226, 114.71618082, 121.75632435, 131.39424631, + 144.96107287, 163.54312836, 185.67160477, 206.4608716, + 222.464149, 233.86642744, 242.05129037, 248.1579107, + 252.91794505, 256.78120715, 260.02821949, 262.83947716, + 265.33522383, 267.59832701, 269.68774445, 271.64669179, + 273.50774886, 275.29614781, 277.03195029, 278.73152558, + 280.35936035] ) ) diff --git a/pvlib/test/test_tools.py b/pvlib/test/test_tools.py index b88372fa32..e8b713df3a 100644 --- a/pvlib/test/test_tools.py +++ b/pvlib/test/test_tools.py @@ -1,10 +1,14 @@ +import datetime + import numpy as np import pandas as pd import pytest import pvlib.tools from pvlib import tools -from pvlib.test.test_solarposition import times + +times = pd.date_range(start=datetime.datetime(2014, 6, 24), + end=datetime.datetime(2014, 6, 26), freq='15Min') @pytest.mark.parametrize('keys, input_dict, expected', [ From 71dc8c58882b5d80cd6d1d5450994230fce406db Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Wed, 22 May 2019 18:43:25 -0700 Subject: [PATCH 12/13] Avoid duplicate import of pvlib tools --- pvlib/test/test_tools.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/pvlib/test/test_tools.py b/pvlib/test/test_tools.py index e8b713df3a..4a8fe73eb8 100644 --- a/pvlib/test/test_tools.py +++ b/pvlib/test/test_tools.py @@ -4,7 +4,6 @@ import pandas as pd import pytest -import pvlib.tools from pvlib import tools times = pd.date_range(start=datetime.datetime(2014, 6, 24), @@ -24,7 +23,7 @@ def test_build_kwargs(keys, input_dict, expected): def test_datetime_to_julian(): """ test transformation from datetime to julians """ - julians = pvlib.tools.datetime_to_julian(pd.to_datetime(times)) + julians = tools.datetime_to_julian(pd.to_datetime(times)) np.testing.assert_array_almost_equal(np.array(julians[:10]), np.array([ 2456832.5, From 03bb14384d896b45b485ea88d70521f878d1e23a Mon Sep 17 00:00:00 2001 From: Daniel Lassahn Date: Thu, 13 Jun 2019 14:21:17 +0200 Subject: [PATCH 13/13] Handle timezones in timeRange --- .stickler.yml | 2 +- pvlib/solarposition.py | 5 +- pvlib/test/test_solarposition.py | 133 +++++-------------------------- pvlib/tools.py | 19 +++-- 4 files changed, 33 insertions(+), 126 deletions(-) diff --git a/.stickler.yml b/.stickler.yml index 41aa301c06..dc360a6a7b 100644 --- a/.stickler.yml +++ b/.stickler.yml @@ -2,7 +2,7 @@ linters: flake8: python: 3 max-line-length: 79 - ignore: E201,E241,E226 + ignore: E201,E241,E226,W503,W504 files: ignore: - 'pvlib/_version.py' diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index 413fc7a212..ef2ca06eaa 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -28,7 +28,6 @@ from pvlib._deprecation import deprecated NS_PER_HR = 1.e9 * 3600. # nanoseconds per hour -JULIAN_2000 = 2451544.5 JULIAN_YEARS = 365.2425 @@ -1488,8 +1487,8 @@ def spencer(times, latitude, longitude): Springer Science & Business Media, 2008. """ - julians = datetime_to_julian(times) - julians_2000 = np.asarray(julians, dtype=np.float) - JULIAN_2000 + julians, julian_2000 = datetime_to_julian(times) + julians_2000 = np.asarray(julians, dtype=np.float) - julian_2000 latitude_radians = np.radians(latitude) day_time = (julians_2000 % 1) * 24 diff --git a/pvlib/test/test_solarposition.py b/pvlib/test/test_solarposition.py index d93df2046f..9052071573 100644 --- a/pvlib/test/test_solarposition.py +++ b/pvlib/test/test_solarposition.py @@ -754,122 +754,25 @@ def test_sun_rise_set_transit_geometric(expected_rise_set_spa, golden_mst): atol=np.abs(expected_transit_error).max()) -def test_spencer(): +def test_get_solar_position_spencer(golden_mst): """ test for the calculation based on spencer 1972 """ - latitude = 32.2 - longitude = -111 - to_test = solarposition.spencer(times, - latitude, - longitude) - np.testing.assert_array_almost_equal( - to_test.zenith.values, - np.array( - [60.22983081, 63.34160549, 66.43402974, 69.50442757, - 72.54998112, 75.56768635, 78.55430954, 81.50634428, - 84.41996557, 87.29098172, 90.11478465, 92.88629624, - 95.5999124, 98.24944658, 100.8280723, 103.32826814, - 105.74177, 108.05953355, 110.2717148, 112.36767828, - 114.33604136, 116.16476744, 117.84132062, 119.35289073, - 120.68669583, 121.83036181, 122.77236599, 123.50252116, - 124.01246252, 124.29608957, 124.34991302, 124.17326279, - 123.76832862, 123.14002805, 122.29571959, 121.24479895, - 119.99822708, 118.56803917, 116.96687825, 115.2075839, - 113.30285269, 111.26497738, 109.10566098, 106.83589588, - 104.46589822, 102.00508468, 99.46207977, 96.84474584, - 94.16022733, 91.4150024, 88.61493961, 85.76535568, - 82.87107126, 79.93646571, 76.96552925, 73.9619112, - 70.92896677, 67.86980202, 64.78731633, 61.68424611, - 58.56321009, 55.42675712, 52.27742195, 49.11779228, - 45.95059256, 42.77879745, 39.60579061, 36.43559544, - 33.27322843, 30.12525958, 27.00073787, 23.9127918, - 20.88152421, 17.93949588, 15.1425391, 12.59129497, - 10.47023949, 9.08837272, 8.80250072, 9.71016099, - 11.53222951, 13.91268049, 16.61145941, 19.49478556, - 22.48953441, 25.55411892, 28.66354138, 31.80189434, - 34.95848314, 38.12572775, 41.29797352, 44.47078752, - 47.64052254, 50.80403616, 53.95850349, 57.10128521, - 60.19189057, 63.30411711, 66.39705174, 69.46801919, - 72.51420324, 75.53260234, 78.51998605, 81.472852, - 84.38737999, 87.25938384, 90.08426172, 92.85694249, - 95.5718297, 98.22274508, 100.80287095, 103.30469507, - 105.71996264, 108.03963849, 110.2538871, 112.35208031, - 114.32284084, 116.15413465, 117.83342459, 119.34789444, - 120.68475046, 121.83160023, 122.77689598, 123.51041877, - 124.02376646, 124.31079722, 124.36797868, 124.19459854, - 123.79280774, 123.16749074, 122.32598014, 121.27765375, - 120.03346257, 118.60543913, 117.0062299, 115.2486827, - 113.34550579, 111.30900578, 109.15090068, 106.88219818, - 104.51312938, 102.0531253, 99.51082375, 96.8940993, - 94.21010748, 91.46533641, 88.6656635, 85.81641338, - 82.92241365, 79.98804982, 77.01731753, 74.01387089, - 70.98106934, 67.92202265, 64.8396335, 61.73664121, - 58.6156671, 55.47926225, 52.3299633, 49.17035952, - 46.00317641, 42.83138906, 39.65838052, 36.48817188, - 33.32577437, 30.17774741, 27.05311932, 23.9649783, - 20.93334666, 17.99061869, 15.19226962, 12.63816039, - 10.51116169, 9.11789559, 8.81495244, 9.70524228, - 11.51457879, 13.88698965, 16.5808083, 19.46100386, - 22.45371889, 25.51695489, 28.62547834, 31.76324145, - 34.91946275, 38.08650657, 41.25868161, 44.43153022, - 47.60138827, 50.76510173, 53.91983782, 57.06295199, - 60.15714089] - ) - ) - - np.testing.assert_array_almost_equal( - to_test.azimuth.values, - np.array( - [280.45072709, 282.11633942, 283.7817931, 285.45643927, - 287.14893238, 288.86746085, 290.61992811, 292.41409693, - 294.2577037, 296.15854802, 298.12456058, 300.16384793, - 302.28471289, 304.49564761, 306.80529163, 309.22234814, - 311.75544988, 314.41296211, 317.20271279, 320.13164164, - 323.20536119, 326.42763345, 329.79977946, 333.32005237, - 336.98302936, 340.77909977, 344.69413863, 348.7094615, - 352.80214024, 356.94571294, 1.11126554, 5.26879407, - 9.38869489, 13.44320688, 17.40764042, 21.26126792, - 24.98782151, 28.57560944, 32.01730949, 35.30953094, - 38.45224149, 41.44814098, 44.3020479, 47.02034252, - 49.61048749, 52.08063548, 54.43932309, 56.69524105, - 58.85707173, 60.93338359, 62.93257078, 64.86283035, - 66.73217087, 68.54844619, 70.31941227, 72.05280693, - 73.75645209, 75.438383, 77.10701127, 78.7713305, - 80.44118034, 82.12759231, 83.84324978, 85.60311369, - 87.42529319, 89.33228209, 91.3527583, 93.52426951, - 95.8973478, 98.54199782, 101.55824322, 105.09379965, - 109.37452354, 114.75768924, 121.82299993, 131.5046235, - 145.14346316, 163.81906282, 186.0060095, 206.77205045, - 222.71272671, 234.05797936, 242.20138813, 248.2792225, - 253.019093, 256.86789112, 260.10425777, 262.90749472, - 265.39707784, 267.65536873, 269.74098333, 271.69690061, - 273.5555329, 275.34199145, 277.07624914, 278.77460915, - 280.40857909, 282.07488513, 283.74081962, 285.41575589, - 287.10836546, 288.8268498, 290.57912232, 292.37295337, - 294.21608515, 296.11632179, 298.08159771, 300.12002294, - 302.239904, 304.44973753, 306.75816911, 309.1739102, - 311.70560495, 314.36163417, 317.1498465, 320.07720857, - 323.14936689, 326.37012512, 329.7408537, 333.25986209, - 336.9217888, 340.71708646, 344.63169085, 348.64697055, - 352.74003696, 356.88444935, 1.05129299, 5.21054047, - 9.33254359, 13.38947986, 17.35658738, 21.21306193, - 24.94256075, 28.53332377, 31.9779703, 35.2730625, - 38.41853254, 41.41705548, 44.27343437, 46.99404174, - 49.58633845, 52.05847988, 54.4190084, 56.6766226, - 58.84001384, 60.91776009, 62.91826507, 64.84973506, - 66.72018726, 68.5374834, 70.30938637, 72.04363982, - 73.74807028, 75.43071613, 77.09999039, 78.7648859, - 80.43523879, 82.12207338, 83.83806078, 85.59814255, - 87.420398, 89.32727554, 91.34738372, 93.51816308, - 95.88997704, 98.53255847, 101.54548029, 105.07568588, - 109.34766226, 114.71618082, 121.75632435, 131.39424631, - 144.96107287, 163.54312836, 185.67160477, 206.4608716, - 222.464149, 233.86642744, 242.05129037, 248.1579107, - 252.91794505, 256.78120715, 260.02821949, 262.83947716, - 265.33522383, 267.59832701, 269.68774445, 271.64669179, - 273.50774886, 275.29614781, 277.03195029, 278.73152558, - 280.35936035] - ) - ) + times = pd.date_range(datetime.datetime(2003, 10, 17, 12, 30, 30), + periods=1, freq='D') + spencer_data = solarposition.spencer(times, + golden_mst.latitude, + golden_mst.longitude) + + expected_solpos = golden_mst.get_solarposition(times) + + np.testing.assert_array_almost_equal(spencer_data.zenith.values, + expected_solpos.zenith.values, + decimal=0) + np.testing.assert_array_almost_equal(spencer_data.azimuth.values, + expected_solpos.azimuth.values, + decimal=0) + np.testing.assert_array_almost_equal(spencer_data.elevation.values, + expected_solpos.elevation.values, + decimal=0) # put numba tests at end of file to minimize reloading diff --git a/pvlib/tools.py b/pvlib/tools.py index 998055c36d..065965a044 100644 --- a/pvlib/tools.py +++ b/pvlib/tools.py @@ -9,9 +9,9 @@ import pandas as pd import pytz -JULIAN_2000 = 2451544.5 -DT_2000 = dt.datetime(2000, 1, 1) +JULIAN_2000 = 2451545 DAY_SECONDS = 60 * 60 * 24 +DT_2000 = pd.to_datetime(dt.datetime(2000, 1, 1)).tz_localize('UTC') def cosd(angle): @@ -437,17 +437,22 @@ def datetime_to_julian(times): Parameters ---------- - times : :class:`pandas.DatetimeIndex` + times : pandas.DatetimeIndex Corresponding timestamps, must be localized to the timezone for the ``latitude`` and ``longitude`` Returns ------- Float64Index The float index contains julian dates + Float + julian 2000 in UTC referenced to local time """ - - delta = times - DT_2000 + dt_2000 = DT_2000.tz_convert(times.tzinfo) + hours_difference = (DT_2000.tz_localize(None) - + dt_2000.tz_localize(None)).seconds / 3600 + julian_2000 = JULIAN_2000 - (hours_difference/24) + delta = times - dt_2000 delta_julians = (delta.seconds + delta.microseconds / 1e6) return ( - JULIAN_2000 + delta.days + delta_julians / DAY_SECONDS - ) + julian_2000 + delta.days + delta_julians / DAY_SECONDS + ), julian_2000