Add N. Martin & J. M. Ruiz spectral response mismatch modifiers model

docs/examples/spectrum/plot_martin_ruiz_mismatch.py

@@ -0,0 +1,219 @@
+"""
+N. Martin & J. M. Ruiz Spectral Mismatch Modifier
+=================================================
+
+How to use this correction factor to adjust the POA global irradiance.
+"""
+
+# %%
+# Effectiveness of a material to convert incident sunlight to current depends
+# on the incident light wavelength. During the day, the spectral distribution
+# of the incident irradiance varies from the standard testing spectra,
+# introducing a small difference between the expected and the real output.
+# In [1]_, N. Martín and J. M. Ruiz propose 3 mismatch factors, one for each
+# irradiance component. These mismatch modifiers are calculated with the help
+# of the airmass, the clearness index and three experimental fitting
+# parameters. In the same paper, these parameters have been obtained for m-Si,
+# p-Si and a-Si modules.
+# With :py:func:`pvlib.spectrum.martin_ruiz` we are able to make use of these
+# already computed values or provide ours.
+#
+# References
+# ----------
+#  .. [1] Martín, N. and Ruiz, J.M. (1999), A new method for the spectral
+#     characterisation of PV modules. Prog. Photovolt: Res. Appl., 7: 299-310.
+#     :doi:`10.1002/(SICI)1099-159X(199907/08)7:4<299::AID-PIP260>3.0.CO;2-0`
+#
+# Calculating the incident and modified global irradiance
+# -------------------------------------------------------
+#
+# Mismatch modifiers are applied to the irradiance components, so first
+# step is to get them. We define an hypothetical POA surface and use TMY to
+# compute sky diffuse, ground reflected and direct irradiance.
+
+import os
+
+import matplotlib.pyplot as plt
+import pvlib
+from scipy import stats
+import pandas as pd
+
+surface_tilt = 40
+surface_azimuth = 180  # Pointing South
+
+# Get TMY data & create location
+datapath = os.path.join(pvlib.__path__[0], 'data',
+                        'tmy_45.000_8.000_2005_2016.csv')
+pvgis_data, _, metadata, _ = pvlib.iotools.read_pvgis_tmy(datapath,
+                                                          map_variables=True)
+site = pvlib.location.Location(metadata['latitude'], metadata['longitude'],
+                               altitude=metadata['elevation'])
+
+# Coerce a year: above function returns typical months of different years
+pvgis_data.index = [ts.replace(year=2022) for ts in pvgis_data.index]
+# Select days to show
+weather_data = pvgis_data['2022-09-03':'2022-09-06']
+
+# Then calculate all we need to get the irradiance components
+solar_pos = site.get_solarposition(weather_data.index)
+
+extra_rad = pvlib.irradiance.get_extra_radiation(weather_data.index)
+
+poa_sky_diffuse = pvlib.irradiance.haydavies(surface_tilt, surface_azimuth,
+                                             weather_data['dhi'],
+                                             weather_data['dni'],
+                                             extra_rad,
+                                             solar_pos['apparent_zenith'],
+                                             solar_pos['azimuth'])
+
+poa_ground_diffuse = pvlib.irradiance.get_ground_diffuse(surface_tilt,
+                                                         weather_data['ghi'])
+
+aoi = pvlib.irradiance.aoi(surface_tilt, surface_azimuth,
+                           solar_pos['apparent_zenith'], solar_pos['azimuth'])
+
+# Get dataframe with all components and global (includes 'poa_direct')
+poa_irrad = pvlib.irradiance.poa_components(aoi, weather_data['dni'],
+                                            poa_sky_diffuse,
+                                            poa_ground_diffuse)
+
+# %%
+# Here come the modifiers. Let's calculate them with the airmass and clearness
+# index.
+# First, let's find the airmass and the clearness index.
+# Little caution: default values for this model were fitted obtaining the
+# airmass through the `'kasten1966'` method, which is not used by default.
+
+airmass = site.get_airmass(solar_position=solar_pos, model='kasten1966')
+airmass_absolute = airmass['airmass_absolute']  # We only use absolute airmass
+clearness_index = pvlib.irradiance.clearness_index(weather_data['ghi'],
+                                                   solar_pos['zenith'],
+                                                   extra_rad)
+
+# Get the spectral mismatch modifiers
+spectral_modifiers = pvlib.spectrum.martin_ruiz(clearness_index,
+                                                airmass_absolute,
+                                                module_type='monosi')
+
+# %%
+# And then we can find the 3 modified components of the POA irradiance
+# by means of a simple multiplication.
+# Note, however, that this does not modify ``poa_global`` nor
+# ``poa_diffuse``, so we should update the dataframe afterwards.
+
+poa_irrad_modified = poa_irrad * spectral_modifiers
+# Above line is equivalent to:
+# poa_irrad_modified = pd.DataFrame()
+# for component in ('poa_direct', 'poa_sky_diffuse', 'poa_ground_diffuse'):
+#     poa_irrad_modified[component] = (poa_irrad[component]
+#                                      * spectral_modifiers[component])
+
+# We want global modified irradiance
+poa_irrad_modified['poa_global'] = (poa_irrad_modified['poa_direct']
+                                    + poa_irrad_modified['poa_sky_diffuse']
+                                    + poa_irrad_modified['poa_ground_diffuse'])
+# Don't forget to update `'poa_diffuse'` if you want to use it
+# poa_irrad_modified['poa_diffuse'] = \
+#     (poa_irrad_modified['poa_sky_diffuse']
+#      + poa_irrad_modified['poa_ground_diffuse'])
+
+# %%
+# Finally, let's plot the incident vs modified global irradiance, and their
+# difference.
+
+poa_irrad_global_diff = (poa_irrad['poa_global']
+                         - poa_irrad_modified['poa_global'])
+plt.figure()
+datetimes = poa_irrad.index  # common to poa_irrad_*
+plt.plot(datetimes, poa_irrad['poa_global'].to_numpy())
+plt.plot(datetimes, poa_irrad_modified['poa_global'].to_numpy())
+plt.plot(datetimes, poa_irrad_global_diff.to_numpy())
+plt.legend(['Incident', 'Modified', 'Difference'])
+plt.ylabel('POA Global irradiance [W/m²]')
+plt.grid()
+plt.show()
+
+# %%
+# Comparison with other models
+# ----------------------------
+# During the addition of this model, a question arose about its trustworthiness
+# so, in order to check the integrity of the implementation, we will
+# compare it against :py:func:`pvlib.pvsystem.sapm_spectral_loss` and
+# :py:func:`pvlib.atmosphere.first_solar`.
+# Former model needs the parameters that characterise a module, but which one?
+# We will take the mean of Sandia parameters `'A0', 'A1', 'A2', 'A3', 'A4'` for
+# the same material type.
+# On the other hand, :py:func:`~pvlib.atmosphere.first_solar` needs the
+# precipitable water. We assume the standard spectrum, `1.42 cm`.
+
+# Retrieve modules and select the subset we want to work with the SAPM model
+module_type = 'mc-Si'  # Equivalent to monosi
+sandia_modules = pvlib.pvsystem.retrieve_sam(name='SandiaMod')
+modules_subset = \
+    sandia_modules.loc[:, sandia_modules.loc['Material'] == module_type]
+
+# Define typical module and get the means of the A0 to A4 parameters
+modules_aggregated = pd.DataFrame(index=('mean', 'std'))
+for param in ('A0', 'A1', 'A2', 'A3', 'A4'):
+    result, _, _ = stats.mvsdist(modules_subset.loc[param])
+    modules_aggregated[param] = result.mean(), result.std()
+
+# Check if 'mean' is a representative value with help of 'std' just in case
+print(modules_aggregated)
+
+# Then apply the SAPM model and calculate introduced difference
+modifier_sapm_f1 = \
+    pvlib.pvsystem.sapm_spectral_loss(airmass_absolute,
+                                      modules_aggregated.loc['mean'])
+poa_irrad_sapm_modified = poa_irrad['poa_global'] * modifier_sapm_f1
+poa_irrad_sapm_difference = (poa_irrad['poa_global']
+                             - poa_irrad_sapm_modified)
+
+# atmosphere.first_solar model
+first_solar_pw = 1.42  # Default for AM1.5 spectrum
+modifier_first_solar = \
+    pvlib.atmosphere.first_solar_spectral_correction(first_solar_pw,
+                                                     airmass_absolute,
+                                                     module_type='monosi')
+poa_irrad_first_solar_mod = poa_irrad['poa_global'] * modifier_first_solar
+poa_irrad_first_solar_diff = (poa_irrad['poa_global']
+                              - poa_irrad_first_solar_mod)
+
+# %%
+# Plot global irradiance difference over time
+# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+datetimes = poa_irrad_global_diff.index  # common to poa_irrad_*_diff*
+plt.figure()
+plt.plot(datetimes, poa_irrad_global_diff.to_numpy(),
+         label='spectrum.martin_ruiz')
+plt.plot(datetimes, poa_irrad_sapm_difference.to_numpy(),
+         label='atmosphere.first_solar')
+plt.plot(datetimes, poa_irrad_first_solar_diff.to_numpy(),
+         label='pvsystem.sapm_spectral_loss')
+plt.legend()
+plt.title('Introduced difference comparison of different models')
+plt.ylabel('POA Global Irradiance Difference [W/m²]')
+plt.grid()
+plt.show()
+
+# %%
+# Plot modifier vs absolute airmass
+# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+ama = airmass_absolute.to_numpy()
+# spectrum.martin_ruiz has 3 modifiers, so we only calculate one as
+# M = S_eff / S_incident that takes into account the global effect
+martin_ruiz_agg_modifier = (poa_irrad_modified['poa_global']
+                            / poa_irrad['poa_global'])
+plt.figure()
+plt.scatter(ama, martin_ruiz_agg_modifier.to_numpy(),
+            label='spectrum.martin_ruiz')
+plt.scatter(ama, modifier_sapm_f1.to_numpy(),
+            label='pvsystem.sapm_spectral_loss')
+plt.scatter(ama, modifier_first_solar.to_numpy(),
+            label='atmosphere.first_solar')
+plt.legend()
+plt.title('Introduced difference comparison of different models')
+plt.xlabel('Absolute airmass')
+plt.ylabel(r'Modifier $M = \frac{S_{effective}}{S_{incident}}$')
+plt.grid()
+plt.show()

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

