DISC model and zenith threshold

[anomam]

Hi,

I tried to look for previous issues that address this question ( #24 and #6 ) but I don't think they did.
I noticed that a threshold was used on zenith values in the DISC model above which the calculated DNI was set to 0.
https://github.com/pvlib/pvlib-python/blob/master/pvlib/irradiance.py#L1158

 dni = np.where((zenith > 87) | (ghi < 0) | (dni < 0), 0, dni)

Does the condition (zenith > 87) come from the original model or was it implemented to filter potential outliers? In which case why 87? And can we create a function argument that would allow us to use different threshold values?

Thanks!

[anomam]

PS: I skimmed through the first reference Maxwell [1] and didn't find the condition there, but maybe I went too fast.

[wholmgren]

The filter is (mostly) inherited from the matlab code, which apparently chose 87 degrees to match a FORTRAN implementation of the algorithm:

https://github.com/sandialabs/MATLAB_PV_LIB/blob/master/pvl_disc.m#L134

DNI(Z>87 | GHI<1 | DNI<0) = 0; % 87 degrees chosen to simulate the FORTRAN code in use by SRRL (from Perez)

I say "mostly" because we changed the ghi filter condition from 1 to 0 at some point, which was probably a mistake. The filter should be ghi|dni <= 0, which for floating point numbers, means ghi|dni < 0 + a small number. Maybe 1 is too big, but it shouldn't be 0.

I'm guessing that we can make a more permissive filter without needing to add an optional argument.

Is the zenith filter actually necessary if we already have the ghi and dni filters? It's not clear to me how the algorithm would behave without it. Maybe filtering on some other quantity (e.g. dni = 0 if Kn > 1?) would be preferable.

The zenith filter is probably preventing nans in the airmass calculation, but you could handle that differently.

[wholmgren]

This is also related to #250

[cwhanse]

The 87 threshold is in the Matlab code to mimic FORTRAN which we got from SRRL. When we started PVLib we tried to implement models faithfully with the published references. In this case, because the publication isn't clear on a few points, we followed code which appeared to be authoritative.

From the comments in the code, it appears at one point in the code development, true zenith > 87 was causing numerical issues. I'm not sure that's still a problem. In the matlab code, one does have issues if zenith > 93.885, because the air mass function will return imaginary numbers, which casts the valid values as also imaginary.

I note that NREL distributes an Excel sheet implementing the DISC model which sets the zenith limit to 80. I think we can change the value if we choose.

We could add a parameter zenith_limit to the function signature with a default of 87.

[adriesse]

Coming at the end of the function, this seems to be an optional step. It would be great if someone could look at what the algorithm produces in this region, to determine whether and how the output becomes unreasonable. I would prefer (again) to see a gradual transition from the model to the constraint--if necessary.

[anomam]

Thank you all for your replies. I just opened a PR (#318) to implement @cwhanse 's recommendation to use an optional argument in the function signature and put 87 as the default value.
I hope this solution will satisfy everybody!

[wholmgren]

This should have been closed by #400