Skip to content

Add analytical azimuth method to solarposition.py. #349

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 11 commits into from
Aug 1, 2017
Merged

Add analytical azimuth method to solarposition.py. #349

Show file tree
Hide file tree
Changes from 7 commits
Commits
Show all changes
11 commits
Select commit Hold shift + click to select a range
52cfb3b
Add analytical azimuth method.
veronicaguo Jul 14, 2017
5395164
Number linked references.
veronicaguo Jul 24, 2017
4e56576
Add test for solar_azimuth_analytical.
veronicaguo Jul 24, 2017
c0698d4
Add an entry in the api.rst file.
veronicaguo Jul 24, 2017
e3c7530
Correct for different azimuth convention.
veronicaguo Jul 24, 2017
9c2b89b
Update test tolerance.
veronicaguo Jul 24, 2017
b6d6fb9
Added an entry in the v0.4.6 whatsnew file.
veronicaguo Jul 24, 2017
d1b38b5
Update assert statements.
veronicaguo Jul 25, 2017
2c346c4
Update threshold for `atol=0.55` assert statement.
veronicaguo Jul 26, 2017
3728f15
Merge branch 'master' of github.com:pvlib/pvlib-python
veronicaguo Jul 27, 2017
31bf30d
Update test, only when the sun is above horizon.
veronicaguo Jul 27, 2017
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions docs/sphinx/source/api.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -73,6 +73,7 @@ calculations.
:toctree: generated/

solarposition.solar_zenith_analytical
solarposition.solar_azimuth_analytical
solarposition.declination_spencer71
solarposition.declination_cooper69
solarposition.equation_of_time_spencer71
Expand Down
3 changes: 3 additions & 0 deletions docs/sphinx/source/whatsnew/v0.4.6.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@ Bug fixes
Enhancements
~~~~~~~~~~~~
* Added default values to docstrings of all functions (:issue:`336`)
* Added analytical method that calculates solar azimuth angle (:issue:`291`)

API Changes
~~~~~~~~~~~
Expand All @@ -28,6 +29,7 @@ Testing
~~~~~~~

* Added explicit tests for aoi and aoi_projection functions.
* Added a test for solar_azimuth_analytical function.


Contributors
Expand All @@ -37,3 +39,4 @@ Contributors
* Uwe Krien
* Alaina Kafkes
* Birgit Schachler
* Siyan (Veronica) Guo
47 changes: 47 additions & 0 deletions pvlib/solarposition.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1052,6 +1052,53 @@ def declination_cooper69(dayofyear):
return np.deg2rad(23.45 * np.sin(day_angle + (2.0 * np.pi / 365.0) * 285.0))


def solar_azimuth_analytical(latitude, hour_angle, declination, zenith):
"""
Analytical expression of solar azimuth angle based on spherical
trigonometry.

Parameters
----------
latitude : numeric
Latitude of location in radians.
hour_angle : numeric
Hour angle in the local solar time in radians.
declination : numeric
Declination of the sun in radians.
zenith : numeric
Solar zenith angle in radians.

Returns
-------
azimuth : numeric
Solar azimuth angle in radians.

References
----------
[1] J. A. Duffie and W. A. Beckman, "Solar Engineering of Thermal
Processes, 3rd Edition" pp. 14, J. Wiley and Sons, New York (2006)

[2] J. H. Seinfeld and S. N. Pandis, "Atmospheric Chemistry and Physics"
p. 132, J. Wiley (1998)

[3] `Wikipedia: Solar Azimuth Angle
<https://en.wikipedia.org/wiki/Solar_azimuth_angle>`_

[4] `PVCDROM: Azimuth Angle <http://www.pveducation.org/pvcdrom/2-properties
-sunlight/azimuth-angle>`_

See Also
--------
declination_spencer71
declination_cooper69
hour_angle
solar_zenith_analytical
"""
return np.sign(hour_angle) * np.abs(np.arccos((np.cos(zenith) * np.sin(
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@veronicaguo and @wholmgren sorry for super late review, and just a tiny nitpick, but I don't think np.abs is necessary here because according to numpy documentation on np.arccos it is only defined from on [0, pi]

The angle of the ray intersecting the unit circle at the given x-coordinate in radians [0, pi].

Therefore, it will never return a negative number, so no need for absolute value.

latitude) - np.sin(declination)) / (np.sin(zenith) * np.cos(
latitude)))) + np.pi


