Skip to content

Commit

Permalink
Fix blockwise coregistration direction
Browse files Browse the repository at this point in the history
  • Loading branch information
rhugonnet committed Dec 6, 2024
1 parent 71c4214 commit 0ac2800
Show file tree
Hide file tree
Showing 4 changed files with 16 additions and 10 deletions.
3 changes: 2 additions & 1 deletion .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -153,5 +153,6 @@ doc/source/gen_modules/
doc/source/sg_execution_times.rst
examples/basic/temp.tif

# Directory where myst_nb executes jupyter code
# Directory where myst_nb executes jupyter code and cache
doc/jupyter_execute/
doc/.jupyter_cache/
15 changes: 10 additions & 5 deletions doc/source/coregistration.md
Original file line number Diff line number Diff line change
Expand Up @@ -495,17 +495,22 @@ For instance, for {class}`xdem.coreg.Deramp`, an input `poly_order` to define th
for {class}`xdem.coreg.DirectionalBias`, an input `angle` to define the angle at which to do the directional correction.

## Dividing coregistration in blocks

### The {class}`~xdem.coreg.BlockwiseCoreg` object

Sometimes, we want to split a coregistration across different spatial subsets of an elevation dataset, running that

```{caution}
The {class}`~xdem.coreg.BlockwiseCoreg` feature is still experimental: it might not support all coregistration
methods, and create edge artefacts.
```

Sometimes, we want to split a coregistration across different spatial subsets of an elevation dataset, running that
method independently in each subset. A {class}`~xdem.coreg.BlockwiseCoreg` can be constructed for this:

```{code-cell} ipython3
blockwise = xdem.coreg.BlockwiseCoreg(xdem.coreg.NuthKaab(), subdivision=16)
```

The subdivision corresponds to an equal-length block division across the extent of the elevation dataset. It needs
The subdivision corresponds to an equal-length block division across the extent of the elevation dataset. It needs
to be a number of the form 2{sup}`n` (such as 4 or 256).

It is run the same way as other coregistrations:
Expand All @@ -528,4 +533,4 @@ ax[0].set_title("Before block NK")
ax[1].set_title("After block NK")
(aligned_dem - ref_dem).plot(cmap='RdYlBu', vmin=-30, vmax=30, ax=ax[1], cbar_title="Elevation differences (m)")
_ = ax[1].set_yticklabels([])
```
```
2 changes: 1 addition & 1 deletion tests/test_coreg/test_base.py
Original file line number Diff line number Diff line change
Expand Up @@ -950,7 +950,7 @@ def test_blockwise_coreg_large_gaps(self) -> None:
ddem_post = (aligned - self.ref).data.compressed()
ddem_pre = (tba - self.ref).data.compressed()
assert abs(np.nanmedian(ddem_pre)) > abs(np.nanmedian(ddem_post))
# assert np.nanstd(ddem_pre) > np.nanstd(ddem_post)
assert np.nanstd(ddem_pre) > np.nanstd(ddem_post)


class TestAffineManipulation:
Expand Down
6 changes: 3 additions & 3 deletions xdem/coreg/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -3393,8 +3393,8 @@ def _apply_rst(
warped_dem = warp_dem(
dem=elev,
transform=transform,
source_coords=all_points[:, :, 0],
destination_coords=all_points[:, :, 1],
source_coords=all_points[:, :, 1],
destination_coords=all_points[:, :, 0],
resampling="linear",
)

Expand Down Expand Up @@ -3538,7 +3538,7 @@ def warp_dem(
if not no_vertical:
grid_offsets = scipy.interpolate.griddata(
points=destination_coords_scaled[:, :2],
values=destination_coords_scaled[:, 2] - source_coords_scaled[:, 2],
values=source_coords_scaled[:, 2] - destination_coords_scaled[:, 2],
xi=(grid_x, grid_y),
method=resampling,
fill_value=np.nan,
Expand Down

0 comments on commit 0ac2800

Please sign in to comment.