Skip to content
101 changes: 101 additions & 0 deletions pvlib/snow.py
Original file line number Diff line number Diff line change
Expand Up @@ -185,3 +185,104 @@ def dc_loss_nrel(snow_coverage, num_strings):
Available at https://www.nrel.gov/docs/fy18osti/67399.pdf
'''
return np.ceil(snow_coverage * num_strings) / num_strings

def townsend_Se(S, N):
'''
Calculates effective snow for a given month based upon the total snowfall
received in a month in inches and the number of events where snowfall is greater
than 1 inch

Parameters
----------
S : numeric
Snowfall in inches received in a month

N: numeric
Number of snowfall events with snowfall > 1"

Returns
-------
effective_snowfall : numeric
Effective snowfall as defined in the townsend model

References
----------
.. [1] Townsend, Tim & Powers, Loren. (2011). Photovoltaics and snow: An
update from two winters of measurements in the SIERRA. Conference
Record of the IEEE Photovoltaic Specialists Conference.
003231-003236. 10.1109/PVSC.2011.6186627.
Available at https://www.researchgate.net/publication/261042016_Photovoltaics_and_snow_An_update_from_two_winters_of_measurements_in_the_SIERRA

'''
return(np.where(N>0, 0.5 * S * (1 + 1/N), 0))

def townsend_snow_loss_model(S, N, tilt, RH, T_air, POA, R, H, P=40):
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Prefer loss_townsend, and would prefer more descriptive parameter names e.g., snow_total instead of S, num_events instead of N. It's a good idea to maintain traceability to the terms used in the reference but we can do that in the description of each parameter.

'''
Calculates monthly snow loss based on a generalized monthly snow loss model discussed in [1]_.

Parameters
----------
S : numeric
Inches of snow received in the current month

N : numeric
Number of snowfall events with snowfall > 1"

tilt : numeric
Array tilt in degrees

RH : numeric
Relative humidity in percentage

T_air : numeric
Ambient temperature in celcius

POA : numeric
Plane of array irradiance in kWh/m2/month

R : numeric
Row length in the slanted plane of array dimension in inches

H : numeric
Drop height from array edge to ground in inches

P : numeric
piled snow angle, assumed to stabilize at 40° , the midpoint of
25°-55° avalanching slope angles

S_prev : numeric
Inches of snow received in the previous month

N_prev : numeric
Number of 1" or greater snow events in the previous month

Returns
-------
loss : numeric
Average monthly DC capacity loss in percentage due to snow coverage

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Let's add a Notes section here pointing out that this model has not been validated for tracking arrays, but the reference suggests using the maximum rotation angle in place of surface_tilt.

References
----------
.. [1] Townsend, Tim & Powers, Loren. (2011). Photovoltaics and snow: An
update from two winters of measurements in the SIERRA. Conference
Record of the IEEE Photovoltaic Specialists Conference.
003231-003236. 10.1109/PVSC.2011.6186627.
Available at https://www.researchgate.net/publication/261042016_Photovoltaics_and_snow_An_update_from_two_winters_of_measurements_in_the_SIERRA
'''

C1 = 5.7e04
C2 = 0.51

S_prev = np.roll(S,1)
N_prev = np.roll(N,1)

Se = townsend_Se(S, N)
Se_prev = townsend_Se(S_prev, N_prev)

Se_weighted = 1/3 * Se_prev + 2/3 * Se
gamma = (R * Se_weighted * np.cos(np.deg2rad(tilt)))/(np.clip((H**2 - Se_weighted**2),a_min=0.01,a_max=None)/2/np.tan(np.deg2rad(P)))

GIT = 1 - C2 * np.exp(-gamma)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'd also change this to ground_interference_term and use more lines in the equation below.

loss = C1 * Se_weighted * (np.cos(np.deg2rad(tilt)))**2 * GIT * RH / (T_air+273.15)**2 / POA**0.67

return (np.round(loss,2))
23 changes: 23 additions & 0 deletions pvlib/tests/test_snow.py
Original file line number Diff line number Diff line change
Expand Up @@ -95,3 +95,26 @@ def test_dc_loss_nrel():
expected = pd.Series([1, 1, .5, .625, .25, .5, 0])
actual = snow.dc_loss_nrel(snow_coverage, num_strings)
assert_series_equal(expected, actual)

def test_townsend_Se():
S = np.array([10, 10, 5, 1, 0, 0, 0, 0, 0, 0, 5, 10])
N = np.array([2, 2, 1, 0, 0, 0, 0, 0, 0, 0, 2, 3])
expected = np.array([7.5, 7.5, 5, 0, 0, 0, 0, 0, 0, 0, 3.75, 6.66666667])
actual = snow.townsend_Se(S, N)
np.testing.assert_allclose(expected, actual, rtol=1e-07)

def test_townsend_snow_loss_model():
S = np.array([10, 10, 5, 1, 0, 0, 0, 0, 0, 0, 5, 10])
N = np.array([2, 2, 1, 0, 0, 0, 0, 0, 0, 0, 2, 3])
tilt = 20
RH = np.array([80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80])
T_air = np.array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0])
POA = np.array([350, 350, 350, 350, 350, 350, 350, 350, 350, 350, 350, 350])
P = 40
R = 100
H = 10
expected = np.array([7.7, 7.99, 6.22, 1.72, 0, 0, 0, 0, 0, 0, 2.64, 6.07])
actual = snow.townsend_snow_loss_model(S,N,tilt,RH,T_air,POA,R,
H,P)
np.testing.assert_allclose(expected, actual, rtol=1e-07)