mdbartos / pysheds

:earth_americas: Simple and fast watershed delineation in python.
GNU General Public License v3.0
707 stars 191 forks source link

implements transport efficiency in flow accumulation (#70) #73

Closed itati01 closed 5 years ago

itati01 commented 5 years ago

As suggested with separate loops to avoid performance hits.

mdbartos commented 5 years ago

Hi @itati01

Let me know when you're ready to submit and I will review.

Thanks! MDB

itati01 commented 5 years ago

Hi @mdbartos , Sorry, I am still not familiar with git. Are my commits / pull requests not available to you?

mdbartos commented 5 years ago

Hi @itati01

I can see your changes. I just wanted to make sure you had everything ready to go before I reviewed and pulled the changes in (there were a couple commented out lines).

mdbartos commented 5 years ago

A couple things:

1) It looks like you fixed #68 in this PR. Can you go ahead and fix issue #74 and include it in this PR as well? 2) Can you add a test for accumulation with efficiency to tests/test_grid.py? The test that you referenced in #70 should be fine.

itati01 commented 5 years ago

Ok, will fix #74 and have look into test_grid.py.

Feel free to remove any comment. I was unsure about the flatten() in efficiency.flatten() / 100. to convert percentages to relative values. Here I simply followed your steps with the weight raster.

mdbartos commented 5 years ago

Hi @itati01

I believe using flatten to create a flat copy is the correct approach (otherwise you can end up modifying the original array inside the function). For example:

import numpy as np

def test_modify(a):
    a[0] = 0

def test_ravel(a):
    a = a.ravel()
    a[0] = 0

def test_flatten(a):
    a = a.flatten()
    a[0] = 0
a = np.ones(3, dtype=int)
test_modify(a)
print(a)
>>> [0 1 1]

a = np.ones(3, dtype=int)
test_ravel(a)
print(a)
>>> [0 1 1]

a = np.ones(3, dtype=int)
test_flatten(a)
print(a)
>>> [1 1 1]
itati01 commented 5 years ago

Hi @mdbartos,

fixed #74, created a new GeoTIFF, extended test_grid.py to cover a basic testing of the transport efficiency, and minor issues with transport efficiency. So, the PR should be ready.

mdbartos commented 5 years ago

Looks like coveralls is complaining because there's no test for efficiency with dinf accumulation. Might be worth adding to the tests.

itati01 commented 5 years ago

This would be grid.accumulation(data='dinf_dir', dirmap=dirmap, out_name='dinf_acc', routing='dinf', efficiency=grid.eff) in test_accumulation but I have to figure out the correct solution first.

However, there is a problem which should also affect the (still untested) weighting: grid.clip_to('patch') does not apply to the efficiency raster. Shall I create and use an unclipped dinf flow direction to test this function?

mdbartos commented 5 years ago

This would be grid.accumulation(data='dinf_dir', dirmap=dirmap, out_name='dinf_acc', routing='dinf', efficiency=grid.eff) in test_accumulation but I have to figure out the correct solution first.

However, there is a problem which should also affect the (still untested) weighting: grid.clip_to('patch') does not apply to the efficiency raster. Shall I create and use an unclipped dinf flow direction to test this function?

I think you can get around this by using clip_to('dir') to return the view to the original raster dimensions. As a side note, I want to expand the clip_to function so that it takes any Raster as an input (instead of just named datasets). See #65 e.g.

itati01 commented 5 years ago

I addressed the clipping issue by using _input_handler() for weights and transport efficiency in accumulation, similar to flow direction. So, we can new pass their names. Flattening is still required, though. Shall I open an issue?

coveralls commented 5 years ago

Pull Request Test Coverage Report for Build 235


Changes Missing Coverage Covered Lines Changed/Added Lines %
pysheds/grid.py 72 74 97.3%
<!-- Total: 72 74 97.3% -->
Totals Coverage Status
Change from base Build 212: 0.4%
Covered Lines: 1559
Relevant Lines: 1880

💛 - Coveralls
itati01 commented 5 years ago

Have opened an issue regarding dinf accumulation in general (#76). Nonetheless, added a simple test case for the efficiency.

mdbartos commented 5 years ago

Hi @itati01

Thanks for all your hard work! I think it's almost ready to pull in.

I addressed the clipping issue by using _input_handler() for weights and transport efficiency in accumulation, similar to flow direction. So, we can new pass their names. Flattening is still required, though. Shall I open an issue?

I like the idea of using _input_handler() for the weights/efficiency. It would make things more consistent overall. There are a couple issues though:

1) Wouldn't passing in an ndarray for weights now raise an error? That would be a breaking change. 2) This same behavior would be needed for the other methods that use weights (such as flow_distance).

It might be best to save this change for later, at least until these two things are resolved.

I think you can get the behavior you want by using grid.view to return a view of the efficiency raster at the grid's current dimensions. So, something like:

grid.clip_to('dir')
eff = grid.view('eff')
grid.accumulation(data='dir', dirmap=dirmap, out_name='acc_eff', efficiency=eff)

See here for more info: https://mdbartos.github.io/pysheds/views.html

itati01 commented 5 years ago

Works, thanks. Why is this not working with grid.view()

grid.eff[grid.eff==grid.eff.nodata] = 1 # change nodata to 1 grid.clip_to("catch") assert(grid.view("eff").min() == grid.eff.min())

but

grid.clip_to("dir") assert(grid.view("eff").min() == grid.eff.min())

mdbartos commented 5 years ago

Works, thanks. Why is this not working with grid.view()

grid.eff[grid.eff==grid.eff.nodata] = 1 # change nodata to 1 grid.clip_to("catch") assert(grid.view("eff").min() == grid.eff.min())

but

grid.clip_to("dir") assert(grid.view("eff").min() == grid.eff.min())

I will take a look later today and let you know.