Diffuse Self-Shading Example: Incorrect transposition of shading.sky_diffuse_passias #2584
Description
Describe the bug
I think the second graph of https://pvlib-python.readthedocs.io/en/stable/gallery/shading/plot_passias_diffuse_shading.html is incorrect, even though it corresponds with the given reference.
In particular, transposing the shading loss in the horizontal plane shading.sky_diffuse_passias(psi) to the plane of array by multiplying it with a transposition ratio seems incorrect to me.
Instead the 'relative diffuse irradiance' (defined as 'ratio of diffuse plane of array irradiance (after accounting for shading) to diffuse horizontal irradiance') should be directly calculated as 1 - shading.sky_diffuse_passias(surface_tilt + psi).
This result also corresponds (more or less) with the result obtained when using the view factor method bifacial.utils.vf_row_sky_2d_integ(surface_tilt, gcr), which I guess is the most accurate method to calculate row-to-row diffuse shading.
The difference between the various methods is illustrated with the following graph:
from pvlib import bifacial, shading, irradiance from cycler import cycler import matplotlib.pyplot as plt import numpy as np surface_tilt = np.arange(0, 90, 0.5) plt.rc('axes', prop_cycle=(cycler('color', ['b', 'orange', 'g', 'r']) * cycler('linestyle', ['-', ':', '--']))) plt.figure() for k in [1, 1.5, 2, 10]: gcr = 1/k psi = shading.masking_angle_passias(surface_tilt, gcr) shading_loss = shading.sky_diffuse_passias(psi) transposition_ratio = irradiance.isotropic(surface_tilt, dhi=1.0) passias_wrong = transposition_ratio * (1-shading_loss) * 100 # % passias_corrected = (1 - shading.sky_diffuse_passias(surface_tilt + psi)) * 100 # % vf_row_sky_2d_integ = bifacial.utils.vf_row_sky_2d_integ(surface_tilt, gcr) * 100 # % plt.plot(surface_tilt, passias_wrong, label=f"passias k={k}") plt.plot(surface_tilt, passias_corrected, label=f"passias corrected k={k}") plt.plot(surface_tilt, vf_row_sky_2d_integ, label=f"vf_row_sky_2d_integ k={k}") plt.xlabel('Inclination angle [degrees]') plt.ylabel('Relative diffuse irradiance [%]') plt.ylim(0, 105) plt.legend() plt.show()