Skip to content

Fixed a half-pixel shift in affine transformation #20

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

Merged
konstntokas merged 7 commits into from
Jan 26, 2026
Merged

Fixed a half-pixel shift in affine transformation #20

Show file tree
Hide file tree
Changes from 6 commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
ec69ecb
test need to be fixed
konstntokas Jan 8, 2026
8ddb8ad
many tests failing
konstntokas Jan 8, 2026
d96df8a
chekout
konstntokas Jan 9, 2026
6b7c10c
downscale is wrong
konstntokas Jan 9, 2026
b03b1ef
read for review
konstntokas Jan 12, 2026
2e9f91e
read for review
konstntokas Jan 12, 2026
edcee1b
address yogesh's review
konstntokas Jan 19, 2026
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
4 changes: 4 additions & 0 deletions CHANGES.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,10 @@
- Improved performance when creating `GridMapping` from 2D coordinates.
Parameter estimation now uses only corner points instead of the full set of
coordinates, significantly reducing runtime for large datasets.
- Bug fix in affine transformations.
Fixed a half-pixel shift that occurred in `xcube_resampling.affine_transform_dataset`
and `xcube_resampling.resample_in_space` when the source and target grid mappings
shared the same CRS and were regular.


## Changes in 0.2.3
Expand Down
102 changes: 55 additions & 47 deletions tests/gridmapping/test_base.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -42,7 +42,7 @@ def __init__(self, **kwargs):
tile_size=self.tile_size,
is_j_axis_up=self.is_j_axis_up,
xy_res=self.xy_res,
xy_min=(self.xy_bbox[0], self.xy_bbox[1]),
xy_min=(self.xy_bbox[0] + self.x_res / 2, self.xy_bbox[1] + self.x_res / 2),
crs=self.crs,
)

Expand Down Expand Up @@ -83,8 +83,13 @@ def kwargs(self, **kwargs):
except TypeError:
x_res, y_res = 2 * (orig_kwargs["xy_res"],)
x_min, y_min = orig_kwargs.pop("xy_min")
x_max, y_max = x_min + x_res * width, y_min + y_res * height
orig_kwargs["xy_bbox"] = x_min, y_min, x_max, y_max
x_max, y_max = x_min + x_res * (width - 1), y_min + y_res * (height - 1)
orig_kwargs["xy_bbox"] = (
x_min - x_res / 2,
y_min - y_res / 2,
x_max + x_res / 2,
y_max + y_res / 2,
)
return orig_kwargs

def test_valid(self):
Expand All @@ -98,10 +103,10 @@ def test_valid(self):
self.assertEqual(180, gm.tile_height)
self.assertEqual((0, 0, 720, 360), gm.ij_bbox)
self.assertEqual((-180.0, -90.0, 180.0, 90.0), gm.xy_bbox)
self.assertEqual(-180.0, gm.x_min)
self.assertEqual(-90.0, gm.y_min)
self.assertEqual(180.0, gm.x_max)
self.assertEqual(90.0, gm.y_max)
self.assertEqual(-179.75, gm.x_min)
self.assertEqual(-89.75, gm.y_min)
self.assertEqual(179.75, gm.x_max)
self.assertEqual(89.75, gm.y_max)
self.assertEqual((0.5, 0.5), gm.xy_res)
self.assertEqual(0.5, gm.x_res)
self.assertEqual(0.5, gm.y_res)
Expand All @@ -113,6 +118,7 @@ def test_valid(self):

self.assertIsInstance(gm.xy_coords, xr.DataArray)
np.testing.assert_equal(
gm.ij_bboxes,
np.array(
[
[0, 0, 360, 180],
Expand All @@ -121,9 +127,9 @@ def test_valid(self):
[360, 180, 720, 360],
]
),
gm.ij_bboxes,
)
np.testing.assert_equal(
gm.xy_bboxes,
np.array(
[
[-180.0, 0.0, 0.0, 90.0],
Expand All @@ -132,7 +138,6 @@ def test_valid(self):
[0.0, -90.0, 180.0, 0.0],
]
),
gm.xy_bboxes,
)

