locationtech / geotrellis

GeoTrellis is a geographic data processing engine for high performance applications.
http://geotrellis.io
Other
1.34k stars 361 forks source link

Anti-meridian crossing polygons reproject incorrectly from Sinusoidal to LatLng (EPSG:4326) #2913

Open philvarner opened 5 years ago

philvarner commented 5 years ago

Reproduction case is in https://github.com/philvarner/geotrellis-sinusoidaltolatlngreprojection

This case is specifically when working with the geometry of a MODIS scene that crosses the antimeridian. When doing the reprojection using code like

polygon.densify(10000).reproject(Sinusoidal, LatLng).toJson

the resulting geojson vector looks like:

Screen Shot 2019-05-01 at 6 46 48 PM

notice the zero-width part of the polygon on the right edge. This is invalid, as there are many duplicated points, which is visually apparent from looking at the geojson text. My guess is that what's happening is that all the parts of the polygon that should on the other side of the antimeridian are getting collapsed into a line along the antimeridian, resulting in an invalid polygon.

In the example code, reading the geojson output with a rigorous parser (spatial4j and JTS) results in the error:

org.locationtech.spatial4j.exception.InvalidShapeException: Self-intersection at or near point (180.0, 67.14285713896875, NaN)

When doing this in gdal with:

aws s3 cp s3://modis-pds/MCD43A4.006/25/02/2019111/MCD43A4.A2019111.h25v02.006.2019120034054_B01.TIF .
gdaltindex MCD43A4.A2019111.h25v02.006.2019120034054_B01.shp MCD43A4.A2019111.h25v02.006.2019120034054_B01.TIF
ogr2ogr -segmentize 25000 -t_srs EPSG:4326 -wrapdateline MCD43A4.A2019111.h25v02.006.2019120034054_B01-clipped.shp MCD43A4.A2019111.h25v02.006.2019120034054_B01.shp
ogr2ogr -t_srs EPSG:4326 -lco RFC7946=YES -lco WRITE_BBOX=YES MCD43A4.A2019111.h25v02.006.2019120034054_B01-clipped.geojson MCD43A4.A2019111.h25v02.006.2019120034054_B01-clipped.shp

I get two polygons split across the antimeridian:

Screen Shot 2019-05-01 at 7 26 53 PM
pomadchin commented 5 years ago

This can definitely be a proj4j bug.