gdalwarp on MSG-SEVIRI fails

[hrobjartur]

I am having some problems reprojecting the resultant MSG geotiff images into polar stereographic.

I have an example geotiff file of MSG channel 9 that I have produced as follows. gdal_translate "MSG($PWD,201112071530,(9),Y,T)" msg_channel9.tif,

The file can be downloaded here:

​http://brunnur.vedur.is/pub/hroi/msg_channel9.tif

A simple reprojection request such as $ gdalwarp -t_srs "+proj=stere" msg_channel9.tif stereo.tif always results in some sort of bounding box error,

---

ERROR 1: Too many points (5624 out of 5625) failed to transform, unable to compute output bounds.

---

I have tried also to limit the output grid to something sensible using -tr and -te but I still get such errors.

These errors occur for all projections, however polar stereographic projection fails with above error. Gdalinfo on the file also suggest some problems with the bounding-box.

I am not clear if this is a bug with gdalwarp, or with the MSG SEVIRI driver.

[Even Rouault]

The issue is not with gdalwarp. The issue is that GDAL doesn't know what PROJECTIONGeostationary_Satellite is ...

$ gdalinfo --debug off /vsicurl/http://brunnur.vedur.is/pub/hroi/msg_channel9.tif
Driver: GTiff/GeoTIFF
Files: /vsicurl/http://brunnur.vedur.is/pub/hroi/msg_channel9.tif
Size is 3712, 3712
Coordinate System is:
PROJCS["Geostationary_Satellite",
    GEOGCS["GCS_WGS_1984",
        DATUM["D_WGS_1984",
            SPHEROID["WGS_1984",6378137,298.257223563]],
        PRIMEM["Greenwich",0],
        UNIT["Degree",0.017453292519943295]],
    PROJECTION["Geostationary_Satellite"],
    PARAMETER["central_meridian",0],
    PARAMETER["satellite_height",35785831],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["Meter",1]]
Origin = (-5568748.275756835937500,5568748.275756835937500)
Pixel Size = (3000.403165817260742,-3000.403165817260742)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
ERROR 1: tolerance condition error
Upper Left  (-5568748.276, 5568748.276) 
ERROR 1: tolerance condition error
Lower Left  (-5568748.276,-5568748.276) 
ERROR 1: tolerance condition error
Upper Right ( 5568748.276, 5568748.276) 
ERROR 1: tolerance condition error
Lower Right ( 5568748.276,-5568748.276) 
Center      (   0.0000000,   0.0000000) (  0d 0' 0.01"E,  0d 0' 0.01"N)
Band 1 Block=3712x1 Type=Float32, ColorInterp=Gray

[retsios]

Hopefully some hints:

When using gdalwarp from gdal 1.6 without the -s_srs parameter, it gives the following error message, but after that gdalwarp proceeds and makes an exact copy of the original image.

ERROR 1: No PROJ.4 translation for source SRS, coordinate
transformation initialization has failed.

When using gdalwarp from gdal 1.6 and forcing the -s_srs parameter, I do get a result:

gdalwarp -s_srs "+proj=geos" -t_srs "+proj=stere" msg_channel9.tif stereo.tif
Creating output file that is 3690P x 3734L.
Processing input file msg_channel9.tif.
0...10...20...30...40...50...60...70...80...90...100 - done.

When using gdal 1.7 to do the same (with the -s_srs parameter), it fails:

ERROR 1: Too many points (440 out of 441) failed to transform,
unable to compute output bounds.

The GEOS projection was fairly new when I added the MSG driver. The OGRSpatialReference::SetGEOS() member function was added specifically to cover the projection of these images, at a time that it was still unknown to EPSG.

[hrobjartur]

Thanks, this is very useful to know.

It turns out that if I set a more complete proj string, as I guess I should have,

-t_srs "+proj=stere +lat_ts=60.0 +lat_0=40.0"

It seems that this is the reason why my warping failed in gdal 1.8.1. The image now gets projected, but with many ERROR 1: tolerance condition error messages. These tolerance error messages actually occur when I project to any other projection it seems.

So I guess the only question that remains is what are the

"ERROR 1: tolerance condition error"

and should I worry about them?

I have managed to avoid these errors by warping onto an already existing geotiff...

btw. thank you for the helpful comments - I am now well on the way to use gdal as a very fast and efficient processor for MSG data here at the Icelandic Met. Office.

[Even Rouault]

Please ignore my first comment about the Geostationary_Satellite projection not being supported. Seems it is well supported by the whole GDAL, libgeotiff and proj4 stack !

As far as the "tolerance condition error", it is generally the sign that a coordinate transformation occured outside the validity area of the projection. Which is the likely case of the 4 corner of the source raster, which are in the black area of the geotiff. When you do a gdalinfo, gdal tries to reproject them in geodetic coordinates, which fails. The same occurs when gdalwarp tries to compute the extent of the output file, when it doesn't yet exist. It tries different approach. First, it takes 21 sample points on each of the 4 edges of the source image and then reproject them in the output SRS. If it fails, it tries on a grid on the whole raster (which is at the minimum of 21 x 21 = 441 points, but since recent GDAL version, it is (MAX(20, MIN(xsize, ysize)/50))+1)^{2, so in your case 75*75 = 5625 points). If the number of points of that grid that fail to be reprojected to the target SRS is close to the total number of points, then the warping will abort. If there are at least a few points that are correctly reprojected, then it will goes on.}

Since recent version, GDAL also check that the reverse transforming gives coordinates close to the original. This improves the detection of it you are in the area of validity or not. This can be disabled by specifying --config CHECK_WITH_INVERT_PROJ FALSE, but this is generally a bad idea because the result will likely be junk because GDAL could not guess correctly the extent of the target raster (I've tried with --config CHECK_WITH_INVERT_PROJ FALSE with -t_srs "+proj=stere", and gdalwarps doesn't abort, but gives a full black dataset). In that case, you need to explicitely define the target extent with the -te parameter.

So, to sum it up, you must a target SRS that is consistant with the spatial extent of the source raster. -t_srs "+proj=stere" is not, but -t_srs "+proj=stere +lat_ts=60.0 +lat_0=40.0" is.

All in all, I don't think there's a bug. It is just that gdalwarp cannot warp automatically in any projection automatically, or it must be guided.