search

GlacioHack / xdem

Analysis of digital elevation models (DEMs)
https://xdem.readthedocs.io/
MIT License
144 stars 41 forks source link

Nitpick; unexpected behavior in resolution term of `xdem.terrain.get_quadric_coefficients()` #589

Open erikmannerfelt opened 1 month ago

erikmannerfelt commented 1 month ago

I'm using the core part of the terrain module for a machine learning exercise, and found the xdem.terrain.get_quadric_coefficients() to be slightly unexpected in how it treated the resolution term. At first, I added the DEM class as an argument, which works, but a resolution term still has to be added (opposed to other tools of ours). So I added dem.res, and got an (at first) strange exception:

ValueError: Resolution must be the same for X and Y directions

It is weird because I expected (and validated) that the X/Y-direction res was the same. Turns out, it's just from how the resolution term is validated! It checks if the given resolution is a Sized object (tuple, list, ndarray), and then flips out if true.

https://github.com/GlacioHack/xdem/blob/34e5049ae76a858b9c95c8763d9763b298c56e3c/xdem/terrain.py#L330-L332

Solutions to this would be:

  1. First check if it's Sized, and then check if resolution[0] != resolution[1].
  2. Clarify the error message to show what's actually wrong (more like "resolution needs to be one value, not a list/tuple")
  3. Allow a Raster (sub-)class and automatically parse the resolution from it.

I'd vouch for either 1), 3), or both.

rhugonnet commented 1 month ago

Since recently we can also do DEM.slope(), DEM.curvature() or more generically DEM.get_terrain_attribute() which adequately handles the resolution. It's actually going to be the preferred way showed in the documentation once it's up, to avoid passing volative metadata. And the reasoning is that it'll be the same with the Xarray accessor ds.dem.slope(), ds.dem.curvature(), etc...

So another question is should we allow a DEM as input of terrain function when it is already possible to call the method from a DEM? Right now I left the terrain module as public so that a user that wants to just pass a NumPy array + a resolution still can without having to rely on our own objects, but we could also discuss if this is really necessary.

But to go back to the issue at hand: in any case I think we should do 1/ you propose. And additionally, if we decide to leave the terrain functions public for NumPy array support, then we could add 3/ as well.

erikmannerfelt commented 1 month ago

Since recently we can also do DEM.slope(), DEM.curvature() or more generically DEM.get_terrain_attribute() which adequately handles the resolution. It's actually going to be the preferred way showed in the documentation once it's up, to avoid passing volative metadata. And the reasoning is that it'll be the same with the Xarray accessor ds.dem.slope(), ds.dem.curvature(), etc...

So another question is should we allow a DEM as input of terrain function when it is already possible to call the method from a DEM? Right now I left the terrain module as public so that a user that wants to just pass a NumPy array + a resolution still can without having to rely on our own objects, but we could also discuss if this is really necessary.

But to go back to the issue at hand: in any case I think we should do 1/ you propose. And additionally, if we decide to leave the terrain functions public for NumPy array support, then we could add 3/ as well.

Thanks for the context! My main concern was more the apparent inconsistency in how it works compared to other functions rather than the functionality in general!

I like the idea of keeping "numpy support" so I'd vouch for implementing 3., as it shouldn't be that much work!

P.S. It's super nice that DEM.slope() etc. works now!!