From 7c6fefd3c20dbe5fd1ed9073bb9dba408f727ef0 Mon Sep 17 00:00:00 2001 From: Brendan Collins Date: Tue, 12 May 2026 16:46:48 -0700 Subject: [PATCH] geotiff: honour ModelTiepointTag when ModelPixelScaleTag is absent (#1750) `_extract_transform` previously returned the default `GeoTransform()` (origin (0, 0), unit pixel size) with `has_georef=True` whenever a ModelTiepointTag was present but ModelPixelScaleTag was missing. The tiepoint's tp_x / tp_y values were silently dropped, so any downstream code trusting `has_georef` constructed coordinates from origin (0, 0) and relocated the raster. Use the tiepoint's X / Y to build the origin (matching the tiepoint-with-scale branch) and fall back to the spec-documented unit pixel scale (1.0, -1.0). The alternative of returning has_georef=False would still mislabel the raster as having no georeferencing even though a real model X / Y is encoded in the tiepoint. Added `TestTiepointWithoutScale_1750` in test_geotags.py with a helper that emits a TIFF carrying ModelTiepointTag but no ModelPixelScaleTag. Verified the new tests fail on main (origin_x = 0.0 instead of the tiepoint X) and pass with the fix. Closes #1750. --- xrspatial/geotiff/_geotags.py | 24 ++++- xrspatial/geotiff/tests/test_geotags.py | 131 ++++++++++++++++++++++++ 2 files changed, 153 insertions(+), 2 deletions(-) diff --git a/xrspatial/geotiff/_geotags.py b/xrspatial/geotiff/_geotags.py index b096dc32..fe56eba6 100644 --- a/xrspatial/geotiff/_geotags.py +++ b/xrspatial/geotiff/_geotags.py @@ -463,9 +463,29 @@ def _extract_transform(ifd: IFD) -> tuple[GeoTransform, bool]: return GeoTransform(pixel_width=sx, pixel_height=-sy), True - # Tiepoint without scale: still flag as georeferenced (origin known) + # Tiepoint without scale: honour the tiepoint origin and fall back to + # unit pixel size. Per the GeoTIFF spec a ModelTiepointTag encodes a + # real-world (X, Y) for pixel (I, J); dropping it would silently relocate + # the raster to (0, 0). Unit scale (1.0, -1.0) is the documented fallback + # when ModelPixelScaleTag is absent. if tiepoint is not None: - return GeoTransform(), True + if not isinstance(tiepoint, tuple): + tiepoint = (tiepoint,) + tp_i = tiepoint[0] if len(tiepoint) > 0 else 0.0 + tp_j = tiepoint[1] if len(tiepoint) > 1 else 0.0 + tp_x = tiepoint[3] if len(tiepoint) > 3 else 0.0 + tp_y = tiepoint[4] if len(tiepoint) > 4 else 0.0 + + # Unit scale: pixel_width = 1.0, pixel_height = -1.0 + origin_x = tp_x - tp_i * 1.0 + origin_y = tp_y + tp_j * 1.0 + + return GeoTransform( + origin_x=origin_x, + origin_y=origin_y, + pixel_width=1.0, + pixel_height=-1.0, + ), True return GeoTransform(), False diff --git a/xrspatial/geotiff/tests/test_geotags.py b/xrspatial/geotiff/tests/test_geotags.py index 3f666bdb..42e68987 100644 --- a/xrspatial/geotiff/tests/test_geotags.py +++ b/xrspatial/geotiff/tests/test_geotags.py @@ -253,3 +253,134 @@ def test_z_coupling_raises(self, tmp_path): from xrspatial.geotiff._reader import read_to_array with pytest.raises(NotImplementedError): read_to_array(str(path)) + + +def _build_tiff_with_tiepoint_only(tiepoint_6: tuple) -> bytes: + """Build a tiny single-strip TIFF carrying ModelTiepointTag but + no ModelPixelScaleTag or ModelTransformationTag. + + The GeoTIFF spec permits this configuration; the tiepoint encodes a + real-world (X, Y) origin for pixel (I, J) and the pixel scale defaults + to (1.0, 1.0). + """ + import struct + + bo = '<' + width, height = 2, 2 + pixels = np.zeros((height, width), dtype=np.uint8) + + tag_list = [] + + def add_short(tag, val): + tag_list.append((tag, 3, 1, struct.pack(f'{bo}H', val))) + + def add_long(tag, val): + tag_list.append((tag, 4, 1, struct.pack(f'{bo}I', val))) + + def add_doubles(tag, vals): + tag_list.append( + (tag, 12, len(vals), struct.pack(f'{bo}{len(vals)}d', *vals))) + + add_short(256, width) # ImageWidth + add_short(257, height) # ImageLength + add_short(258, 8) # BitsPerSample + add_short(259, 1) # Compression: none + add_short(262, 1) # PhotometricInterpretation + add_short(277, 1) # SamplesPerPixel + add_short(278, height) # RowsPerStrip + add_long(273, 0) # StripOffsets (placeholder) + add_long(279, len(pixels.tobytes())) # StripByteCounts + add_short(339, 1) # SampleFormat + # ModelTiepointTag (33922): 6 doubles (I, J, K, X, Y, Z). + # Deliberately no ModelPixelScaleTag (33550). + add_doubles(33922, list(tiepoint_6)) + + tag_list.sort(key=lambda t: t[0]) + + num_entries = len(tag_list) + ifd_start = 8 + ifd_size = 2 + 12 * num_entries + 4 + + overflow = bytearray() + overflow_offsets = {} + for tag, _typ, _count, raw in tag_list: + if len(raw) > 4: + overflow_offsets[tag] = ifd_start + ifd_size + len(overflow) + overflow.extend(raw) + if len(overflow) % 2: + overflow.append(0) + + pixel_start = ifd_start + ifd_size + len(overflow) + + patched = [] + for tag, typ, count, raw in tag_list: + if tag == 273: + patched.append((tag, typ, count, struct.pack(f'{bo}I', pixel_start))) + else: + patched.append((tag, typ, count, raw)) + tag_list = patched + + out = bytearray() + out.extend(b'II') + out.extend(struct.pack(f'{bo}H', 42)) + out.extend(struct.pack(f'{bo}I', ifd_start)) + out.extend(struct.pack(f'{bo}H', num_entries)) + for tag, typ, count, raw in tag_list: + out.extend(struct.pack(f'{bo}HHI', tag, typ, count)) + if len(raw) <= 4: + out.extend(raw.ljust(4, b'\x00')) + else: + out.extend(struct.pack(f'{bo}I', overflow_offsets[tag])) + out.extend(struct.pack(f'{bo}I', 0)) + out.extend(overflow) + out.extend(pixels.tobytes()) + return bytes(out) + + +class TestTiepointWithoutScale_1750: + """Issue #1750: ModelTiepointTag present, ModelPixelScaleTag absent. + + Previously the reader returned a default GeoTransform with origin (0, 0) + while still flagging the raster as georeferenced, silently relocating it. + The fix honours the tiepoint's X/Y and falls back to unit pixel scale + (the GeoTIFF spec convention when ModelPixelScale is absent). + """ + + def test_tiepoint_origin_preserved(self, tmp_path): + # I=0, J=0 maps to (X=500000, Y=4500000) + tiepoint = (0.0, 0.0, 0.0, 500000.0, 4500000.0, 0.0) + data = _build_tiff_with_tiepoint_only(tiepoint) + path = tmp_path / 'tiepoint_no_scale_1750.tif' + path.write_bytes(data) + + header = parse_header(data) + ifds = parse_all_ifds(data, header) + from xrspatial.geotiff._geotags import _extract_transform + transform, has_georef = _extract_transform(ifds[0]) + + assert has_georef is True + assert transform.origin_x == pytest.approx(500000.0) + assert transform.origin_y == pytest.approx(4500000.0) + # Unit pixel scale fallback per GeoTIFF spec + assert transform.pixel_width == pytest.approx(1.0) + assert transform.pixel_height == pytest.approx(-1.0) + + def test_tiepoint_with_nonzero_ij(self, tmp_path): + # I=2, J=3 maps to (X=100, Y=200) -- origin shifts to (98, 203) + tiepoint = (2.0, 3.0, 0.0, 100.0, 200.0, 0.0) + data = _build_tiff_with_tiepoint_only(tiepoint) + path = tmp_path / 'tiepoint_no_scale_ij_1750.tif' + path.write_bytes(data) + + header = parse_header(data) + ifds = parse_all_ifds(data, header) + from xrspatial.geotiff._geotags import _extract_transform + transform, has_georef = _extract_transform(ifds[0]) + + assert has_georef is True + # origin_x = tp_x - tp_i * 1.0 = 100 - 2 = 98 + # origin_y = tp_y + tp_j * 1.0 = 200 + 3 = 203 + assert transform.origin_x == pytest.approx(98.0) + assert transform.origin_y == pytest.approx(203.0) + assert transform.pixel_width == pytest.approx(1.0) + assert transform.pixel_height == pytest.approx(-1.0)