def solar_zenith_analytical(latitude, hour_angle, declination):
"""
Analytical expression of solar zenith angle based on spherical trigonometry.
Expand Down
37 changes: 37 additions & 0 deletions pvlib/test/test_solarposition.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -420,3 +420,40 @@ def test_analytical_zenith():
zenith_2 = solarposition.solar_zenith_analytical(lat_rad, hour_angle, decl)
assert np.allclose(zenith_1, solar_zenith, atol=0.015)
assert np.allclose(zenith_2, solar_zenith, atol=0.025)


def test_analytical_azimuth():
times = pd.DatetimeIndex(start="1/1/2015 0:00", end="12/31/2015 23:00",
freq="H").tz_localize('Etc/GMT+8')
lat, lon = 37.8, -122.25
lat_rad = np.deg2rad(lat)
output = solarposition.spa_python(times, lat, lon, 100)
solar_azimuth = np.deg2rad(output['azimuth']) # spa
# spencer
eot = solarposition.equation_of_time_spencer71(times.dayofyear)
hour_angle = np.deg2rad(solarposition.hour_angle(times, lon, eot))
decl = solarposition.declination_spencer71(times.dayofyear)
zenith = solarposition.solar_zenith_analytical(lat_rad, hour_angle, decl)
azimuth_1 = solarposition.solar_azimuth_analytical(lat_rad, hour_angle,
decl, zenith)
# pvcdrom and cooper
eot = solarposition.equation_of_time_pvcdrom(times.dayofyear)
hour_angle = np.deg2rad(solarposition.hour_angle(times, lon, eot))
decl = solarposition.declination_cooper69(times.dayofyear)
zenith = solarposition.solar_zenith_analytical(lat_rad, hour_angle, decl)
azimuth_2 = solarposition.solar_azimuth_analytical(lat_rad, hour_angle,
decl, zenith)
for idx, a in enumerate(azimuth_1):
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can you add some comments as to why you've structured the for loop and the assert statements this way? Please also replace a with something more descriptive such as azi or azimuth.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sure. The analytical azimuth formula is quite accurate in every direction, with an error or 1-2 degrees, except when the sun is in the north. So I separated those cases out through the for loops and if statements and gave them different atol values. Also, since the sun is only in the north during midnight when no power is produced, it shouldn't impact solar performance calculations. The test code is not in the most concise form right now. I will try slim it down a bit. Thank you for the feedback.

if a < 0.7:
assert np.allclose(a, solar_azimuth[idx], atol=0.025) or \
np.allclose(a + np.pi * 2, solar_azimuth[idx], atol=0.55)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

atol=0.55 implies an accuracy of 32 degrees, correct? That seems quite bad. Not sure if the test needs to change or if that's a limitation of the algorithm under some input parameters.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

atol=0.55 only when the sun is in the north, which is every night from 12-1 AM. It is otherwise quite accurate during the rest of the day. I assume it's probably a limitation with the analytical method that most references share.

Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In case a system in the Southern hemisphere is analyzed, the sun will be in the North during the day. Shouldn't the code also be able to correctly handle this case?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks @trichrt for pointing it out. The analytical formula actually takes care of it inherently. When in the Southern hemisphere, the calculation is less accurate when the sun is in the south during midnight. The current test is for a Northern hemisphere location. I can put together a test for Southern hemisphere if desired.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In that case we should restrict the atol=0.55 tests to azimuth angles that are close to 0/360 (as defined by the more accurate models).

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Agreed. The discrepancies start to appear within a 30-degree range around 0/360, i.e. from 345 to 15 degree. Have updated the azimuth threshold in the latest commit.

else:
assert np.allclose(a, solar_azimuth[idx], atol=0.025)

for idx, a in enumerate(azimuth_2):
if a < 0.7:
assert np.allclose(a, solar_azimuth[idx], atol=0.025) or \
np.allclose(a + np.pi * 2, solar_azimuth[idx], atol=0.55)
else:
assert np.allclose(a, solar_azimuth[idx], atol=0.025)