def test_invalids(self):
Expand Down Expand Up @@ -176,58 +181,61 @@ def test_ij_to_xy_transform(self):
**self.kwargs(size=(1200, 1200), xy_min=(0, 0), xy_res=1, crs=NOT_A_GEO_CRS)
)
i2crs = image_geom.ij_to_xy_transform
self.assertMatrixPoint((0, 0), i2crs, (0, 1200))
self.assertMatrixPoint((1024, 0), i2crs, (1024, 1200))
self.assertMatrixPoint((0, 1024), i2crs, (0, 1200 - 1024))
self.assertMatrixPoint((1024, 1024), i2crs, (1024, 1200 - 1024))
self.assertEqual(((1, 0, 0), (0.0, -1, 1200)), i2crs)
self.assertMatrixPoint((0, 0), i2crs, (0, 1199))
self.assertMatrixPoint((1024, 0), i2crs, (1024, 1199))
self.assertMatrixPoint((0, 1024), i2crs, (0, 1199 - 1024))
self.assertMatrixPoint((1024, 1024), i2crs, (1024, 1199 - 1024))
self.assertEqual(((1, 0, 0), (0.0, -1, 1199)), i2crs)

image_geom = _TestGridMapping(
**self.kwargs(size=(1440, 720), xy_min=(-180, -90), xy_res=0.25)
**self.kwargs(size=(1440, 720), xy_min=(-179.875, -89.875), xy_res=0.25)
)
i2crs = image_geom.ij_to_xy_transform
self.assertMatrixPoint((-180, 90), i2crs, (0, 0))
self.assertMatrixPoint((0, 0), i2crs, (720, 360))
self.assertMatrixPoint((180, -90), i2crs, (1440, 720))
self.assertEqual(((0.25, 0.0, -180.0), (0.0, -0.25, 90.0)), i2crs)
self.assertMatrixPoint((-179.875, 89.875), i2crs, (0, 0))
self.assertMatrixPoint((0.125, -0.125), i2crs, (720, 360))
self.assertMatrixPoint((179.875, -89.875), i2crs, (1439, 719))
self.assertEqual(((0.25, 0.0, -179.875), (0.0, -0.25, 89.875)), i2crs)

image_geom = _TestGridMapping(
**self.kwargs(
size=(1440, 720), xy_min=(-180, -90), xy_res=0.25, is_j_axis_up=True
size=(1440, 720),
xy_min=(-179.875, -89.875),
xy_res=0.25,
is_j_axis_up=True,
)
)
i2crs = image_geom.ij_to_xy_transform
self.assertMatrixPoint((-180, -90), i2crs, (0, 0))
self.assertMatrixPoint((0, 0), i2crs, (720, 360))
self.assertMatrixPoint((180, 90), i2crs, (1440, 720))
self.assertEqual(((0.25, 0.0, -180.0), (0.0, 0.25, -90.0)), i2crs)
self.assertMatrixPoint((-179.875, -89.875), i2crs, (0, 0))
self.assertMatrixPoint((0.125, 0.125), i2crs, (720, 360))
self.assertMatrixPoint((179.875, 89.875), i2crs, (1439, 719))
self.assertEqual(((0.25, 0.0, -179.875), (0.0, 0.25, -89.875)), i2crs)

def test_xy_to_ij_transform(self):
image_geom = _TestGridMapping(
**self.kwargs(size=(1200, 1200), xy_min=(0, 0), xy_res=1, crs=NOT_A_GEO_CRS)
)
crs2i = image_geom.xy_to_ij_transform
self.assertMatrixPoint((0, 0), crs2i, (0, 1200))
self.assertMatrixPoint((1024, 0), crs2i, (1024, 1200))
self.assertMatrixPoint((0, 1024), crs2i, (0, 1200 - 1024))
self.assertMatrixPoint((1024, 1024), crs2i, (1024, 1200 - 1024))
self.assertEqual(((1, 0, 0), (0.0, -1, 1200)), crs2i)
self.assertMatrixPoint((0, 0), crs2i, (0, 1199))
self.assertMatrixPoint((1024, 0), crs2i, (1024, 1199))
self.assertMatrixPoint((0, 1024), crs2i, (0, 1199 - 1024))
self.assertMatrixPoint((1024, 1024), crs2i, (1024, 1199 - 1024))
self.assertEqual(((1, 0, 0), (0.0, -1, 1199)), crs2i)

