>From what I could tell, you are attempting to re-project a bounding box across the antimeridian.
Here is one way to handle that ( https://pyproj4.github.io/pyproj/stable/api/transformer.html#pyproj.transformer.Transformer.transform_bounds ): >>> import shapely.geometry >>> from pyproj import Transformer >>> transformer = Transformer.from_crs(32601, 4326, always_xy=True) >>> def bounding_polygon(left, bottom, right, top): ... if right < left: ... return shapely.geometry.MultiPolygon( ... [ ... shapely.geometry.box(left, bottom, 180, top), ... shapely.geometry.box(-180, bottom, right, top), ... ] ... ) ... return shapely.geometry.box(left, bottom, right, top) ... >>> shapely.geometry.mapping(bounding_polygon(*transformer.transform_bounds(377120, 7577600, 418080, 7618560))) {'type': 'MultiPolygon', 'coordinates': [(((180.0, 68.28494231541386), (180.0, 68.66685731991039), (179.97430985556775, 68.66685731991039), (179.97430985556775, 68.28494231541386), (180.0, 68.28494231541386)),), (((-178.98557986959838, 68.28494231541386), (-178.98557986959838, 68.66685731991039), (-180.0, 68.66685731991039), (-180.0, 68.28494231541386), (-178.98557986959838, 68.28494231541386)),)]} Alternatively, you could use rasterio.warp.transform_bounds (rasterio 1.3+ and GDAL 3.4+). https://github.com/rasterio/rasterio/issues/2313 On Tue, Mar 29, 2022 at 3:30 PM Robin Wilson <[email protected]> wrote: > Hi, > > I’ve been running into problems creating STAC items for images that cross > the antimeridian. The geometries for the images are coming out as spanning > the width of the whole world (from -180 to +180 degrees). With the help of > the author of the tool I was using to create the STAC items (rio-stac), > I’ve narrowed it down to a simple test case of reprojecting geometries > (which is what the tool is doing internally). > > If you put the following GeoJSON for an antimeridian-crossing polygon > in WGS 84 / UTM zone 1N in a file called test_32601.geojson: > > ``` > > {"type":"Polygon","coordinates":[[[377120,7577600],[418080,7577600],[418080,7618560],[377120,7618560],[377120,7577600]]]} > ``` > > and then run the following: > > ``` > ogr2ogr -f GeoJSON -s_srs epsg:32601 -t_srs epsg:4326 out_4326.geojson > test_32601.geojson > ``` > > The contents of the output file will be: > > ``` > { "type": "Polygon", "coordinates": [ [ [ -179.018099960087909, > 68.666857319910392 ], [ -178.985579869598382, 68.299719485120988 ], [ > -179.976982644969041, 68.284942315413858 ], [ 179.974309855567753, > 68.651801088224246 ], [ 180.0, 68.652259945574556 ], [ 180.0, > 68.652259945541758 ], [ 180.0, 68.652259944541754 ], [ -180.0, > 68.652259944541754 ], [ -180.0, 68.652259945541758 ], [ -180.0, > 68.652259945610126 ], [ -179.018099960087909, 68.666857319910392 ] ] ] } > ``` > > In that output geometry there are longitudes of 180 and -180, so when > plotted on a map it spans the whole width of the world. > > The author of the rio-stac library did a PR to try and fix this problem ( > https://github.com/developmentseed/rio-stac/pull/30), but it doesn’t seem > to have fixed it for this example, as the underlying GDAL code used by his > library (through rasterio, I believe) is giving this result. > > I find this result unexpected, but I have no knowledge of how GDAL deals > with the antimeridian. Is this to be expected, or is this potentially a > bug? If it isn’t a bug, is there any way around this? > > Best regards, > > Robin > > Dr Robin Wilson > www.rtwilson.com > > > _______________________________________________ > gdal-dev mailing list > [email protected] > https://lists.osgeo.org/mailman/listinfo/gdal-dev > -- Alan Snow
_______________________________________________ gdal-dev mailing list [email protected] https://lists.osgeo.org/mailman/listinfo/gdal-dev
