more mode -1 changes

docs/examples/snow-detection/snow-mode.py

@@ -125,8 +125,11 @@
 data[dc_current_cols] = data[dc_current_cols].replace({np.nan: 0, None: 0})
 data.loc[:, ac_power_cols] = data[ac_power_cols].replace({np.nan: 0, None: 0})
 
-# %% Plot DC voltage for each combiner input relative to inverter nameplate
+
+# %% 
+# Plot DC voltage for each combiner input relative to inverter nameplate
 # limits
+
 fig, ax = plt.subplots(figsize=(10, 10))
 date_form = DateFormatter("%m/%d")
 ax.xaxis.set_major_formatter(date_form)
@@ -140,6 +143,7 @@
 ax.set_xlabel('Date', fontsize='large')
 ax.set_ylabel('Voltage [V]', fontsize='large')
 ax.legend(loc='lower left')
+plt.show()
 
 # %% Plot AC power relative to inverter nameplate limits
 
@@ -154,6 +158,7 @@
 ax.set_xlabel('Date', fontsize='large')
 ax.set_ylabel('AC Power [kW]', fontsize='large')
 ax.legend(loc='upper left')
+plt.show()
 
 # %% Filter data.
 # Identify periods where the system is operating off of its maximum power
@@ -209,6 +214,7 @@
 ax.set_xlabel('Date', fontsize='large')
 ax.set_ylabel('Voltage [V]', fontsize='large')
 ax.legend(loc='lower left')
+plt.show()
 
 # %% We want to exclude periods where array voltage is affected by horizon
 # shading
@@ -235,6 +241,7 @@
 ax.legend()
 ax.set_xlabel(r'Azimuth [$\degree$]')
 ax.set_ylabel(r'Elevation [$\degree$]')
+plt.show()
 
 # %% Exclude data collected while the sun is below the horizon
 data = data[data['Horizon Mask']]
@@ -267,6 +274,7 @@
 ax.plot(data['Cell Temp [C]'], alpha=0.3, c='b')
 ax.set_ylabel('Cell Temp [C]', c='b', fontsize='xx-large')
 ax.set_xlabel('Date', fontsize='xx-large')
+plt.show()
 
 # %% For one combiner, demonstrate the transmission calculation using two
 # different approaches to modeling effective irradiance from measured Imp.
@@ -315,6 +323,7 @@
 ax.legend()
 ax.set_ylabel('Transmission', fontsize='xx-large')
 ax.set_xlabel('Date + Time', fontsize='large')
+plt.show()
 
 # %%
 # Model voltage using calculated transmission (two different approaches)
@@ -369,6 +378,7 @@
 ax.legend(fontsize='xx-large')
 ax.set_ylabel('Voltage [V]', fontsize='xx-large')
 ax.set_xlabel('Date', fontsize='large')
+plt.show()
 
 
 # %% Function to do analysis so we can loop over combiner boxes
