[Gdal-dev] gdalwarp shift problem

[Markus Neteler]

Dear list members,

at time I struggle to reproject a LANDSAT TM 5 subscene with
gdalwarp (from UTM to Gauss-Boaga). The resulting map is
shifted around 920km to West and also slightly distorted.

The input map 'lsat1_utm.tif' was created with gdal_translate from
the original LANDSAT TM 5 scene (from Geotiff to Geotiff). The
smaller subscene is located in Nothern Italy and looks like this:

gdalinfo lsat1_utm.tif
Driver: GTiff/GeoTIFF
Size is 3816, 3461
Coordinate System is:
PROJCS["WGS 84 / UTM zone 32N",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.2572235629972,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",9],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AUTHORITY["EPSG","32632"]]
Origin = (636861.000000,5152686.000000)
Pixel Size = (28.500000,-28.500000)
Corner Coordinates:
Upper Left  (  636861.000, 5152686.000) ( 10d47'2.74"E, 46d30'49.78"N)
Lower Left  (  636861.000, 5054047.500) ( 10d45'20.60"E, 45d37'35.33"N)
Upper Right (  745617.000, 5152686.000) ( 12d12'0.22"E, 46d28'58.64"N)
Lower Right (  745617.000, 5054047.500) ( 12d 8'57.27"E, 45d35'47.57"N)
Center      (  691239.000, 5103366.750) ( 11d28'20.39"E, 46d 3'25.67"N)
Band 1 Block=3816x2 Type=Byte, ColorInterp=Gray


The Lat/Long corners look o.k. But then, after gdalwarp:

#26591: Monte Mario (Rome) / Italy zone 1
gdalwarp -t_srs '+init=epsg:26591 +towgs84=-104.1,-49.1,-9.9,0.971,-2.917,0.714,-11.68' \
  lsat1_utm.tif lsat1_gb.tif

The datum comes from the official site:
http://crs.ifag.de/

gdalinfo lsat1_gb2.tif
Driver: GTiff/GeoTIFF
Size is 4311, 4015
Coordinate System is:
PROJCS["Monte Mario (Rome) / Italy zone 1",
    GEOGCS["Monte Mario (Rome)",
        DATUM["Monte_Mario_Rome",
            SPHEROID["International 1924",6378388,297.000000000005,
                AUTHORITY["EPSG","7022"]],
            AUTHORITY["EPSG","6806"]],
        PRIMEM["Rome",12.45233333333333],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4806"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",-3.45233333333333],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",1500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AUTHORITY["EPSG","26591"]]
Origin = (666266.468468,5206815.783256)
Pixel Size = (28.697167,-28.697167)
Corner Coordinates:
Upper Left  (  666266.468, 5206815.783) ( 14d19'12.01"W, 46d29'41.38"N)
Lower Left  (  666266.468, 5091596.659) ( 14d 7'12.46"W, 45d28'34.07"N)
Upper Right (  789979.954, 5206815.783) ( 12d43'50.79"W, 46d38'10.47"N)
Lower Right (  789979.954, 5091596.659) ( 12d33'32.53"W, 45d36'45.44"N)
Center      (  728123.211, 5149206.221) ( 13d25'57.64"W, 46d 3'27.81"N)
Band 1 Block=4311x1 Type=Byte, ColorInterp=Gray

Looking at LatLong here, it is heavily misplaced (especially shifted to
West).

To me "Central Meridian" in the Gauss-Boaga file looks wrong.

Looking into
/usr/local/share/proj/epsg
I am not sure if
+lon_0=21.45233333333333
is right. Shouldn't it be 
+lon_0=9.0
?

Thanks for any hint,

 Markus Neteler