image_geom = _TestGridMapping(**self.kwargs(size=(1440, 720), xy_res=0.25))
crs2i = image_geom.xy_to_ij_transform
self.assertMatrixPoint((0, 720), crs2i, (-180, -90))
self.assertMatrixPoint((720, 360), crs2i, (0, 0))
self.assertMatrixPoint((1440, 0), crs2i, (180, 90))
self.assertEqual(((4.0, 0.0, 720.0), (0.0, -4.0, 360.0)), crs2i)
self.assertMatrixPoint((0, 719), crs2i, (-179.875, -89.875))
self.assertMatrixPoint((720, 360), crs2i, (0.125, -0.125))
self.assertMatrixPoint((1439, 0), crs2i, (179.875, 89.875))
self.assertEqual(((4.0, 0.0, 719.5), (0.0, -4.0, 359.5)), crs2i)

image_geom = _TestGridMapping(
**self.kwargs(size=(1440, 720), xy_res=0.25, is_j_axis_up=True)
)
crs2i = image_geom.xy_to_ij_transform
self.assertMatrixPoint((0, 0), crs2i, (-180, -90))
self.assertMatrixPoint((720, 360), crs2i, (0, 0))
self.assertMatrixPoint((1440, 720), crs2i, (180, 90))
self.assertEqual(((4.0, 0.0, 720.0), (0.0, 4.0, 360.0)), crs2i)
self.assertMatrixPoint((0, 0), crs2i, (-179.875, -89.875))
self.assertMatrixPoint((720, 360), crs2i, (0.125, 0.125))
self.assertMatrixPoint((1439, 719), crs2i, (179.875, 89.875))
self.assertEqual(((4.0, -0.0, 719.5), (-0.0, 4.0, 359.5)), crs2i)

def test_ij_transform_to_and_from(self):
gm1 = _TestGridMapping(
Expand All @@ -239,16 +247,16 @@ def test_ij_transform_to_and_from(self):
)
)
self.assertEqual(
((10.0, 0.0, -7600.0), (0.0, 10.0, -5600.0)), gm1.ij_transform_to(gm2)
((10.0, 0.0, -7595.0), (0.0, 10.0, -5595.0)), gm1.ij_transform_to(gm2)
)
self.assertEqual(
((10.0, 0.0, -7600.0), (0.0, 10.0, -5600.0)), gm2.ij_transform_from(gm1)
((10.0, -0.0, -7595.0), (-0.0, 10.0, -5595.0)), gm2.ij_transform_from(gm1)
)
self.assertEqual(
((0.1, 0.0, 760.0), (0.0, 0.1, 560.0)), gm2.ij_transform_to(gm1)
((0.1, 0.0, 759.5), (0.0, 0.1, 559.5)), gm2.ij_transform_to(gm1)
)
self.assertEqual(
((0.1, 0.0, 760.0), (0.0, 0.1, 560.0)), gm1.ij_transform_from(gm2)
((0.1, 0.0, 759.5), (0.0, 0.1, 559.5)), gm1.ij_transform_from(gm2)
)

