Fix GeoTIFF localization grid when tie points are rectilinear but unevenly spaced

[JF3Env]

Problem

Localization.nonLinear(…) cannot build a localization grid when the model tie points are on a
rectilinear grid — every combination of the distinct x and y pixel coordinates present
exactly once — but the spacing between those coordinates is not an integer amount of pixels
repeated exactly.

LocalizationGridBuilder infers the grid size from the greatest common divisor of the tie point
coordinates (LocalizationGridBuilder.infer(…)).
That divisor is much smaller than the actual step when the step is fractional, and it collapses to
1 when a rounding to integers makes one step differ from the others by one pixel. Either way the
inferred grid becomes as large as the image, and the constructor throws ArithmeticException.

The recovery heuristic then splits the tie points in four parts around a threshold computed per
axis. It was written for Sentinel 1, whose last step is genuinely 18 pixels shorter than the
others, and it behaves badly on a rectilinear grid:

1. Crash. If the spacing is irregular on a single axis, threshold(…) returns NaN for the
   other axis, the comparisons against it are always false, two of the four parts receive no
   point, and Vector.pick(int[]) is called with an empty index array:

   IndexOutOfBoundsException: Range [0 … -1] is not valid.
   at org.apache.sis.math.Vector.subSampling(Vector.java:932)
   at org.apache.sis.math.Vector.pick(Vector.java:1158)
   at org.apache.sis.storage.geotiff.reader.Localization.localizationGrid(Localization.java:184)
   

2. Silently wrong coordinates. If a part is itself irregular it is split again, but the
   transform computed by the recursive call is added to the map of specializations only when it
   could be built without splitting — in the catch branch the result of
   MathTransforms.specialize(…) is returned and never put into addTo. The sub-area then falls
   back to an extrapolation of the largest part, far from the tie points it was supposed to honor.

Impact

Every ICEYE product takes this path: they carry no ModelPixelScaleTag and no
ModelTransformationTag, and their tie points are at k × (size-1) / (n-1) pixels, sometimes
rounded to integers. ICEYE GRD and SLC products stay in native SAR geometry and are
"tagged with ground control points (GCP) and rapid positioning capability polynomial coefficients (RPC's)";
SIS has no RPC support, so the tie points are the only georeferencing available to it.

Worst distance between a tie point as declared in the file and as computed by the transform,
measured on five products. Before is this same code with only the empty parts skipped, because
without that the products cannot be opened at all:

product	size	tie points	spacing	before	after
ScanSAR (public, see below)	19250 × 19510	39 × 40	fractional	2161.51 m	4.55 m
GRD	20000 × 20000	10 × 10	rounded, both axes	2.09 m	0.33 m
GRD	20000 × 20000	10 × 10	rounded, both axes	2.64 m	0.33 m
SLC	34484 × 15342	10 × 10	rounded, alternating	17.09 m	0.34 m
SLC	114644 × 16714	10 × 10	rounded, one axis	2.23 m	0.16 m

About 0.2 m of the remaining error on the GRD and SLC products, and about 4 m on the ScanSAR
product, is the half pixel between the tie point convention of those files
(GTRasterTypeGeoKey = RasterPixelIsArea) and PixelInCell.CELL_CORNER. It is a constant offset,
identical before and after, and out of scope here.

The residual is not evenly distributed. On the 34484 × 15342 SLC, with only the empty parts
skipped, only the diagonal of the tie point grid is honored (residual in metres):

 0 3831 7663 11494 15326 19157 22989 26820 30652 34483
0 0.18 0.18 0.12 0.15 0.43 0.75 1.19 1.68 2.26 2.90
5114 1.33 1.25 1.19 0.18 1.40 1.56 1.83 2.22 2.66 3.18
10227 7.07 6.86 6.67 6.62 6.60 6.66 0.18 7.01 7.32 7.69
15341 17.09 16.70 16.41 16.18 16.02 15.95 15.94 15.99 16.19 0.18

Fix

rectilinearGrid(…) is attempted before the split. When the tie points are a complete rectilinear
grid whose distinct coordinates are evenly spaced up to a rounding to integers, the grid is built
without splitting: the ranks of the distinct pixel coordinates are the grid indices, and the
linear relationship between pixel coordinates and ranks is applied before the localization grid.
The transform is then continuous and honors every tie point.

Grids whose spacing is irregular for a real reason — Sentinel 1 — deviate from an evenly spaced
sequence by much more than one pixel, are rejected by isUniformAfterRounding(…), and are still
handled by the existing split, unchanged. The empty parts of that split are skipped, so failure
mode 1 cannot happen even for the grids that are still split.

Failure mode 2 is not fixed here. It is reachable only by grids that are still split, for which
I have no sample product, so fixing it would ship less well validated than the rest. It is
described above so it is not lost.

Tests

LocalizationTest covers the spacings of the products above, the Sentinel 1 spacing, and a regular
grid as a control. The tie points are built from a non-linear model — with a linear model every
candidate transform reproduces the tie points exactly and the tests only verify the absence of
exception, which is why failure mode 2 is invisible to such a test.

