Projected zenith convenience function

docs/sphinx/source/reference/effects_on_pv_system_output/shading.rst

@@ -9,4 +9,5 @@ Shading
    shading.ground_angle
    shading.masking_angle
    shading.masking_angle_passias
-   shading.sky_diffuse_passias
+   shading.sky_diffuse_passias
+   shading.projected_solar_zenith_angle

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

@@ -12,6 +12,8 @@ Enhancements
 * :py:func:`pvlib.bifacial.infinite_sheds.get_irradiance` and
   :py:func:`pvlib.bifacial.infinite_sheds.get_irradiance_poa` now include
   shaded fraction in returned variables. (:pull:`1871`)
+* Added function :py:func:`pvlib.shading.projected_solar_zenith_angle`,
+  a common calculation in shading and tracking. (:issue:`1734`, :pull:`1904`)
 
 Bug fixes
 ~~~~~~~~~
@@ -33,3 +35,4 @@ Contributors
 * Will Hobbs (:ghuser:`williamhobbs`)
 * Anton Driesse (:ghuser:`adriesse`)
 * :ghuser:`matsuobasho`
+* Echedey Luis (:ghuser:`echedey-ls`)

pvlib/shading.py

@@ -232,3 +232,64 @@ def sky_diffuse_passias(masking_angle):
        Available at https://www.nrel.gov/docs/fy18osti/67399.pdf
     """
     return 1 - cosd(masking_angle/2)**2
+
+
+def projected_solar_zenith_angle(surface_tilt, surface_azimuth,
+                                 solar_apparent_elevation, solar_azimuth):
+    r"""
+    Calculate projected solar zenith angle in degrees.
+
+    This is common in track and shadow computation [1]_ [2]_ [3]_.
+
+    Parameters
+    ----------
+    surface_tilt : numeric
+        Axis tilt angle in degrees. From horizontal plane to array plane.
+    surface_azimuth : numeric
+        Axis azimuth angle in degrees.
+        North = 0°; East = 90°; South = 180°; West = 270°
+    solar_apparent_elevation : numeric
+        Sun's apparent elevation in degrees.
+    solar_azimuth : numeric
+        Sun's azimuth in degrees.
+
+    Returns
+    -------
+    Projected_solar_zenith : numeric
+        In degrees.
+
+    References
+    ----------
+    .. [1] K. Anderson and M. Mikofski, ‘Slope-Aware Backtracking for
+       Single-Axis Trackers’, National Renewable Energy Lab. (NREL), Golden,
+       CO (United States); Det Norske Veritas Group, Oslo (Norway),
+       NREL/TP-5K00-76626, Jul. 2020. :doi:`10.2172/1660126`.
+    .. [2] W. F. Marion and A. P. Dobos, ‘Rotation Angle for the Optimum
+       Tracking of One-Axis Trackers’, National Renewable Energy Lab. (NREL),
+       Golden, CO (United States), NREL/TP-6A20-58891, Jul. 2013.
+       :doi:`10.2172/1089596`.
+    .. [3] E. Lorenzo, L. Narvarte, and J. Muñoz, ‘Tracking and back-tracking’,
+       Progress in Photovoltaics: Research and Applications, vol. 19, no. 6,
+       pp. 747–753, 2011, :doi:`10.1002/pip.1085`.
+    """
+    # Avoid recalculating these values
+    cosd_solar_apparent_elevation = cosd(solar_apparent_elevation)
+    cosd_surface_azimuth = cosd(surface_azimuth)
+    sind_surface_azimuth = sind(surface_azimuth)
+    sind_surface_tilt = sind(surface_tilt)
+
+    # Notation from [1]
+    # Sun's x, y, z coords
+    sx = cosd_solar_apparent_elevation * sind(solar_azimuth)
+    sy = cosd_solar_apparent_elevation * cosd(solar_azimuth)
+    sz = sind(solar_apparent_elevation)
+    # Eq. (4); sx', sz' values from sun coordinates projected onto surface
+    sx_prime = sx * cosd_surface_azimuth - sy * sind_surface_azimuth
+    sz_prime = (
+        sx * sind_surface_azimuth * sind_surface_tilt
+        + sy * sind_surface_tilt * cosd_surface_azimuth
+        + sz * cosd(surface_tilt)
+    )
+    # Eq. (5); angle between sun's beam and surface
+    theta_T = np.degrees(np.arctan2(sx_prime, sz_prime))
+    return theta_T

pvlib/tests/test_shading.py

@@ -2,8 +2,11 @@
 import pandas as pd
 
 from pandas.testing import assert_series_equal
+from numpy.testing import assert_allclose
 import pytest
+from datetime import timezone, timedelta
 
+import pvlib
 from pvlib import shading
 
 
@@ -37,7 +40,7 @@ def test__ground_angle_zero_gcr():
 
 @pytest.fixture
 def surface_tilt():
-    idx = pd.date_range('2019-01-01', freq='h', periods=3)
+    idx = pd.date_range("2019-01-01", freq="h", periods=3)
     return pd.Series([0, 20, 90], index=idx)
 
 
@@ -104,3 +107,87 @@ def test_sky_diffuse_passias_scalar(average_masking_angle, shading_loss):
     for angle, loss in zip(average_masking_angle, shading_loss):
         actual_loss = shading.sky_diffuse_passias(angle)
         assert np.isclose(loss, actual_loss)
+
+
+@pytest.fixture
+def true_tracking_angle_and_inputs():
+    # data retrieved from NREL Slope-Aware Backtracking for Single-Axis
+    # Trackers
+    # doi.org/10.2172/1660126 ; Accessed on 2023-11-06.
+    tzinfo = timezone(timedelta(hours=-5))
+    axis_tilt_angle = 9.666  # deg
+    axis_azimuth_angle = 195.0  # deg
+    timedata = pd.DataFrame(
+        columns=("Apparent Elevation", "Solar Azimuth", "True-Tracking"),
+        data=(
+            (2.404287, 122.791770, -84.440),
+            (11.263058, 133.288729, -72.604),
+            (18.733558, 145.285552, -59.861),
+            (24.109076, 158.939435, -45.578),
+            (26.810735, 173.931802, -28.764),
+            (26.482495, 189.371536, -8.475),
+            (23.170447, 204.136810, 15.120),
+            (17.296785, 217.446538, 39.562),
+            (9.461862, 229.102218, 61.587),
+            (0.524817, 239.330401, 79.530),
+        ),
+    )
+    timedata.index = pd.date_range(
+        "2019-01-01T08", "2019-01-01T17", freq="1H", tz=tzinfo
+    )
+    timedata["Apparent Zenith"] = 90.0 - timedata["Apparent Elevation"]
+    return (axis_tilt_angle, axis_azimuth_angle, timedata)
+
+
+def test_projected_solar_zenith_angle_numeric(true_tracking_angle_and_inputs):
+    psz_func = shading.projected_solar_zenith_angle
+    axis_tilt, axis_azimuth, timedata = true_tracking_angle_and_inputs
+    psz = psz_func(
+        axis_tilt,
+        axis_azimuth,
+        timedata["Apparent Elevation"],
+        timedata["Solar Azimuth"],
+    )
+    assert_allclose(psz, timedata["True-Tracking"], atol=1e-3)
+    # test equivalence against pvlib.tracking.singleaxis
+    singleaxis = pvlib.tracking.singleaxis(timedata["Apparent Zenith"],
+                                           timedata["Solar Azimuth"],
+                                           axis_tilt, axis_azimuth,
+                                           backtrack=False)
+    assert_allclose(psz, singleaxis["tracker_theta"])
+    # test by changing axis azimuth and tilt
+    psz = psz_func(
+        -axis_tilt,
+        axis_azimuth-180,
+        timedata["Apparent Elevation"],
+        timedata["Solar Azimuth"],
+    )
+    assert_allclose(psz, -timedata["True-Tracking"], atol=1e-3)
+
+
+@pytest.mark.parametrize(
+    "cast_type, cast_func",
+    [
+        (float, float),
+        (np.ndarray, lambda x: np.array([x])),
+        (pd.Series, lambda x: pd.Series(data=[x])),
+    ],
+)
+def test_projected_solar_zenith_angle_datatypes(
+    cast_type, cast_func, true_tracking_angle_and_inputs
+):
+    psz_func = shading.projected_solar_zenith_angle
+    axis_tilt, axis_azimuth, timedata = true_tracking_angle_and_inputs
+    sun_apparent_zenith = timedata["Apparent Zenith"].iloc[0]
+    sun_azimuth = timedata["Solar Azimuth"].iloc[0]
+
+    axis_tilt, axis_azimuth, sun_apparent_zenith, sun_azimuth = (
+        cast_func(axis_tilt),
+        cast_func(axis_azimuth),
+        cast_func(sun_apparent_zenith),
+        cast_func(sun_azimuth),
+    )
+    psz = psz_func(
+        axis_tilt, axis_azimuth, sun_apparent_zenith, axis_azimuth
+    )
+    assert isinstance(psz, cast_type)