|
| 1 | +""" |
| 2 | +Calculating daily diffuse PAR using Spitter's relationship |
| 3 | +========================================================== |
| 4 | +
|
| 5 | +This example demonstrates how to calculate the diffuse photosynthetically |
| 6 | +active radiation (PAR) from diffuse fraction of broadband insolation. |
| 7 | +""" |
| 8 | + |
| 9 | +# %% |
| 10 | +# The photosynthetically active radiation (PAR) is a key metric in quantifying |
| 11 | +# the photosynthesis process of plants. As with broadband irradiance, PAR can |
| 12 | +# be divided into direct and diffuse components. The diffuse fraction of PAR |
| 13 | +# with respect to the total PAR is important in agrivoltaic systems, where |
| 14 | +# crops are grown under solar panels. The diffuse fraction of PAR can be |
| 15 | +# calculated using the Spitter's relationship [1]_ implemented in |
| 16 | +# :py:func:`~pvlib.irradiance.diffuse_par_spitters`. |
| 17 | +# This model requires the average daily solar zenith angle and the |
| 18 | +# daily fraction of the broadband insolation that is diffuse as inputs. |
| 19 | +# |
| 20 | +# .. note:: |
| 21 | +# Understanding the distinction between the broadband insolation and the PAR |
| 22 | +# is a key concept. Broadband insolation is the total amount of solar |
| 23 | +# energy that gets to a surface, often used in PV applications, while PAR |
| 24 | +# is a measurement of a narrower spectrum of wavelengths that are involved |
| 25 | +# in photosynthesis. See section on *Photosynthetically Active insolation* |
| 26 | +# in pp. 222-223 of [1]_. |
| 27 | +# |
| 28 | +# References |
| 29 | +# ---------- |
| 30 | +# .. [1] C. J. T. Spitters, H. A. J. M. Toussaint, and J. Goudriaan, |
| 31 | +# 'Separating the diffuse and direct component of global radiation and its |
| 32 | +# implications for modeling canopy photosynthesis Part I. Components of |
| 33 | +# incoming radiation', Agricultural and Forest Meteorology, vol. 38, |
| 34 | +# no. 1, pp. 217-229, Oct. 1986, :doi:`10.1016/0168-1923(86)90060-2`. |
| 35 | +# |
| 36 | +# Read some example data |
| 37 | +# ^^^^^^^^^^^^^^^^^^^^^^ |
| 38 | +# Let's read some weather data from a TMY3 file and calculate the solar |
| 39 | +# position. |
| 40 | + |
| 41 | +import pvlib |
| 42 | +import pandas as pd |
| 43 | +import matplotlib.pyplot as plt |
| 44 | +from matplotlib.dates import AutoDateLocator, ConciseDateFormatter |
| 45 | +from pathlib import Path |
| 46 | + |
| 47 | +# Datafile found in the pvlib distribution |
| 48 | +DATA_FILE = Path(pvlib.__path__[0]).joinpath("data", "723170TYA.CSV") |
| 49 | + |
| 50 | +tmy, metadata = pvlib.iotools.read_tmy3( |
| 51 | + DATA_FILE, coerce_year=2002, map_variables=True |
| 52 | +) |
| 53 | +tmy = tmy.filter( |
| 54 | + ["ghi", "dhi", "dni", "pressure", "temp_air"] |
| 55 | +) # remaining columns are not needed |
| 56 | +tmy = tmy["2002-09-06":"2002-09-21"] # select some days |
| 57 | + |
| 58 | +solar_position = pvlib.solarposition.get_solarposition( |
| 59 | + # TMY timestamp is at end of hour, so shift to center of interval |
| 60 | + tmy.index.shift(freq="-30T"), |
| 61 | + latitude=metadata["latitude"], |
| 62 | + longitude=metadata["longitude"], |
| 63 | + altitude=metadata["altitude"], |
| 64 | + pressure=tmy["pressure"] * 100, # convert from millibar to Pa |
| 65 | + temperature=tmy["temp_air"], |
| 66 | +) |
| 67 | +solar_position.index = tmy.index # reset index to end of the hour |
| 68 | + |
| 69 | +# %% |
| 70 | +# Calculate daily values |
| 71 | +# ^^^^^^^^^^^^^^^^^^^^^^ |
| 72 | +# The daily average solar zenith angle and the daily diffuse fraction of |
| 73 | +# broadband insolation are calculated as follows: |
| 74 | + |
| 75 | +daily_solar_zenith = solar_position["zenith"].resample("D").mean() |
| 76 | +# integration over the day with a time step of 1 hour |
| 77 | +daily_tmy = tmy[["ghi", "dhi"]].resample("D").sum() * 1 |
| 78 | +daily_tmy["diffuse_fraction"] = daily_tmy["dhi"] / daily_tmy["ghi"] |
| 79 | + |
| 80 | +# %% |
| 81 | +# Calculate Photosynthetically Active Radiation |
| 82 | +# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ |
| 83 | +# The total PAR can be approximated as 0.50 times the broadband horizontal |
| 84 | +# insolation (integral of GHI) for an average solar elevation higher that 10°. |
| 85 | +# See section on *Photosynthetically Active Radiation* in pp. 222-223 of [1]_. |
| 86 | + |
| 87 | +par = pd.DataFrame({"total": 0.50 * daily_tmy["ghi"]}, index=daily_tmy.index) |
| 88 | +if daily_solar_zenith.min() < 10: |
| 89 | + raise ValueError( |
| 90 | + "The total PAR can't be assumed to be half the broadband insolation " |
| 91 | + + "for average zenith angles lower than 10°." |
| 92 | + ) |
| 93 | + |
| 94 | +# Calculate broadband insolation diffuse fraction, input of the Spitter's model |
| 95 | +daily_tmy["diffuse_fraction"] = daily_tmy["dhi"] / daily_tmy["ghi"] |
| 96 | + |
| 97 | +# Calculate diffuse PAR fraction using Spitter's relationship |
| 98 | +par["diffuse_fraction"] = pvlib.irradiance.diffuse_par_spitters( |
| 99 | + solar_position["zenith"], daily_tmy["diffuse_fraction"] |
| 100 | +) |
| 101 | + |
| 102 | +# Finally, calculate the diffuse PAR |
| 103 | +par["diffuse"] = par["total"] * par["diffuse_fraction"] |
| 104 | + |
| 105 | +# %% |
| 106 | +# Plot the results |
| 107 | +# ^^^^^^^^^^^^^^^^ |
| 108 | +# Insolation on left axis, diffuse fraction on right axis |
| 109 | + |
| 110 | +fig, ax_l = plt.subplots(figsize=(12, 6)) |
| 111 | +ax_l.set( |
| 112 | + xlabel="Time", |
| 113 | + ylabel="Daily insolation $[Wh/m^2/day]$", |
| 114 | + title="Diffuse PAR using Spitter's relationship", |
| 115 | +) |
| 116 | +ax_l.xaxis.set_major_formatter( |
| 117 | + ConciseDateFormatter(AutoDateLocator(), tz=daily_tmy.index.tz) |
| 118 | +) |
| 119 | +ax_l.plot( |
| 120 | + daily_tmy.index, |
| 121 | + daily_tmy["ghi"], |
| 122 | + label="Broadband: total", |
| 123 | + color="deepskyblue", |
| 124 | +) |
| 125 | +ax_l.plot( |
| 126 | + daily_tmy.index, |
| 127 | + daily_tmy["dhi"], |
| 128 | + label="Broadband: diffuse", |
| 129 | + color="skyblue", |
| 130 | + linestyle="-.", |
| 131 | +) |
| 132 | +ax_l.plot(daily_tmy.index, par["total"], label="PAR: total", color="orangered") |
| 133 | +ax_l.plot( |
| 134 | + daily_tmy.index, |
| 135 | + par["diffuse"], |
| 136 | + label="PAR: diffuse", |
| 137 | + color="coral", |
| 138 | + linestyle="-.", |
| 139 | +) |
| 140 | +ax_l.grid() |
| 141 | +ax_l.legend(loc="upper left") |
| 142 | + |
| 143 | +ax_r = ax_l.twinx() |
| 144 | +ax_r.set(ylabel="Diffuse fraction") |
| 145 | +ax_r.plot( |
| 146 | + daily_tmy.index, |
| 147 | + daily_tmy["diffuse_fraction"], |
| 148 | + label="Broadband diffuse fraction", |
| 149 | + color="plum", |
| 150 | + linestyle=":", |
| 151 | +) |
| 152 | +ax_r.plot( |
| 153 | + daily_tmy.index, |
| 154 | + par["diffuse_fraction"], |
| 155 | + label="PAR diffuse fraction", |
| 156 | + color="chocolate", |
| 157 | + linestyle=":", |
| 158 | +) |
| 159 | +ax_r.legend(loc="upper right") |
| 160 | + |
| 161 | +plt.show() |
0 commit comments