test	grid	main	empty parts skipped	this PR
testRegularGrid	evenly spaced, integer step	pass	pass	pass
testGridWithFractionalSpacing	39 × 40, step 506.55 / 500.23	IndexOutOfBoundsException	wrong coordinates	pass
testGridWithIrregularStepOnOneAxis	rounded, one axis	IndexOutOfBoundsException	pass	pass
testGridWithIrregularStepOnBothAxes	rounded, both axes	IndexOutOfBoundsException	wrong coordinates	pass
testGridWithAlternatingSteps	rounded, alternating	IndexOutOfBoundsException	wrong coordinates	pass
testGridWithShorterLastStep	Sentinel 1	pass	pass	pass
testGridWithShorterLastStepOnOneAxis	Sentinel 1, one axis	IndexOutOfBoundsException	pass	pass

Full :endorsed:test on this branch: 3786 tests, 0 failures, 275 skipped.

How to verify

Without any data. The tests construct the tie points directly:

git fetch https://github.com/JF3Env/sis fix/SIS-000-localization-grid-single-axis-split
git checkout FETCH_HEAD
./gradlew :endorsed:test --tests "org.apache.sis.storage.geotiff.reader.LocalizationTest"# same tests against the unmodified reader, to see the failures of the middle column:
git checkout main -- endorsed/src/org.apache.sis.storage.geotiff/main/org/apache/sis/storage/geotiff/reader/Localization.java
./gradlew :endorsed:test --tests "org.apache.sis.storage.geotiff.reader.LocalizationTest"

Use JDK 24: building SIS requires Java 22 or later, and the pinned Gradle 8.14.4 runs only on
JDK 17–24.

With a public product. ESA publishes ICEYE sample products, no licence or account required.
The ScanSAR one contains a GeoTIFF and is the first row of the table above:

curl -O https://earth.esa.int/eogateway/ftp/missions/sample-data/third-party-missions/iceye/ICEYE-Scan-mode.zip
unzip ICEYE-Scan-mode.zip
# 2631255/ICEYE_X11_GRD_SC_2631255_20230902T005432.tif# sha256 31f7cd8bf7c8f5dc800a7583b763ecc3d776784274561a2c46c1e4d6305fa32d

That file is 19250 × 19510 pixels. Its ModelTiepointTag holds 1560 tie points on a 39 × 40 grid:
39 distinct x coordinates spaced by exactly 506.5526315789… pixels and 40 distinct y
coordinates spaced by exactly 500.2307692307… pixels. ModelPixelScaleTag and
ModelTransformationTag are absent, RPCCoefficientTag is present. GDAL exposes the tie points as
GCPs, so gdalinfo on it is a quick independent check of those numbers.

Opening it with GeoTiffStore and calling getGridGeometry():

	result
main	IndexOutOfBoundsException: Range [0 … -1] is not valid.
empty parts skipped only	opens; tie point residual min 4.07 m, median 421.67 m, max 2161.51 m
this PR	opens; tie point residual min 4.01 m, median 4.23 m, max 4.55 m

The ~4 m floor in the last two rows is the half pixel of anchor convention mentioned above; on this
product one pixel is about 8 m.

Checking that the tie points are honored. This transforms each tie point of a GeoTIFF through
the transform SIS built for it, and prints the worst disagreement:

vartiff = java.nio.file.Path.of(args[0]);
try (varstore = org.apache.sis.storage.DataStores.open(tiff)) {
varresource = (org.apache.sis.storage.GridCoverageResource)
((org.apache.sis.storage.Aggregate) store).components().iterator().next();
vargridToCRS = resource.getGridGeometry().getGridToCRS(
org.apache.sis.coverage.grid.PixelInCell.CELL_CORNER);
// tiePoints = the ModelTiepointTag records (I,J,K,X,Y,Z) read from the same filedoubleworst = 0;
varsource = newdouble[2];
vartarget = newdouble[2];
for (inti = 0; i < tiePoints.length; i += 6) {
source[0] = tiePoints[i];
source[1] = tiePoints[i+1];
gridToCRS.transform(source, 0, target, 0, 1);
worst = Math.max(worst, Math.hypot(target[0] - tiePoints[i+3],
target[1] - tiePoints[i+4]));
}
System.out.println(worst);
}

References

• ICEYE product documentation — Data Products, Geospatial Considerations
• ESA Earth Online — ICEYE sample data
• Changed file — Localization.java
• Branch — fix/SIS-000-localization-grid-single-axis-split (the branch name predates the wider fix; it is kept so this PR is not closed by a rename)

If a Jira issue is wanted before this is applied, I will file one and rename the branch accordingly.

---

Disclosure: this contribution was prepared with the assistance of an AI coding tool (Claude Code).
The defect was found while integrating Apache SIS into a downstream SAR application; the fix, the
tests, and every measurement above were produced and verified against the project's own test suite
and against real products, and a human maintainer of that downstream project reviewed and
authorised this submission. Please review it as you would any other outside contribution.