@@ -433,35 +443,25 @@ def wrapper(voltage, current, temp_cell, effective_irradiance,
     T = snow.get_transmission(effective_irradiance, modeled_e_e,
                               current/config['num_str_per_cb'])
 
-    name_T = inv_cb + ' Transmission'
-    data[name_T] = T
-
     # Model voltage for a single module, scale up to array
-    modeled_vmp = pvlib.pvsystem.sapm(effective_irradiance*T, temp_cell,
-                                      coeffs)['v_mp']
-    modeled_vmp *= config['num_mods_per_str']
-
-    # Voltage is modeled as NaN if T = 0, but V = 0 makes more sense
-    modeled_vmp[T == 0] = 0
-
-    # Identify periods where modeled voltage is outside of MPPT range,
-    # and correct values
-    modeled_vmp[modeled_vmp > config['max_dcv']] = np.nan
-    modeled_vmp[modeled_vmp < config['min_dcv']] = 0
-
-    # Calculate voltage ratio
-    with np.errstate(divide='ignore'):
-        vmp_ratio = voltage / modeled_vmp
-
-    # take care of divide by zero
-    vmp_ratio[modeled_vmp == 0] = np.nan
+    modeled_voltage_with_calculated_transmission =\
+        pvlib.pvsystem.sapm(effective_irradiance*T, temp_cell,
+                            coeffs)['v_mp'] * config['num_mods_per_str']
+    modeled_voltage_with_ideal_transmission =\
+        pvlib.pvsystem.sapm(effective_irradiance, temp_cell,
+                            coeffs)['v_mp'] * config['num_mods_per_str']
+
+    mode, vmp_ratio = snow.categorize(
+        T, voltage, modeled_voltage_with_calculated_transmission,
+        modeled_voltage_with_ideal_transmission, config['min_dcv'],
+        config['max_dcv'], config['threshold_vratio'],
+        config['threshold_transmission'])
 
-    mode = snow.categorize(vmp_ratio, T, voltage, modeled_vmp,
-                           config['min_dcv'],
-                           config['threshold_vratio'],
-                           config['threshold_transmission'])
     my_dict = {'transmission': T,
-               'modeled_vmp': modeled_vmp,
+               'modeled_voltage_with_calculated_transmission':
+               modeled_voltage_with_calculated_transmission,
+               'modeled_voltage_with_ideal_transmission':
+               modeled_voltage_with_ideal_transmission,
                'vmp_ratio': vmp_ratio,
                'mode': mode}
 
@@ -544,7 +544,8 @@ def wrapper(voltage, current, temp_cell, effective_irradiance,
 # %%
 # Look at transmission for all DC inputs
 
-transmission_cols = [c for c in data.columns if 'transmission' in c]
+transmission_cols = [c for c in data.columns if 'transmission' in c and
+                     'voltage' not in c]
 fig, ax = plt.subplots(figsize=(10, 10))
 date_form = DateFormatter("%m/%d")
 ax.xaxis.set_major_formatter(date_form)
@@ -554,6 +555,7 @@ def wrapper(voltage, current, temp_cell, effective_irradiance,
     ax.scatter(temp.index, temp[c], s=0.5, label=c)
 ax.set_xlabel('Date', fontsize='xx-large')
 ax.legend()
+plt.show()
 
 # %%
 # Look at voltage ratios for all DC inputs
@@ -571,6 +573,7 @@ def wrapper(voltage, current, temp_cell, effective_irradiance,
 ax.set_ylabel('Voltage Ratio (measured/modeled)', fontsize='xx-large')
 ax.axhline(1, c='k', alpha=0.1, ls='--')
 ax.legend()
+plt.show()
 
 # %% Calculate all power losses - snow and non-snow
 
@@ -664,6 +667,7 @@ def wrapper(voltage, current, temp_cell, effective_irradiance,
 ax.legend(handles=handles, fontsize='xx-large', loc='upper left')
 ax.set_title('Measured and modeled production for INV1 CB2',
              fontsize='xx-large')
+plt.show()
 
 # %%
 
@@ -714,5 +718,6 @@ def wrapper(voltage, current, temp_cell, effective_irradiance,
 ax.set_xticks(xvals, days)
 ax.xaxis.set_major_formatter(date_form)
 ax.set_title('Losses incurred in modes -1, 0, 1, 2, 3', fontsize='xx-large')
+plt.show()
 
 # %%

pvanalytics/features/snow.py

@@ -147,8 +147,10 @@ def get_transmission(measured_e_e, modeled_e_e, i_mp):
     return T
 
 
-def categorize(vmp_ratio, transmission, measured_voltage, modeled_voltage,
-               min_dcv, threshold_vratio, threshold_transmission):
+def categorize(transmission, measured_voltage,
+               modeled_voltage_with_calculated_transmission,
+               modeled_voltage_with_ideal_transmission,
+               min_dcv, max_dcv, threshold_vratio, threshold_transmission):
 
     """
     Categorizes electrical behavior into a snow-related mode.
@@ -173,9 +175,14 @@ def categorize(vmp_ratio, transmission, measured_voltage, modeled_voltage,
 
     An additional mode was added to cover a case that was not addressed in
     [1]_:
-    * Mode -1: Indicates periods where voltage modeled with measured irradiance
-      is below the inverter's turn-on voltage. Under Mode -1, it is unknown if
-      and how snow impacts power output.
+    * Mode -1: Indicates periods where it is unknown if and how snow impacts
+    power output. Mode -1 includes periods when
+      (a) voltage modeled with measured irradiance is below the inverter's
+      turn-on voltage OR
+      (b) voltage modeled with measured irradiance exceeds the upper bound of
+      the inverter's MPPT algorithm OR
+      (c) measured voltage exceeds the upper bound of the inverter's MPPT
+      algorithm
 
     Parameters
     ----------
@@ -187,12 +194,18 @@ def categorize(vmp_ratio, transmission, measured_voltage, modeled_voltage,
         through the snow cover. [dimensionless]
     measured_voltage : array-like
         Measured DC voltage. [V]
-    modeled_voltage : array-like
+    modeled_voltage_with_calculated_transmission : array-like
+        DC voltage modeled using measured plane-of-array irradiance,
+        back-of-module temperature, and calculated transmission. [V]
+    modeled_voltage_with_ideal_transmission : array-like
         DC voltage modeled using measured plane-of-array irradiance and
-        back-of-module temperature. [V]
+        back-of-module temperature. Assumes transmission equals 1. [V]
     min_dcv : float
         The lower voltage bound on the inverter's maximum power point
         tracking (MPPT) algorithm. [V]
+    max_dcv : float
+        The upper voltage bound on the inverter's maximum power point
+        tracking (MPPT) algorithm. [V]
     threshold_vratio : float
         The lower bound for vratio that is representative of snow-free
         conditions. Determined empirically. Depends on system configuration and
@@ -211,28 +224,42 @@ def categorize(vmp_ratio, transmission, measured_voltage, modeled_voltage,
        2023, pp. 1-5, :doi:`10.1109/PVSC48320.2023.10360065`.
     """
     umin_meas = measured_voltage >= min_dcv
-    umin_model = modeled_voltage >= min_dcv
+    umax_meas = measured_voltage < max_dcv
+
+    umin_model = modeled_voltage_with_ideal_transmission >= min_dcv
+    umax_model = modeled_voltage_with_ideal_transmission < max_dcv
+
+    # Voltage is modeled as NaN if T = 0, but V = 0 makes more sense
+    modeled_voltage_with_calculated_transmission[transmission == 0] = 0
+
+    with np.errstate(divide='ignore'):
+        vmp_ratio =\
+            measured_voltage / modeled_voltage_with_calculated_transmission
 
-    # offline if model says that voltage is too low.
-    # if modeled is above the minimum, then system is
-    # possibly generating
-    # offline = ~umin_meas & ~umin_model
+    # take care of divide by zero
+    vmp_ratio[modeled_voltage_with_calculated_transmission == 0] = 1
 
     # vmp_ratio discrimates between states (1,2) and (3,4)
     uvr = np.where(vmp_ratio >= threshold_vratio, 3, 1)
 
     # transmission discrimates within (1,2) and (3,4)
     utrans = np.where(transmission >= threshold_transmission, 1, 0)
 
-    # None if nan
-    # if offline:
-    #   - -1 if umin_model is 0
-    # if not offline:
-    #   - 0 if umin_meas is 0, i.e., measurement indicate no power but
-    #     it must be that umin_model is 1
-    #   - state 1, 2, 3, 4 defined by uvr + utrans
-    mode = np.where(np.isnan(vmp_ratio) | np.isnan(transmission),
-                    None, umin_meas * (uvr + utrans))
-    mode = np.where(~umin_model, -1, mode)
-
-    return mode
+    # None if transmission, vmp_ratio, modeled_voltage_with_ideal_transmission,
+    # modeled_voltage_with_calculated_transmission, or measured_voltage is nan
+    # -1 if umin_model is 0
+    # 0 if umin_meas is 0, i.e., measurement indicate no power but it must be
+    # that umin_model is 1
+    # state 1, 2, 3, 4 defined by uvr + utrans
+    # -1 if umax_model is 0
+    # -1 if umax_meas is 0
+
+    mode = umin_meas * (uvr + utrans)
+    mode = np.where(~umin_model | ~umax_model | ~umax_meas, -1, mode)
+    mode = np.where(np.isnan(vmp_ratio) | np.isnan(transmission) |
+                    np.isnan(measured_voltage) |
+                    np.isnan(modeled_voltage_with_calculated_transmission) |
+                    np.isnan(modeled_voltage_with_ideal_transmission),
+                    None, mode)
+
+    return mode, vmp_ratio

pvanalytics/tests/features/test_snow.py

@@ -66,25 +66,27 @@ def test_get_transmission():
 
 
 def test_categorize():
-    vmp_ratio = np.array([np.nan, 0.9, 0.1, 0.6, 0.7, 0.9, 0.9, 0.5])
     measured_voltage = np.array([400., 450., 400., 420., 420., 495., 270.,
-                                 275.])
-    modeled_voltage = np.array([np.nan, 500, 4000, 700, 600, 550, 300, 550])
-    transmission = np.array([0.5, np.nan, 0.5, 0.9, 0.5, 0.9, 0.9, 0.9])
+                                 200., 850., 0., 700.])
+    modeled_voltage_with_calculated_transmission = np.array([np.nan, 500, 4000,
+                                                             700, 600, 550,
+                                                             300, 200, 600,
+                                                             0, 700])
+    modeled_voltage_with_ideal_transmission = np.array([700, 600, 6000, 750,
+                                                        700, 700, 350, 250,
+                                                        600, 400, 850])
+    transmission = np.array([0.5, np.nan, 0.5, 0.9, 0.5, 0.9, 0.9, 0.9, 1, 0,
+                             0.9])
+
     min_dcv = 300
+    max_dcv = 800
     threshold_vratio = 0.7
     threshold_transmission = 0.6
-    # vmp_ratio: np.nan, >thres, <thres, <thres, =thres, >thres, >thres, <thres
-    # measured_voltage: >thres, >thres, >thres, >thres, >thres, >thres, <thres,
-    #                   <thres
-    # modeled_voltage: np.nan, >thres, >thres, >thres, >thres, >thres, <thres,
-    #                  >thres
-    # transmission: <thres, np.nan, <thres, >thres, <thres, >thres, >thres,
-    #               >thres, > thres
-    # None (vmp_ratio, modeled_voltage), None (transmission), 1, 2, 3, 4, -1,
-    # 0
-    expected = np.array([None, None, 1, 2, 3, 4, -1, 0])
-    result = snow.categorize(vmp_ratio, transmission, measured_voltage,
-                             modeled_voltage, min_dcv,
-                             threshold_vratio, threshold_transmission)
+
+    expected = np.array([None, None, 1, 2, 3, 4, 0, -1, -1, 0, -1])
+    result, _ = snow.categorize(transmission, measured_voltage,
+                                modeled_voltage_with_calculated_transmission,
+                                modeled_voltage_with_ideal_transmission,
+                                min_dcv, max_dcv, threshold_vratio,
+                                threshold_transmission)
     assert_array_equal(result, expected)