Repository: pvlib/pvlib-python (python)

Issue:
incorrect horizon test in tracker.singleaxis
**Describe the bug**
`    wid[zp <= 0] = np.nan` sets the tracker rotation angle to `np.nan` when the sun is behind the array plane, but the array plane is not yet defined. At this line, the array plane is defined only by the tracker tilt and azimuth, but not its rotation.

**To Reproduce**
```
import pvlib
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

loc = pvlib.location.Location(latitude=30, longitude=110)
sat = pvlib.tracking.SingleAxisTracker(axis_tilt=90, axis_azimuth=180, backtrack=False, max_angle=180)
sat.localize(location=loc)

solar_azimuth = pd.Series(data=np.arange(0.0, 359.0))
solar_zenith = 45.0 * solar_azimuth / solar_azimuth

df = sat.singleaxis(apparent_zenith=solar_zenith, apparent_azimuth=solar_azimuth)

p1 = plt.plot(df)
plt.legend(df.columns)
```
The plot shows that the tracker rotation `tracker_theta` is limited to between -90 and 90.

**Expected behavior**
In the above example, the sun rotates in azimuth at 45 zenith. The vertical tracker without rotation limit should follow the sun around the complete circle.

**Screenshots**
![test_sat](https://user-images.githubusercontent.com/5393711/45325941-88428d80-b50f-11e8-918a-2cfc3e8aa460.png)

**Versions:**
 - ``pvlib.__version__``: 0.5.2 

**Additional context**

One fix would be to substitute a test for solar zenith, so that the tracker rotation is not set when the sun is below the visible horizon. I also wonder if we should set the rotation to 0 (stow position)  instead of `np.nan` in this case.

`    wid[z <= 0] = np.nan`
`    wid[z <= 0] = 0.0`


Pull request description:
fix incorrect nan assignments in singleaxis, allow scalar and 1d array input
 - [x] Closes #569
 - [x] I am familiar with the [contributing guidelines](http://pvlib-python.readthedocs.io/en/latest/contributing.html).
 - [x] Fully tested. Added and/or modified tests to ensure correct behavior for all reasonable inputs. Tests (usually) must pass on the TravisCI and Appveyor testing services.
 - [x] Updates entries to `docs/sphinx/source/api.rst` for API changes.
 - [x] Adds description and name entries in the appropriate `docs/sphinx/source/whatsnew` file for all changes.
 - [x] Code quality and style is sufficient. Passes LGTM and SticklerCI checks.
 - [x] New code is fully documented. Includes sphinx/numpydoc compliant docstrings and comments in the code where necessary.
 - [x] Pull request is nearly complete and ready for detailed review.

I punted on the question of nan/0. 

~Would also be nice to remove the requirement that input data be series. Yuck.~ Now accepts scalar/1d input

Patch:
diff --git a/docs/sphinx/source/whatsnew/v0.6.0.rst b/docs/sphinx/source/whatsnew/v0.6.0.rst
index 2a5bf2d..83af55a 100644
--- a/docs/sphinx/source/whatsnew/v0.6.0.rst
+++ b/docs/sphinx/source/whatsnew/v0.6.0.rst
@@ -125,6 +125,7 @@ Enhancements
   `doc` (requirements for minimal documentation build), `test` (requirements
   for testing), and `all` (optional + doc + test). (:issue:`553`, :issue:`483`)
 * Set default alpha to 1.14 in :func:`~pvlib.atmosphere.angstrom_aod_at_lambda` (:issue:`563`)
+* tracking.singleaxis now accepts scalar and 1D-array input.
 
 
 Bug fixes
@@ -152,6 +153,9 @@ Bug fixes
 * Fix bug in get_relative_airmass(model='youngirvine1967'). (:issue:`545`)
 * Fix bug in variable names returned by forecast.py's HRRR_ESRL model.
   (:issue:`557`)
+* Fixed bug in tracking.singleaxis that mistakenly assigned nan values when
+  the Sun was still above the horizon. No effect on systems with axis_tilt=0.
+  (:issue:`569`)
 
 
 Documentation
diff --git a/pvlib/tracking.py b/pvlib/tracking.py
index c1d504b..a2d9700 100644
--- a/pvlib/tracking.py
+++ b/pvlib/tracking.py
@@ -260,10 +260,10 @@ def singleaxis(apparent_zenith, apparent_azimuth,
 
     Parameters
     ----------
-    apparent_zenith : Series
+    apparent_zenith : float, 1d array, or Series
         Solar apparent zenith angles in decimal degrees.
 
-    apparent_azimuth : Series
+    apparent_azimuth : float, 1d array, or Series
         Solar apparent azimuth angles in decimal degrees.
 
     axis_tilt : float, default 0
@@ -296,7 +296,7 @@ def singleaxis(apparent_zenith, apparent_azimuth,
 
     Returns
     -------
-    DataFrame with the following columns:
+    dict or DataFrame with the following columns:
 
     * tracker_theta: The rotation angle of the tracker.
         tracker_theta = 0 is horizontal, and positive rotation angles are
@@ -318,6 +318,18 @@ def singleaxis(apparent_zenith, apparent_azimuth,
     # MATLAB to Python conversion by
     # Will Holmgren (@wholmgren), U. Arizona. March, 2015.
 
+    if isinstance(apparent_zenith, pd.Series):
+        index = apparent_zenith.index
+    else:
+        index = None
+
+    # convert scalars to arrays
+    apparent_azimuth = np.atleast_1d(apparent_azimuth)
+    apparent_zenith = np.atleast_1d(apparent_zenith)
+
+    if apparent_azimuth.ndim > 1 or apparent_zenith.ndim > 1:
+        raise ValueError('Input dimensions must not exceed 1')
+
     # Calculate sun position x, y, z using coordinate system as in [1], Eq 2.
 
     # Positive y axis is oriented parallel to earth surface along tracking axis
@@ -334,15 +346,6 @@ def singleaxis(apparent_zenith, apparent_azimuth,
     # Rotate sun azimuth to coordinate system as in [1]
     # to calculate sun position.
 
-    try:
-        pd.util.testing.assert_index_equal(apparent_azimuth.index,
-                                           apparent_zenith.index)
-    except AssertionError:
-        raise ValueError('apparent_azimuth.index and '
-                         'apparent_zenith.index must match.')
-
-    times = apparent_azimuth.index
-
     az = apparent_azimuth - 180
     apparent_elevation = 90 - apparent_zenith
     x = cosd(apparent_elevation) * sind(az)
@@ -408,10 +411,11 @@ def singleaxis(apparent_zenith, apparent_azimuth,
 
     # Calculate angle from x-y plane to projection of sun vector onto x-z plane
     # and then obtain wid by translating tmp to convention for rotation angles.
-    wid = pd.Series(90 - np.degrees(np.arctan2(zp, xp)), index=times)
+    wid = 90 - np.degrees(np.arctan2(zp, xp))
 
     # filter for sun above panel horizon
-    wid[zp <= 0] = np.nan
+    zen_gt_90 = apparent_zenith > 90
+    wid[zen_gt_90] = np.nan
 
     # Account for backtracking; modified from [1] to account for rotation
     # angle convention being used here.
@@ -423,14 +427,11 @@ def singleaxis(apparent_zenith, apparent_azimuth,
         # (always positive b/c acosd returns values between 0 and 180)
         wc = np.degrees(np.arccos(temp))
 
-        v = wid < 0
-        widc = pd.Series(index=times)
-        widc[~v] = wid[~v] - wc[~v]  # Eq 4 applied when wid in QI
-        widc[v] = wid[v] + wc[v]     # Eq 4 applied when wid in QIV
+        # Eq 4 applied when wid in QIV (wid < 0 evalulates True), QI
+        tracker_theta = np.where(wid < 0, wid + wc, wid - wc)
     else:
-        widc = wid
+        tracker_theta = wid
 
-    tracker_theta = widc.copy()
     tracker_theta[tracker_theta > max_angle] = max_angle
     tracker_theta[tracker_theta < -max_angle] = -max_angle
 
@@ -447,7 +448,6 @@ def singleaxis(apparent_zenith, apparent_azimuth,
 
     # calculate angle-of-incidence on panel
     aoi = np.degrees(np.arccos(np.abs(np.sum(sun_vec*panel_norm, axis=0))))
-    aoi = pd.Series(aoi, index=times)
 
     # calculate panel tilt and azimuth
     # in a coordinate system where the panel tilt is the
@@ -491,9 +491,8 @@ def singleaxis(apparent_zenith, apparent_azimuth,
 #     surface_azimuth = pd.Series(
 #         np.degrees(np.arctan(projected_normal[:,1]/projected_normal[:,0])),
 #                                 index=times)
-    surface_azimuth = pd.Series(
-        np.degrees(np.arctan2(projected_normal[:, 1], projected_normal[:, 0])),
-        index=times)
+    surface_azimuth = \
+        np.degrees(np.arctan2(projected_normal[:, 1], projected_normal[:, 0]))
 
     # 2. Clean up atan when x-coord or y-coord is zero
 #     surface_azimuth[(projected_normal[:,0]==0) & (projected_normal[:,1]>0)] =  90
@@ -545,18 +544,17 @@ def singleaxis(apparent_zenith, apparent_azimuth,
     surface_azimuth[surface_azimuth >= 360] -= 360
 
     # Calculate surface_tilt
-    # Use pandas to calculate the sum because it handles nan values better.
-    surface_tilt = (90 - np.degrees(np.arccos(
-                            pd.DataFrame(panel_norm_earth * projected_normal,
-                                         index=times).sum(axis=1))))
+    dotproduct = (panel_norm_earth * projected_normal).sum(axis=1)
+    surface_tilt = 90 - np.degrees(np.arccos(dotproduct))
 
     # Bundle DataFrame for return values and filter for sun below horizon.
-    df_out = pd.DataFrame({'tracker_theta': tracker_theta, 'aoi': aoi,
-                           'surface_azimuth': surface_azimuth,
-                           'surface_tilt': surface_tilt},
-                          index=times)
-    df_out = df_out[['tracker_theta', 'aoi',
-                     'surface_azimuth', 'surface_tilt']]
-    df_out[apparent_zenith > 90] = np.nan
-
-    return df_out
+    out = {'tracker_theta': tracker_theta, 'aoi': aoi,
+           'surface_azimuth': surface_azimuth, 'surface_tilt': surface_tilt}
+    if index is not None:
+        out = pd.DataFrame(out, index=index)
+        out = out[['tracker_theta', 'aoi', 'surface_azimuth', 'surface_tilt']]
+        out[zen_gt_90] = np.nan
+    else:
+        out = {k: np.where(zen_gt_90, np.nan, v) for k, v in out.items()}
+
+    return out


Test patch:
diff --git a/pvlib/test/test_tracking.py b/pvlib/test/test_tracking.py
index c8c6c09..d54e85d 100644
--- a/pvlib/test/test_tracking.py
+++ b/pvlib/test/test_tracking.py
@@ -14,9 +14,11 @@ from pvlib import tracking
 SINGLEAXIS_COL_ORDER = ['tracker_theta', 'aoi',
                         'surface_azimuth', 'surface_tilt']
 
+
 def test_solar_noon():
-    apparent_zenith = pd.Series([10])
-    apparent_azimuth = pd.Series([180])
+    index = pd.DatetimeIndex(start='20180701T1200', freq='1s', periods=1)
+    apparent_zenith = pd.Series([10], index=index)
+    apparent_azimuth = pd.Series([180], index=index)
     tracker_data = tracking.singleaxis(apparent_zenith, apparent_azimuth,
                                        axis_tilt=0, axis_azimuth=0,
                                        max_angle=90, backtrack=True,
@@ -24,12 +26,90 @@ def test_solar_noon():
 
     expect = pd.DataFrame({'tracker_theta': 0, 'aoi': 10,
                            'surface_azimuth': 90, 'surface_tilt': 0},
-                           index=[0], dtype=np.float64)
+                          index=index, dtype=np.float64)
     expect = expect[SINGLEAXIS_COL_ORDER]
 
     assert_frame_equal(expect, tracker_data)
 
 
+def test_scalars():
+    apparent_zenith = 10
+    apparent_azimuth = 180
+    tracker_data = tracking.singleaxis(apparent_zenith, apparent_azimuth,
+                                       axis_tilt=0, axis_azimuth=0,
+                                       max_angle=90, backtrack=True,
+                                       gcr=2.0/7.0)
+    assert isinstance(tracker_data, dict)
+    expect = {'tracker_theta': 0, 'aoi': 10, 'surface_azimuth': 90,
+              'surface_tilt': 0}
+    for k, v in expect.items():
+        assert_allclose(tracker_data[k], v)
+
+
+def test_arrays():
+    apparent_zenith = np.array([10])
+    apparent_azimuth = np.array([180])
+    tracker_data = tracking.singleaxis(apparent_zenith, apparent_azimuth,
+                                       axis_tilt=0, axis_azimuth=0,
+                                       max_angle=90, backtrack=True,
+                                       gcr=2.0/7.0)
+    assert isinstance(tracker_data, dict)
+    expect = {'tracker_theta': 0, 'aoi': 10, 'surface_azimuth': 90,
+              'surface_tilt': 0}
+    for k, v in expect.items():
+        assert_allclose(tracker_data[k], v)
+
+
+def test_nans():
+    apparent_zenith = np.array([10, np.nan, 10])
+    apparent_azimuth = np.array([180, 180, np.nan])
+    with np.errstate(invalid='ignore'):
+        tracker_data = tracking.singleaxis(apparent_zenith, apparent_azimuth,
+                                           axis_tilt=0, axis_azimuth=0,
+                                           max_angle=90, backtrack=True,
+                                           gcr=2.0/7.0)
+    expect = {'tracker_theta': np.array([0, nan, nan]),
+              'aoi': np.array([10, nan, nan]),
+              'surface_azimuth': np.array([90, nan, nan]),
+              'surface_tilt': np.array([0, nan, nan])}
+    for k, v in expect.items():
+        assert_allclose(tracker_data[k], v)
+
+    # repeat with Series because nans can differ
+    apparent_zenith = pd.Series(apparent_zenith)
+    apparent_azimuth = pd.Series(apparent_azimuth)
+    with np.errstate(invalid='ignore'):
+        tracker_data = tracking.singleaxis(apparent_zenith, apparent_azimuth,
+                                           axis_tilt=0, axis_azimuth=0,
+                                           max_angle=90, backtrack=True,
+                                           gcr=2.0/7.0)
+    expect = pd.DataFrame(np.array(
+        [[ 0., 10., 90.,  0.],
+         [nan, nan, nan, nan],
+         [nan, nan, nan, nan]]),
+        columns=['tracker_theta', 'aoi', 'surface_azimuth', 'surface_tilt'])
+    assert_frame_equal(tracker_data, expect)
+
+
+def test_arrays_multi():
+    apparent_zenith = np.array([[10, 10], [10, 10]])
+    apparent_azimuth = np.array([[180, 180], [180, 180]])
+    # singleaxis should fail for num dim > 1
+    with pytest.raises(ValueError):
+        tracker_data = tracking.singleaxis(apparent_zenith, apparent_azimuth,
+                                           axis_tilt=0, axis_azimuth=0,
+                                           max_angle=90, backtrack=True,
+                                           gcr=2.0/7.0)
+    # uncomment if we ever get singleaxis to support num dim > 1 arrays
+    # assert isinstance(tracker_data, dict)
+    # expect = {'tracker_theta': np.full_like(apparent_zenith, 0),
+    #           'aoi': np.full_like(apparent_zenith, 10),
+    #           'surface_azimuth': np.full_like(apparent_zenith, 90),
+    #           'surface_tilt': np.full_like(apparent_zenith, 0)}
+    # for k, v in expect.items():
+    #     assert_allclose(tracker_data[k], v)
+
+
 def test_azimuth_north_south():
     apparent_zenith = pd.Series([60])
     apparent_azimuth = pd.Series([90])
@@ -163,14 +243,38 @@ def test_axis_azimuth():
     assert_frame_equal(expect, tracker_data)
 
 
-def test_index_mismatch():
-    apparent_zenith = pd.Series([30])
-    apparent_azimuth = pd.Series([90,180])
-    with pytest.raises(ValueError):
-        tracker_data = tracking.singleaxis(apparent_zenith, apparent_azimuth,
-                                           axis_tilt=0, axis_azimuth=90,
-                                           max_angle=90, backtrack=True,
-                                           gcr=2.0/7.0)
+def test_horizon_flat():
+    # GH 569
+    solar_azimuth = np.array([0, 180, 359])
+    solar_zenith = np.array([100, 45, 100])
+    solar_azimuth = pd.Series(solar_azimuth)
+    solar_zenith = pd.Series(solar_zenith)
+    # depending on platform and numpy versions this will generate
+    # RuntimeWarning: invalid value encountered in > < >=
+    out = tracking.singleaxis(solar_zenith, solar_azimuth, axis_tilt=0,
+                              axis_azimuth=180, backtrack=False, max_angle=180)
+    expected = pd.DataFrame(np.array(
+        [[ nan,  nan,  nan,  nan],
+         [  0.,  45., 270.,   0.],
+         [ nan,  nan,  nan,  nan]]),
+        columns=['tracker_theta', 'aoi', 'surface_azimuth', 'surface_tilt'])
+    assert_frame_equal(out, expected)
+
+
+def test_horizon_tilted():
+    # GH 569
+    solar_azimuth = np.array([0, 180, 359])
+    solar_zenith = np.full_like(solar_azimuth, 45)
+    solar_azimuth = pd.Series(solar_azimuth)
+    solar_zenith = pd.Series(solar_zenith)
+    out = tracking.singleaxis(solar_zenith, solar_azimuth, axis_tilt=90,
+                              axis_azimuth=180, backtrack=False, max_angle=180)
+    expected = pd.DataFrame(np.array(
+        [[ 180.,  45.,   0.,  90.],
+         [   0.,  45., 180.,  90.],
+         [ 179.,  45., 359.,  90.]]),
+        columns=['tracker_theta', 'aoi', 'surface_azimuth', 'surface_tilt'])
+    assert_frame_equal(out, expected)
 
 
 def test_SingleAxisTracker_creation():
@@ -285,19 +389,25 @@ def test_get_irradiance():
                              end='20160101 1800-0700', freq='6H')
     location = Location(latitude=32, longitude=-111)
     solar_position = location.get_solarposition(times)
-    irrads = pd.DataFrame({'dni':[900,0], 'ghi':[600,0], 'dhi':[100,0]},
+    irrads = pd.DataFrame({'dni': [900, 0], 'ghi': [600, 0], 'dhi': [100, 0]},
                           index=times)
     solar_zenith = solar_position['apparent_zenith']
     solar_azimuth = solar_position['azimuth']
-    tracker_data = system.singleaxis(solar_zenith, solar_azimuth)
-
-    irradiance = system.get_irradiance(tracker_data['surface_tilt'],
-                                       tracker_data['surface_azimuth'],
-                                       solar_zenith,
-                                       solar_azimuth,
-                                       irrads['dni'],
-                                       irrads['ghi'],
-                                       irrads['dhi'])
+
+    # invalid warnings already generated in horizon test above,
+    # no need to clutter test output here
+    with np.errstate(invalid='ignore'):
+        tracker_data = system.singleaxis(solar_zenith, solar_azimuth)
+
+    # some invalid values in irradiance.py. not our problem here
+    with np.errstate(invalid='ignore'):
+        irradiance = system.get_irradiance(tracker_data['surface_tilt'],
+                                           tracker_data['surface_azimuth'],
+                                           solar_zenith,
+                                           solar_azimuth,
+                                           irrads['dni'],
+                                           irrads['ghi'],
+                                           irrads['dhi'])
 
     expected = pd.DataFrame(data=np.array(
         [[961.80070,   815.94490,   145.85580,   135.32820, 10.52757492],