Optimizing raster.icrop
#660
Conversation
|
Are we sure that the DEM is not completely load before crop ? |
| @@ -173,7 +173,7 @@ def _outside_image( | |||
| """See description of Raster.outside_image.""" | |||
|
|
|||
| if not index: | |||
| xi, xj = _xy2ij(xi, yj, transform=transform, area_or_point=area_or_point) | |||
| xi, yj = _xy2ij(xi, yj, transform=transform, area_or_point=area_or_point) | |||
| @@ -2464,15 +2480,29 @@ def icrop( | |||
| if len(bbox) != 4: | |||
| raise ValueError("bbox must be a list or tuple of four integers: (xmin, ymin, xmax, ymax).") | |||
|
|
|||
| # Convert pixel coordinates to georeferenced coordinates using the transform | |||
There was a problem hiding this comment.
The argument description above was ambiguous for me:
:param bbox: Pixel-based bounding box as (xmin, ymin, xmax, ymax) in pixel coordinates.
Maybe similar language as in the function description "Bounding box based on indices of the raster array".
We could also add a note in the general description to clarify why this is strikingly different, something like: "Note that, by convention, a raster has the Y axis as first in its array, and that this Y axis is also flipped (first index corresponds to the highest-value coordinate)."
| # Crop the raster data based on the given bounding box (bbox), limit the range to the borders of the raster | ||
| crop_img = self.data[..., bbox[0] : min(self.height, bbox[2]), bbox[1] : min(bbox[3], self.width)] | ||
|
|
||
| # Create a new affine transform that corresponds to the cropped area |
There was a problem hiding this comment.
Took me a little while to understand why this transform looked like a translation!
When we crop within the extent, the new array size determines the new height/width, so it is only a translation for the affine transform!
In this case, let's use our tested _translate function to make it easier to read/debug for future us 😉:
| @@ -1818,28 +1818,33 @@ def test_icrop(self, example: str, bbox: list[int] | tuple[int, ...]) -> None: | |||
| raster = gu.Raster(example) | |||
|
|
|||
| # Call the icrop method with the given bounding box, returning a cropped raster | |||
| raster_cropped = raster.icrop(bbox=bbox) | |||
| raster_icropped = raster.icrop(bbox=bbox) | |||
|
|
|||
| # Compute the expected affine transform for the cropped raster | |||
| transformed_cropped = rio.transform.Affine( | |||
| bottom_left = self.transform * (bbox[0], self.height - bbox[1]) | ||
| top_right = self.transform * (bbox[2], self.height - bbox[3]) | ||
| # Crop the raster data based on the given bounding box (bbox), limit the range to the borders of the raster | ||
| crop_img = self.data[..., bbox[0] : min(self.height, bbox[2]), bbox[1] : min(bbox[3], self.width)] |
There was a problem hiding this comment.
Ah sorry I didn't notice this line and @adebardo's comment, good catch.
Indeed, this will load the entire array, then subset it. To circumvent this, we need to use crop() as before.
There was a problem hiding this comment.
We should add a test checking if Raster.is_loaded is still False after calling icrop().
Resolves #659
Description
icropbounding box to match with bouding box of other objects.cropto avoid useless coordinates conversion.Tests
icroptests.