Calculating horizon profiles and associated shading losses

[JPalakapillyKWH]

• [x] I am familiar with the contributing guidelines.
• [x] Fully tested. Added and/or modified tests to ensure correct behavior for all reasonable inputs. Tests (usually) must pass on the TravisCI and Appveyor testing services.
• [ ] Updates entries to docs/sphinx/source/api.rst for API changes.
• [ ] Adds description and name entries in the appropriate docs/sphinx/source/whatsnew file for all changes.
• [ ] Code quality and style is sufficient. Passes LGTM and SticklerCI checks.
• [x] New code is fully documented. Includes sphinx/numpydoc compliant docstrings and comments in the code where necessary.
• [x] Pull request is nearly complete and ready for detailed review.

Code to get generate horizon profiles using the googlemaps API. Also has functionality to create shading losses due to the horizon profile. Currently only done with an isotropic sky diffuse model and the direct normal. Follows much of the methodology in https://doi.org/10.1016/j.solener.2014.09.037

I've been developing this code in recent weeks. Is this something PVlib wants to incorporate?

[adriesse]

I think this would make a great addition!

[cwhanse]

Is this something PVlib wants to incorporate?

Yes, I think so. horizon.py is OK, I would like to think a bit more about how this module relates to existing modules (irradiance.py, bifacial.py) and other process models not yet in pvlib, e.g., effects of partial shading on an array.

[cwhanse]

@JPalakapillyKWH I'm not familiar with Google maps elevation data. What is the spatial resolution? Do the data represent the near-field or only the far-field horizon? Near-field accounts for trees, buildings, etc., far-field is terrain such as mountains and hills.

[JPalakapillyKWH]

