Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Update the SDR LS Factor #1289

Merged
merged 11 commits into from
Aug 22, 2023
Merged

Conversation

phargogh
Copy link
Member

@phargogh phargogh commented Apr 28, 2023

This PR implements changes to the SDR model's LS Factor as discussed in #915 and reviewed in an ADR in #1288. In summary, we are:

  1. Replacing the average-aspect (derived from MFD flow direction) approach with on-pixel aspect derived from slope.
  2. Replacing upstream contributing area with an estimate for specific catchment area, $\sqrt{upstream\ contributing\ area}$.

A discussion of the differences between SAGA's LS Factor and InVEST's can be found in #1288.

A branch with changes to the user's guide is available at https://github.com/natcap/invest.users-guide/tree/task/update-ls-factor.

Fixes #915

Checklist

  • Updated HISTORY.rst and link to any relevant issue (if these changes are user-facing)
  • Updated the user's guide (if needed)
  • Tested the affected models' UIs (if relevant)

@phargogh
Copy link
Member Author

As mentioned in #915 , @schmittrjp says we're good to go on including this in the next release of InVEST!

@phargogh phargogh marked this pull request as ready for review August 17, 2023 22:18
@phargogh phargogh requested a review from emlys August 21, 2023 20:54
Copy link
Member

@emlys emlys left a comment

Choose a reason for hiding this comment

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

@phargogh this looks good and it's really helpful to have that ADR for context! Just some minor suggestions.

HISTORY.rst Outdated Show resolved Hide resolved
@@ -1020,23 +1009,19 @@ def _what_drains_to_stream(flow_dir_mfd, dist_to_channel):


def _calculate_ls_factor(
flow_accumulation_path, slope_path, avg_aspect_path, l_max,
flow_accumulation_path, slope_path, l_max,
target_ls_prime_factor_path):
"""Calculate LS factor.

Calculates a modified LS factor as Equation 3 from "Extension and
Copy link
Member

Choose a reason for hiding this comment

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

Is this modified from Equation 3, or is Equation 3 the modified version?

Copy link
Member

Choose a reason for hiding this comment

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

Including the LS factor equation in the docstring could be helpful in the future

Copy link
Member Author

Choose a reason for hiding this comment

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

Yeah, this docstring is a bit of a mess from so many changes, both in development and now with this new update. I have updated the docstring to better reflect the current state, including the (very long) equation, which mirrors what's in the UG. 74877a4


flow_accumulation_info = pygeoprocessing.get_raster_info(
flow_accumulation_path)
flow_accumulation_nodata = flow_accumulation_info['nodata'][0]
cell_size = abs(flow_accumulation_info['pixel_size'][0])
cell_area = cell_size ** 2

def ls_factor_function(
percent_slope, flow_accumulation, avg_aspect, l_max):
def ls_factor_function(percent_slope, flow_accumulation, l_max):
"""Calculate the LS' factor.
Copy link
Member

Choose a reason for hiding this comment

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

Extra '?

Copy link
Member Author

Choose a reason for hiding this comment

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

The prime notation (') was intentional at one point, when this output was called "LS Prime", but now we're just referring to it as the LS factor so I've removed this in 74877a4

# matches the SAGA LS Factor option "square root of catchment area".
# See the InVEST ADR-0001 for more information.
contributing_area = numpy.sqrt(
(flow_accumulation[valid_mask]-1) * cell_area)
Copy link
Member

Choose a reason for hiding this comment

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

I'm not sure I understand the purpose of the -1. Is it that flow_accumulation includes each pixel itself as well as its upstream pixels, and we only want the upstream pixels?

Copy link
Member Author

Choose a reason for hiding this comment

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

Exactly! Flow accumulation includes itself in the count of pixels upstream of this pixel, but that's in conflict with the equation's definition where we want only the pixels that are strictly upstream of this pixel. Added a comment about this in 4b5e13e

'sed_dep': 8.80130577087,
'avoid_exp': 57971.87890625,
'avoid_eros': 1458232.5,
'usle_tot': 2.62457418442,
Copy link
Member

Choose a reason for hiding this comment

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

These are some dramatic changes! I know we didn't have a unit test for _calculate_ls_factor before, but might I suggest adding one? It would be great to confirm that the new results are as expected on some fake data where we can easily spot-check the math.

Copy link
Member Author

Choose a reason for hiding this comment

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

I added a test in aac76d2, though the LS factor function is complicated enough where this is nontrivial to spot check by hand! Even so, I think the test covers our bases w/r/t the various conditions that the function has, so any future changes should be caught by the test.

@phargogh phargogh requested a review from emlys August 22, 2023 22:42
Copy link
Member

@emlys emlys left a comment

Choose a reason for hiding this comment

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

wow, thanks for adding that extensive doctring and test!

@emlys emlys merged commit 1ce08ba into natcap:main Aug 22, 2023
20 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

Reported discrepancy between LS Factor in SDR vs SAGA
2 participants