Add analytical azimuth method to solarposition.py. (#349)

* Add analytical azimuth method.

* Number linked references.

* Add test for solar_azimuth_analytical.

* Add an entry in the api.rst file.

* Correct for different azimuth convention.

* Update test tolerance.

Not as accurate when azimuth angle is close to due north.

* Added an entry in the v0.4.6 whatsnew file.

* Update assert statements.

* Update threshold for `atol=0.55` assert statement.

* Update test, only when the sun is above horizon.

Loops no longer needed.

Author: veronicaguo

docs/sphinx/source/api.rst

@@ -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

docs/sphinx/source/whatsnew/v0.4.6.rst

@@ -17,6 +17,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
 ~~~~~~~~~~~
@@ -31,6 +32,7 @@ Testing
 
 * Added explicit tests for aoi and aoi_projection functions.
 * Update test of `ModelChain.__repr__` to take in account :issue:`352`
+* Added a test for solar_azimuth_analytical function.
 
 
 Contributors
@@ -41,3 +43,4 @@ Contributors
 * Alaina Kafkes
 * Birgit Schachler
 * Jonathan Gaffiot
+* Siyan (Veronica) Guo

pvlib/solarposition.py

@@ -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(
+        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.

pvlib/test/test_solarposition.py

@@ -420,3 +420,33 @@ 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
+    solar_zenith = np.deg2rad(output['zenith'])
+    # 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)
+
+    idx = np.where(solar_zenith < np.pi/2)
+    assert np.allclose(azimuth_1[idx], solar_azimuth.as_matrix()[idx],
+                       atol=0.01)
+    assert np.allclose(azimuth_2[idx], solar_azimuth.as_matrix()[idx],
+                       atol=0.017)