gdalwarp odd antimeridian handling

[Plantain]

What is the bug?

I have a space-projected GeoTIFF that spans the antimeridian - ( http://static.skysight.io/raw.tiff )

I want to transform it to EPSG:3857, and clip only a square from the "middle" of the data by defining extents that cross the antimeridian. I am not sure yet if it matters if the resulting image wraps or not, whether it could be 0-360 or must be -180 to 180, the result will have overviews generated and be used for tile serving. I am trying to avoid approaches that might lose fidelity or add processing time as this is to be done in realtime.

The image looks something like this:

And I want to extract the middle:

I tried several approaches, and I am not sure if I am bouncing off bugs, undefined behaviour, or my own incompetence.

Steps to reproduce the issue

Firstly, warping to
gdalwarp -of PNG -t_srs EPSG:3857 raw.tiff out.png
Works as expected

Attempting to crop using extents beyond -180 crashes, I guess to be expected:
gdalwarp -of PNG -t_srs EPSG:3857 raw.tiff out.png -te_srs EPSG:4326 -190 -65 -90 70
ERROR 1: Invalid dfWestLongitudeDeg
Creating output file that is -39965P x 28552L.
ERROR 1: Attempt to create -39965x28552 dataset is illegal,sizes must be larger than zero.

Interesting, the same in EPSG:3857 coordinates does run! But the result is missing everything west of the antimeridian
gdalwarp -of PNG -t_srs EPSG:3857 raw.tiff out.png -te -22094105.797 -11068715.659 -10018754.171 11068715.659
ERROR 1: Point outside of projection domain
Creating output file that is 23198P x 42528L.

I then tried to get clever and specify +over explicitly on the projection, but get the same result, missing everything west of the antimeridian:
gdalwarp -of PNG -t_srs '+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext +no_defs +over' -te -22094105.797 -11068715.659 -10018754.171 11068715.659 raw.tiff over.png

Trying with CENTER_LONG seems to have no effect
gdalwarp -of PNG -t_srs 'EPSG:3857' --config CENTER_LONG 180 raw.tiff over.png

Unless I change the projection to EPSG:4326, in which case there's nothing east of the antimeridian??
gdalwarp -of PNG -t_srs 'EPSG:4326' --config CENTER_LONG 180 raw.tiff over.png

Versions and provenance

GDAL 3.10.0, released 2024/11/01
OSX via Brew

Additional context

No response

[jratike80]

I think you must clip the raster into two halves.

gdalwarp -of PNG -t_srs EPSG:3857 raw.tiff out1.png -te_srs EPSG:4326 -180 -65 -90 70
gdalwarp -of PNG -t_srs EPSG:3857 raw.tiff out2.png -te_srs EPSG:4326  170 -65 180 70

[Plantain]

I think this should not have to be the case with CENTER_LONG or lon_wrap=180?

Because of some downstream decision making my tileserving speed is linearly related to the number of files I have to check.

[jratike80]

I don't know how you could make such single wrapping image in the target CRS, EPSG:3857. The Proj library supports wrapping only with geographic coordinate systems, not with projected ones. See https://proj.org/en/stable/usage/projections.html

The +lon_wrap option can be used to provide an alternative means of doing longitude wrapping. It has only effect with operations that output angular coordinates, such as +proj=longlat

[Plantain]

Should either:

gdalwarp -of PNG -t_srs EPSG:3857 raw.tiff out.png -te_srs '+init=epsg:4326 +proj=longlat +ellps=WGS84
+datum=WGS84 +no_defs +towgs84=0,0,0 +over' -190 -65 -90 70

or

gdalwarp -of PNG -t_srs EPSG:3857 raw.tiff out.png -te_srs '+init=epsg:4326 +proj=longlat +ellps=WGS84
+datum=WGS84 +no_defs +towgs84=0,0,0 +lon_wrap=180' 170 -65 270 70

not work then?

Or is there any other way to achieve this?

[jratike80]

Please clarify what is your aim. Is it a single image like this, in EPSG:3857?

[Plantain]

Yes.

However, I also wanted to try centered on the antimeridian (like below) to compare. I was sure if they had would have the same computational complexity to generate or size to store.

[jratike80]

The projected bounds of EPSG:3857 are

-20037508.34 -20048966.1
20037508.34 20048966.1

I do not believe that there is any supported way to centre EPSG:3857 at any other longitude than 0, or to go beyond the bounds of the valid coordinate range. So either you must split the image into two halves, or to select some other coordinate system and create a custom tiling schema for that CRS.

The two half-images can be merged into one with gdal_merge or gdalbuildvrt+gdal_translate if it is essential to have one image instead of two, but then you will have lots of nodata pixels in the middle.

Note: I have never done anything with a gestatationary CRS like this and I noticed that my own commands do not really work. The CRS is +proj=geos +sweep=x +lon_0=-137 +h=35786023 +x_0=0 +y_0=0 +ellps=GRS80 +units=m +no_defs"

[jratike80]

Below are tested commands which creates two halves in EPSG:3857 which show correctly in QGIS. The images can then be converted into Web Mercator tiles according to the Google schema, if that was your final goal, The Google tiles are never split by the antimeridian so the end result will be the same from these two half-images of from a theoretical single image. Theoretical because I do not know if it is possible to create such.

gdalwarp -of GTiff -t_srs EPSG:3857 raw.tiff out1.tiff -co tiled=yes -co compress=LZW -te_srs EPSG:4326 -te -180 -65 -90 70
gdalwarp -of GTiff -t_srs EPSG:3857 raw.tiff out2.tiff -co tiled=yes -co compress=LZW -te_srs EPSG:4326 -te  170 -65 180 70
gdaladdo out1.tiff
gdaladdo out2.tiff

By creating Cloud Optimized GeoTIFFs the overviews and tiles are generated automatically but it is slower and I do not know if you would have any benefit from COG for your use case.

gdalwarp -of COG -t_srs EPSG:3857 raw.tiff out1_cog.tiff -co compress=LZW -te_srs EPSG:4326 -te -180 -65 -90 70
gdalwarp -of COG -t_srs EPSG:3857 raw.tiff out2_cog.tiff -co compress=LZW -te_srs EPSG:4326 -te  170 -65 180 70

I do not believe that you you have found any bug in GDAL so this issue may be closed soon. The right forum for asking questions about how to use GDAL is the gdal-dev mailing list. I know that in GitHub it is more convenient because it is easy to add screenshots etc.