@cwhanse The spatial resolution varies since googlemaps (as I've gathered) uses a patchwork of different signals but is generally at least 90m. The code will capture the near-field horizon to the extent that the googlemaps AP is affected by buildings/trees etc.

I included the googlemaps elevation calls because it was a simple way to get elevation data for sites distributed across the world, but the bulk of the horizon module is agnostic to the source of the elevation data.

[cwhanse]

Ok, was just thinking how to describe the output of this module, as near-field or far-field shading. Horizon is usually interpreted as far field.

[JPalakapillyKWH]

It is designed for far-field.

[cwhanse]

@JPalakapillyKWH super content for pvlib. I've read the reference and look through the code. Some high-level comments about modules and connecting with the rest of pvlib. would appreciate @wholmgren @mikofski input here.

• horizon.py main purpose is to provide functions that return the horizon profile (list of (azimuth, dip angle) tuples) and sky diffuse loss factor dtf. I'd put calculate_dtf and collection_plane_dip_angle in horizon.py. I'd add the Boolean "in-front" variable as an output if that is being calculated somewhere.
• Some of the helper functions e.g. pol2cart, latitude_to_geocentric could go into pvlib.tools. I'd leave the grid, sample, filter and gmaps functions in horizon.py.
• There's a mix of "to" and "from" function names, I'd prefer "to".
• I am wondering if there's a way to use a modified calculation of dtf in the current isotropic sky diffuse model, rather than making a new model horizon_adjusted. Could we change the denominator in the dtf calculation to sum over only the cells in front of the array? Then sky diffuse = dhi * (1 + cos(tilt))/2 * dtf and dtf could be an optional kwarg for the current isotropic function.
• it's too bad the reference isn't more explicit about the time series simulation.

[wholmgren]

I skimmed the code and it looks like a great start - thanks @JPalakapillyKWH! I've not yet looked at the reference.

In addition to expanding some of the function documentation, I think we're going to want some documentation that synthesizes the proposed code.

There are a handful of places that generate tuples or loop over them in calculations. Can the code make more use of numpy arrays to improve performance and readability?

I would encourage marking some of the helper functions as private with an _ and/or putting some of these functions in pvlib.tools (which is not explicitly private but is also not documented in the api.rst file).

Great that the code is agnostic about elevation data source. I am not so sure about including the googlemaps api wrapper code in pvlib python, but I'm definitely ok using it in example documentation. I'd probably want to use the functions with a gridded data source such as SRTM.

[JPalakapillyKWH]

@cwhanse and @wholmgren I made a lot of the code restructuring changes y'all suggested. Many functions were moved to tools.py.

I am currently working on making the code more numpy-friendly. I'm also working on creating documentation by creating a .rst file in the sphinx source folder and updating api.rst. Is this right way to create documentation?

I did include a few functions for visualizing horizon profiles but they require matplotlib to run. Do these seem worth including in the horizon module itself or is this something that would go in the documentation?

@cwhanse I don't think summing over the front of the panel will work since the horizon could be at a higher elevation angle than the back of the panel (e.g. panel tilt is 20 deg but the horizon is at a 30 deg elevation angle behind it). I just included the horizon adjustment code in the isotropic model with horizon profile as a kwarg. Let me know what you think.

[mikofski]

Sorry to jump in late. Great addition @JPalakapillyKWH.

One caveat that should probably be mentioned somewhere, is that I believe some irradiance sources, such as NSRDB may already include the far field horizon loss, since they are calibrated by ground stations, please correct me if I'm wrong.

I too would like to see examples from other geographical sources like NASA SRTM and ESRI, but I think those should be separate PRs .

I'll try to take a closer look later today.

Thanks!

[cwhanse]

@cwhanse I don't think summing over the front of the panel will work

I was not clear - I meant sum over the sky cells visible from the plane of array. The revised implementation accomplishes the same end for the isotropic sky and is arguably closer to the published reference.

Not in this PR, but it appears to me that we can extend the use of horizon data to other transposition models. That's what I had in mind when I made the comment.

[JPalakapillyKWH]

@mikofski You are correct about ground station data already incorporating far field losses. This is designed for satellite irradiance data (PSM, TWC etc.) @adriesse I took many of your suggestions and ran with them.

There's been a lot of restructuring in the code to move towards numpy. I'm still updating docstrings though.

[mikofski]

ground station data already incorporating far field losses. This is designed for satellite irradiance data (PSM, TWC etc.)

Awesome! Would you consider adding a advisory note in the docstring with this warning on usage?

[JPalakapillyKWH]

Awesome! Would you consider adding a advisory note in the docstring with this warning on usage?

Which docstring do you mean?

[mikofski]

I would recommend adding a warning note in any method/function that uses the new horizon_* args/kwargs:

• pvlib.irradiance.isotropic
• pvlib.irradiance.get_sky_diffuse
• pvlib.irradiance.get_total_irradiance

But this is not a blocker for me. It's just a suggestion, my opinion, take it or leave it :)

[JPalakapillyKWH]

There seem to be a lot of timeout issues in the tests. Anybody have any insight? I pass pytest pvlib locally.

[adriesse]

Is it necessary/desirable to pull these options into the existing irradiance functions? For my purposes I would be happy with two separate functions:

• one to give me an adjustment factor for diffuse irradiance, and
• one to give me an adjustment factor for DNI (binary in its simplest form).

Anyway even if the functionality is integrated somehow, it doesn't all have to be done in one PR.

[JPalakapillyKWH]

@adriesse Fair point. I can include just horizon.py and test_horizon.py in this PR.

[adriesse]

@adriesse Fair point. I can include just horizon.py and test_horizon.py in this PR.

@JPalakapillyKWH Beware though, that my opinion is only one; there may be others...

[JPalakapillyKWH]

@JPalakapillyKWH Beware though, that my opinion is only one; there may be others...

I tend to agree that this code should be split into multiple PRs. Hopefully this will be easier to review now.

[JPalakapillyKWH]

Working my way through the docstrings and code. I am a little confused about the purpose of the "sample" functions - if a grid is being supplied as input, why does the grid need to be sampled? Is the input grid not regularly spaced?

The input grid will be regularly spaced. The sampling functions are to produce a smooth horizon profile. Without sampling, the "smoothness" of the horizon profile is highly dependent on the size of your grid. If you don't have enough points on your grid then you will have a very choppy horizon profile. Some azimuths might not even have elevation angles. So the sampling functions are a way of pseudo-generating more elevation data to make the horizon smoother.

[wholmgren]

I stumbled across earthpy and this hillshade example. Thought people interested in this PR might want to check it out.

[cwhanse]

Working my way through the docstrings and code. I am a little confused about the purpose of the "sample" functions - if a grid is being supplied as input, why does the grid need to be sampled? Is the input grid not regularly spaced?

The input grid will be regularly spaced. The sampling functions are to produce a smooth horizon profile. Without sampling, the "smoothness" of the horizon profile is highly dependent on the size of your grid. If you don't have enough points on your grid then you will have a very choppy horizon profile. Some azimuths might not even have elevation angles. So the sampling functions are a way of pseudo-generating more elevation data to make the horizon smoother.

Makes sense now, thanks.

[adriesse]

I read the paper mentioned at the top of this PR and my impression is that there is not much overlap between the methods discussed there, and those programmed here. Is this correct?

I would very much like to see notes in the docstrings regarding the nature of the algorithms used, their origins, and their validation. If there are new original methods, please also identify and take credit for them! If there is a forthcoming publication, I would love to see it; but if you're still working on one, I would suggest pausing the PR until it is finished and you can share it.

[JPalakapillyKWH]

@adriesse You are right in that there is not a whole lot of overlap. Much of the code in this PR is taken from the geometry used in viewshed analysis. The sampling functions and binning of the horizon profile were ideas that I got from http://sherrytowers.com/2014/04/13/archeoastronomy-calculating-the-horizon-profile-using-online-us-geographic-survey-data/. I don't have any plans for a publication, just working on this for the summer.

For my own validation of horizon profile modeling, I used the algorithms in this module with elevation data taken from the googlemaps elevation API to calculate horizon profiles. These were mostly compared by eye to horizon profiles from other sources and the online tool http://www.heywhatsthat.com/.

I also did validation of shading losses due to the horizon. I found shading loss by comparing the output of a modelchain run with the horizon adjusted irradiance vs without. I then compared this to far field losses outputted by PVSyst for a few modeled projects and found that they agreed.

[adriesse]

Perhaps a good next step would be to add (or expand) Notes sections for the doc strings, where you explain the algorithms used, and give references to the sources you used, where applicable.

If validation is not provided in some formal publication, then perhaps it can be done in a jupyter notebook or something, but I don't know what the best place would be to store such a thing. Of course you could also reconsider whether maybe you could make a publication out of this work...

[cwhanse]

Perhaps a good next step would be to add (or expand) Notes sections for the doc strings, where you explain the algorithms used, and give references to the sources you used, where applicable.

calculate_dtf can reference the Goss et al. Solar Energy paper for the idea of dtf. We need a reference for the summarization formula (lines 598 - 601). I thought the code was implementing the sum over equal area sky patches? I may have missed a decision along the way.

I'd reference the ACM transactions paper instead of the class notes for uniformly_sample_triangle.

Add the reference to elevation_and_azimuth - it is called out as [1] but not listed. The rest of the functions seem straight-forward enough.

If validation is not provided in some formal publication, then perhaps it can be done in a jupyter notebook or something, but I don't know what the best place would be to store such a thing. Of course you could also reconsider whether maybe you could make a publication out of this work...

I vote for a page on readthedocs. Separate PR though.

[cwhanse]

@JPalakapillyKWH could you contact me privately at [email protected]? I have a few a questions I'd appreciate your input.

[mikofski]

I would like to see this topic picked up again, I added an issue #1290 to spur more interest, that issue also adds the need to download the horizon data (or figure out how to generate it form GIS elevation data). I couldn't tell if this PR already encompasses downloading horizon data or just applying it to existing solar resource.

[cwhanse]

Most of this PR involved turning a [lat, long, elev] grid into [azimuth, elevation] for an observer location. I don't see any code to download the topographic data, but I could have missed it.

The original submitter @JPalakapillyKWH isn't available to complete this PR, so feel free (anyone) to pick it up from here.