@@ -10,3 +10,4 @@ Spectrum
    spectrum.get_example_spectral_response
    spectrum.get_am15g
    spectrum.calc_spectral_mismatch_field
+   spectrum.martin_ruiz

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

@@ -20,6 +20,9 @@ Enhancements
   :py:func:`pvlib.snow.loss_townsend` (:issue:`1636`, :pull:`1653`)
 * Added optional ``n_ar`` parameter to :py:func:`pvlib.iam.physical` to
   support an anti-reflective coating. (:issue:`1501`, :pull:`1616`)
+* Added :py:func:`pvlib.spectrum.martin_ruiz`, a spectral
+  mismatch correction factor for POA components, dependant in module type,
+  airmass and clearness index (:pull:`1658`)
 
 Bug fixes
 ~~~~~~~~~
@@ -64,3 +67,4 @@ Contributors
 * Mark Mikofski (:ghuser:`mikofski`)
 * Anton Driesse (:ghuser:`adriesse`)
 * Michael Deceglie (:ghuser:`mdeceglie`)
+* Echedey Luis (:ghuser:`echedey-ls`)

pvlib/atmosphere.py

@@ -148,7 +148,7 @@ def get_relative_airmass(zenith, model='kastenyoung1989'):
         Available models include the following:
 
         * 'simple' - secant(apparent zenith angle) -
