[gdal-dev] gdal_translate - ECW format and non-UTM projections

[Armin Burger]

Hi all

I am using "gdal_translate" to convert GeoTiff or Erdas images to ECW. This works fine when using standard projections like UTM-WGS84. For UTM I guess GDAL even writes the correct projection information into the ECW header. Several problems appear however if I convert images that are in one of the various European national projections.  
    
The way I'm trying to assure GDAL uses the projection information it understands is to call gdal_translate with a flag like

  -a_srs "EPSG:3912"
  
This usually creates an XML file (*.ecw.aux.xml) that contains the projection information in an XML tag <SRS>. But then some quite strange issues arise that I don't understand:

-------------

1. "gdalinfo" displays the projection correctly, but the values for the corner coordinates are just the pixel values, starting with 0/0 for upper-left. This for the above listed projection 3912 gives eg. for an ECW file

Size is 40258, 20016
Coordinate System is:
PROJCS["MGI 1901 / Slovene National Grid",
    GEOGCS["MGI 1901",
        DATUM["MGI_1901",
            SPHEROID["Bessel 1841",6377397.155,299.1528128,
                AUTHORITY["EPSG","7004"]],
            AUTHORITY["EPSG","1031"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.01745329251994328,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","3906"]],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",15],
    PARAMETER["scale_factor",0.9999],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",-5000000],
    AUTHORITY["EPSG","3912"],
    AXIS["Y",EAST],
    AXIS["X",NORTH]]
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0,20016.0)
Upper Right (40258.0,    0.0)
Lower Right (40258.0,20016.0)
Center      (20129.0,10008.0)


In programs like ArcGIS 10 the coordinates are displayed correctly and the projection information is also right. But in QGIS I have just like with gdalinfo the pixel coordinates even though the projection information is correct (displayed in proj4 notation). 

Using such an image as layer in MapServer via WMS, defining the projection like
   "init=EPSG:3912"
the image has more or less the correct projected coordinates, but datum shift parameters defined in the epsg file are not taken into account. If I instead of the ECW file use the GeoTIFF version instead, the returned imagettes from the WMS are correct with datum shift parameters applied.  

2. If I remove the *.ecw.aux.xml file then gdalinfo does not find any more any projection information but at least displays correctly the corner coordinates. The same for QGIS.

Also MapServer WMS then returns the images correctly with datum shift parameters applied.


3. An even more exotic problem appears eg. with an image in European ETM projection (the EU-specific UTM projection using ETRS89 datum, but in practical terms +/- identical to UTM WGS84).

If I convert a GeoTiff with such a projection (EPSG 25832)
PROJCS["UTM Zone 32, Northern Hemisphere",
    GEOGCS["ETRF 1989",
        DATUM["ETRF 1989",
            SPHEROID["WGS 84",6378137,298.2572235630016],
            TOWGS84[0,0,0,0,0,0,0]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    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["Meter",1]]

the resulting ECW shows the same problems mentioned above, just pixel coordinates displayed. But in addition, it does not display any pixels in some programs like ERDAS, ER Viewer or gives errors in ERDAS IWS.... It doesn't matter if I keep the aux.xml file or if I remove it. 

If I only re-write the projection information in the GeoTiff header to the equivalent UTM using GDAL Python bindings, convert the image again to ECW with gdal_translate, then the resulting ECW file is perfectly correct and displays fine in all programs I tested.

-----------

The GDAL version was 1.8.0, using the ECW library 3.3, running on Linux. The same problems appear though on Windows using GDAL 1.7.3 from FWTools  or a newer 1.8.0 compilation from www.gisinternals.com.
 
Did others stumble over similar types of problem as well? Is there any possibility this could be sorted out somehow? If needed I can provide some test images that show the problems via FTP download.

Best regards
Armin



-- 
NEU: FreePhone - kostenlos mobil telefonieren!			
Jetzt informieren: http://www.gmx.net/de/go/freephone