PR #349: Antimeridian longitude wrap + NaN from oversized tile matrices on global COGs #352
Description
Context
PR #349 (merged as 6adddc4) introduces CARTESIAN/EPSG:3857 rendering for WebMercator viewports. This works great for regional COGs but produces rendering artifacts on large global EPSG:4326 datasets where tiles span ±180° longitude.
Test file
gdalinfo /vsicurl/https://s3.us-west-2.amazonaws.com/us-west-2.opendata.source.coop/walkthru-earth/dem-terrain/GEDTM30/gedtm30.tif
1440010×600010 pixels, EPSG:4326, Int32, Deflate, 2048×2048 tiles, 10 overviews, bbox [-180.00125, -65.00125, 180.00125, 85.00125].
Issues found
1. inverseFrom3857 returns +180° at the antimeridian → 360° round-trip error
For tiles touching lon=-180°, the mesh refinement round-trip CRS → 3857 → CRS can flip longitude from -180° to +180° due to floating-point ambiguity in proj4. This causes a round-trip error equal to the full image width (5625 pixels for this file), and the mesh never converges.
Workaround: Normalize longitude in inverseFrom3857/inverseFrom4326:
if (r[0] >= 180) r[0] -= 360;
else if (r[0] < -180) r[0] += 360;2. sampleReferencePointsInEPSG3857 receives NaN from oversized tile matrices
When tileHeight > imageHeight in a COG overview (e.g. 2048px tile for a 2343px overview), the tile matrix bounds extend far beyond the data bbox (lat=-177° for this file). projectTo3857 returns NaN for these out-of-domain coordinates, which propagates through tile traversal → "Invalid number null" crash in TileLayer initialization.
Workaround: Clamp non-finite values in sampleReferencePointsInEPSG3857:
if (!Number.isFinite(projected[0])) projected[0] = 0;
if (!Number.isFinite(projected[1])) projected[1] = geoY > 0
? EPSG_3857_HALF_CIRCUMFERENCE : -EPSG_3857_HALF_CIRCUMFERENCE;3. RasterReprojector.run() has no iteration limit → browser hang
For tiles with extreme Mercator non-linearity (partial tiles spanning -65° to -177° latitude), the Delatin mesh refinement while (getMaxError() > maxError) { refine() } loop never converges and hangs the browser.
Workaround: Add iteration cap in delatin.js:
const MAX_ITERATIONS = 5000;
while (this.getMaxError() > maxError) {
this.refine();
if (++iterations >= MAX_ITERATIONS) break;
}4. Partial edge tiles produce distorted meshes
At z=0, tiles where actual data height is <25% of tile matrix height (e.g. 295px data in 2048px tile) have forwardTransform mapping UV [0,1] to CRS coordinates far beyond the data bbox. The mesh is built for the full UV range, creating distorted geometry.
Workaround: Skip rendering tiles where data.height / tileMatrix.tileHeight < 0.25.
5. Bbox slightly beyond ±180° wraps easting sign
The test file has bbox lon=±180.00125°. forwardTo3857(-180.00125, lat) wraps to positive easting (+20037369 instead of -20037647), placing the tile on the wrong side of the world.
Workaround: Clamp EPSG:4326 bbox to exactly [-180, lat_min, 180, lat_max] before calling generateTileMatrixSet.
Environment
@developmentseed/deck.gl-geotiff@0.3.0+ PR fix: Render mesh from Web Mercator coordinates #349 backport@developmentseed/deck.gl-raster@0.3.0@developmentseed/raster-reproject@0.3.0- deck.gl 9.2, MapLibre 5.20, Firefox/Chrome on macOS