-          Note that this gives -Inf at zenith=90
+          Note that this gives +Inf at zenith=90
         * 'kasten1966' - See reference [1] -
           requires apparent sun zenith
         * 'youngirvine1967' - See reference [2] -

pvlib/ivtools/sdm.py

@@ -81,7 +81,7 @@ def fit_cec_sam(celltype, v_mp, i_mp, v_oc, i_sc, alpha_sc, beta_voc,
 
     Notes
     -----
-    The CEC model and estimation method  are described in [1]_.
+    The CEC model and estimation method are described in [1]_.
     Inputs ``v_mp``, ``i_mp``, ``v_oc`` and ``i_sc`` are assumed to be from a
     single IV curve at constant irradiance and cell temperature. Irradiance is
     not explicitly used by the fitting procedure. The irradiance level at which

pvlib/pvsystem.py

@@ -2497,8 +2497,8 @@ def sapm(effective_irradiance, temp_cell, module):
                        Imp, Vmp, Ix, and Ixx to effective irradiance
     Isco               Short circuit current at reference condition (amps)
     Impo               Maximum power current at reference condition (amps)
-    Voco               Open circuit voltage at reference condition (amps)
-    Vmpo               Maximum power voltage at reference condition (amps)
+    Voco               Open circuit voltage at reference condition (volts)
+    Vmpo               Maximum power voltage at reference condition (volts)
     Aisc               Short circuit current temperature coefficient at
                        reference condition (1/C)
     Aimp               Maximum power current temperature coefficient at
@@ -2674,10 +2674,7 @@ def sapm_effective_irradiance(poa_direct, poa_diffuse, airmass_absolute, aoi,
     and diffuse irradiance on the plane of array to the irradiance absorbed by
     a module's cells.
 
-    The model is
-    .. math::
-
-        `Ee = f_1(AM_a) (E_b f_2(AOI) + f_d E_d)`
+    The model is :math:`Ee = f_1(AM_a) (E_b f_2(AOI) + f_d E_d)`
 
     where :math:`Ee` is effective irradiance (W/m2), :math:`f_1` is a fourth
     degree polynomial in air mass :math:`AM_a`, :math:`E_b` is beam (direct)

pvlib/spectrum/__init__.py

@@ -1,3 +1,4 @@
 from pvlib.spectrum.spectrl2 import spectrl2  # noqa: F401
 from pvlib.spectrum.mismatch import (get_example_spectral_response, get_am15g,
-                                     calc_spectral_mismatch_field)
+                                     calc_spectral_mismatch_field,
+                                     martin_ruiz)