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.

Sign up for GitHub

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

Modular coreg.py revision #71

Merged
merged 17 commits into from
Apr 8, 2021
Merged

Modular coreg.py revision #71

Changes from 1 commit
Commits
Show all changes
17 commits
Select commit Hold shift + click to select a range
36adb81
Added Coreg base class and BiasCorr subclass
erikmannerfelt Mar 26, 2021
921c675
Incremental commit on new Coreg class
erikmannerfelt Mar 27, 2021
7845b56
Added NuthKaab, Deramp, and CoregPipeline classes. Incremental progre…
erikmannerfelt Apr 1, 2021
4841d4a
Added + operator functionality for Coreg classes
erikmannerfelt Apr 5, 2021
038678e
Added CoregPipeline + operator functionality. Improved/added docstrings
erikmannerfelt Apr 5, 2021
60fc38d
Added documentation for new class
erikmannerfelt Apr 5, 2021
7603b8d
Added deprecation warnings for old coregistration functions.
erikmannerfelt Apr 5, 2021
9038292
Fixed bug in tests
erikmannerfelt Apr 5, 2021
987c306
Removed code block from titles
erikmannerfelt Apr 5, 2021
625798b
Small improvement. Dummy commit to trigger new readthedocs version
erikmannerfelt Apr 5, 2021
593ef00
Added subsampling functionality for coreg fitting.
erikmannerfelt Apr 7, 2021
b569902
Improved test documentation
erikmannerfelt Apr 7, 2021
843b7b0
Changed 'mask' argument in Coreg.fit() to 'inlier_mask' and improved …
erikmannerfelt Apr 7, 2021
2ea683b
Changed _trasnform_to_bounds_and_res to take a shape instead of an array
erikmannerfelt Apr 8, 2021
bb4325c
Changed post-ICP resampling from nearest to cubic. Much slower but mu…
erikmannerfelt Apr 8, 2021
56a4d1b
Cleaned code
erikmannerfelt Apr 8, 2021
1f3e735
Added 'verbose' argument and associated callbacks to every coreg.fit(…
erikmannerfelt Apr 8, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
21 changes: 11 additions & 10 deletions xdem/coreg.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1083,11 +1083,11 @@ def coregister(reference_raster: Union[str, gu.georaster.Raster], to_be_aligned_
return aligned_raster, error


def _transform_to_bounds_and_res(dem: np.ndarray,
def _transform_to_bounds_and_res(shape: tuple[int, int],
adehecq marked this conversation as resolved.
Show resolved Hide resolved
transform: rio.transform.Affine) -> tuple[rio.coords.BoundingBox, float]:
bounds = rio.coords.BoundingBox(
*rio.transform.array_bounds(dem.shape[0], dem.shape[1], transform=transform))
resolution = (bounds.right - bounds.left) / dem.shape[1]
*rio.transform.array_bounds(shape[0], shape[1], transform=transform))
resolution = (bounds.right - bounds.left) / shape[1]

return bounds, resolution

Expand Down Expand Up @@ -1275,6 +1275,7 @@ class ICP(Coreg):
Estimates a rigid transform (rotation + translation) between two DEMs.

Requires 'opencv'
See opencv docs for more info: https://docs.opencv.org/master/dc/d9b/classcv_1_1ppf__match__3d_1_1ICP.html
"""

def __init__(self, max_iterations=100, tolerance=0.05, rejection_scale=2.5, num_levels=6): # pylint: disable=super-init-not-called
Expand All @@ -1300,7 +1301,7 @@ def _fit_func(self, ref_dem: np.ndarray, tba_dem: np.ndarray, transform: Optiona
if weights is not None:
warnings.warn("ICP was given weights, but does not support it.")

bounds, resolution = _transform_to_bounds_and_res(ref_dem, transform)
bounds, resolution = _transform_to_bounds_and_res(ref_dem.shape, transform)
points: dict[str, np.ndarray] = {}
# Generate the x and y coordinates for the reference_dem
x_coords, y_coords = np.meshgrid(
Copy link
Member

Choose a reason for hiding this comment

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

See my comment about subsampling. It could be done at this stage for example, to reduce memory usage.

Expand Down Expand Up @@ -1340,7 +1341,7 @@ def _fit_func(self, ref_dem: np.ndarray, tba_dem: np.ndarray, transform: Optiona

def _apply_func(self, dem: np.ndarray, transform: rio.transform.Affine) -> np.ndarray:
"""Apply the coregistration matrix to a DEM."""
bounds, resolution = _transform_to_bounds_and_res(dem, transform)
bounds, resolution = _transform_to_bounds_and_res(dem.shape, transform)
x_coords, y_coords = np.meshgrid(
np.linspace(bounds.left + resolution / 2, bounds.right - resolution / 2, num=dem.shape[1]),
np.linspace(bounds.bottom + resolution / 2, bounds.top - resolution / 2, num=dem.shape[0])[::-1]
Expand All @@ -1362,7 +1363,7 @@ def _apply_func(self, dem: np.ndarray, transform: rio.transform.Affine) -> np.nd
np.linspace(bounds.left, bounds.right, dem.shape[1]) - bounds.left,
np.linspace(bounds.bottom, bounds.top, dem.shape[0])[::-1] - bounds.bottom
)),
method="nearest"
method="linear"
)
aligned_dem[~valid_mask] = np.nan

Expand Down Expand Up @@ -1398,7 +1399,7 @@ def __init__(self, degree: int = 1):
def _fit_func(self, ref_dem: np.ndarray, tba_dem: np.ndarray, transform: Optional[rio.transform.Affine],
weights: Optional[np.ndarray]):
"""Fit the dDEM between the DEMs to a least squares polynomial equation."""
bounds, resolution = _transform_to_bounds_and_res(ref_dem, transform)
bounds, resolution = _transform_to_bounds_and_res(ref_dem.shape, transform)
x_coords, y_coords = np.meshgrid(
np.linspace(bounds.left + resolution / 2, bounds.right - resolution / 2, num=ref_dem.shape[1]),
np.linspace(bounds.bottom + resolution / 2, bounds.top - resolution / 2, num=ref_dem.shape[0])[::-1]
Expand Down Expand Up @@ -1450,7 +1451,7 @@ def residuals(coefs: np.ndarray, x_coords: np.ndarray, y_coords: np.ndarray, tar

def _apply_func(self, dem: np.ndarray, transform: rio.transform.Affine) -> np.ndarray:
"""Apply the deramp function to a DEM."""
bounds, resolution = _transform_to_bounds_and_res(dem, transform)
bounds, resolution = _transform_to_bounds_and_res(dem.shape, transform)
x_coords, y_coords = np.meshgrid(
np.linspace(bounds.left + resolution / 2, bounds.right - resolution / 2, num=dem.shape[1]),
np.linspace(bounds.bottom + resolution / 2, bounds.top - resolution / 2, num=dem.shape[0])[::-1]
Expand Down Expand Up @@ -1585,7 +1586,7 @@ def __init__(self, max_iterations: int = 50, error_threshold: float = 0.05): #
def _fit_func(self, ref_dem: np.ndarray, tba_dem: np.ndarray, transform: Optional[rio.transform.Affine],
weights: Optional[np.ndarray]):
"""Estimate the x/y/z offset between two DEMs."""
bounds, resolution = _transform_to_bounds_and_res(ref_dem, transform)
bounds, resolution = _transform_to_bounds_and_res(ref_dem.shape, transform)
verbose = False # TODO: Make this an argument somewhere.
erikmannerfelt marked this conversation as resolved.
Show resolved Hide resolved
# Make a new DEM which will be modified inplace
aligned_dem = tba_dem.copy()
Expand Down Expand Up @@ -1651,7 +1652,7 @@ def _fit_func(self, ref_dem: np.ndarray, tba_dem: np.ndarray, transform: Optiona

def _apply_func(self, dem: np.ndarray, transform: rio.transform.Affine) -> np.ndarray:
"""Apply the estimated x/y/z offsets to a DEM."""
bounds, resolution = _transform_to_bounds_and_res(dem, transform)
bounds, resolution = _transform_to_bounds_and_res(dem.shape, transform)
scaling_factor = self._meta["resolution"] / resolution

# Make index grids for the east and north dimensions
Expand Down