def assertMatrixPoint(self, expected_point, matrix, point):
Expand Down Expand Up @@ -361,8 +369,8 @@ def test_transform_xy_res(self):
self.assertEqual(pyproj.CRS.from_string("EPSG:32633"), transformed_gm.crs)
self.assertEqual((400, 200), transformed_gm.size)
self.assertEqual((200, 200), transformed_gm.tile_size)
self.assertAlmostEqual(605.482, transformed_gm.xy_res[0])
self.assertAlmostEqual(1109.46, transformed_gm.xy_res[1])
self.assertAlmostEqual(605.566, transformed_gm.xy_res[0])
self.assertAlmostEqual(1109.47, transformed_gm.xy_res[1])
self.assertEqual(False, transformed_gm.is_j_axis_up)
self.assertEqual(
("transformed_x", "transformed_y"), transformed_gm.xy_var_names
Expand All @@ -375,10 +383,10 @@ def test_transform_xy_res(self):
self.assertEqual(
pyproj.CRS.from_string("EPSG:32633"), transformed_gm_regular.crs
)
self.assertEqual((440, 410), transformed_gm_regular.size)
self.assertEqual((439, 223), transformed_gm_regular.size)
self.assertEqual((200, 200), transformed_gm_regular.tile_size)
self.assertAlmostEqual(605.482, transformed_gm.xy_res[0])
self.assertAlmostEqual(1109.46, transformed_gm.xy_res[1])
self.assertAlmostEqual(605.566, transformed_gm_regular.xy_res[0])
self.assertAlmostEqual(1109.47, transformed_gm_regular.xy_res[1])
self.assertEqual(False, transformed_gm_regular.is_j_axis_up)
self.assertEqual(("x", "y"), transformed_gm_regular.xy_var_names)
self.assertEqual(("x", "y"), transformed_gm_regular.xy_dim_names)
Expand All @@ -400,7 +408,7 @@ def test_to_regular(self):
self.assertEqual(
pyproj.CRS.from_string("EPSG:32633"), transformed_gm_regular.crs
)
self.assertEqual((1103, 1551), transformed_gm_regular.size)
self.assertEqual((1102, 1046), transformed_gm_regular.size)
self.assertEqual((1000, 1000), transformed_gm_regular.tile_size)
self.assertEqual(False, transformed_gm_regular.is_j_axis_up)
self.assertEqual(False, transformed_gm_regular.is_lon_360)
Expand Down
12 changes: 6 additions & 6 deletions tests/gridmapping/test_coords.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -271,10 +271,10 @@ def test_2d_regular(self):
self.assertEqual((4, 3), gm.tile_size)
self.assertAlmostEqual(0.1, gm.x_res)
self.assertAlmostEqual(0.2, gm.y_res)
self.assertAlmostEqual(10.15, gm.x_min)
self.assertAlmostEqual(52.3, gm.y_min)
self.assertAlmostEqual(10.55, gm.x_max)
self.assertAlmostEqual(52.9, gm.y_max)
self.assertAlmostEqual(10.2, gm.x_min)
self.assertAlmostEqual(52.4, gm.y_min)
self.assertAlmostEqual(10.5, gm.x_max)
self.assertAlmostEqual(52.8, gm.y_max)
self.assertEqual(GEO_CRS, gm.crs)
self.assertEqual(True, gm.is_regular)
self.assertEqual(True, gm.is_j_axis_up)
Expand Down Expand Up @@ -317,9 +317,9 @@ def test_to_regular(self):
gm_irr = GridMapping.from_coords(lon, lat, GEO_CRS)
gm_reg_actual = gm_irr.to_regular()
gm_reg_expected = GridMapping.regular(
size=(4, 4),
size=(3, 3),
tile_size=(2, 2),
xy_min=(-1.75, 48.25),
xy_min=(0, 50),
xy_res=3.5,
crs=GEO_CRS,
)
Expand Down
2 changes: 1 addition & 1 deletion tests/gridmapping/test_dataset.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -103,7 +103,7 @@ def test_from_real_olci(self):
self.assertEqual(("coord", "y", "x"), gm.xy_coords.dims)

gm = gm.to_regular()
self.assertEqual((2867, 2133), gm.size)
self.assertEqual((1636, 2132), gm.size)
Comment thread
b-yogesh marked this conversation as resolved.

def test_from_sentinel_2(self):
dataset = create_s2plus_dataset()
Expand Down
Loading