[gdal-dev] Handling data raster data in decimal seconds

[Moses.Gone at t-systems.com]

Hello GDAL Experts,
We have a grid dataset that is in WGS84 decimal-second coordinates. The first problem is that when we run gdalinfo against this data, the coordinates systems comes out wrong i.e in decimal degrees rather than in decimal seconds:

bash-2.05$ gdalinfo -noct /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1
Driver: AIG/Arc/Info Binary Grid
Files: /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1
       /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1.aux.xml
       /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1/hdr.adf
       /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1/metadata.xml
       /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1/prj.adf
       /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1/w001001.adf
       /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1/w001001x.adf
       /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1/dblbnd.adf
       /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1/sta.adf
       /users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/aa_ds_1/vat.adf
Size is 1324, 936
Coordinate System is:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        TOWGS84[0,0,0,0,0,0,0],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9108"]],
    AUTHORITY["EPSG","4326"]]
Origin = (5.984444444444445,50.945000000000000)
Pixel Size = (0.000277777777778,-0.000277777777778)
Corner Coordinates:
Upper Left  (   5.9844444,  50.9450000) (  5d59' 4.00"E, 50d56'42.00"N)
Lower Left  (   5.9844444,  50.6850000) (  5d59' 4.00"E, 50d41' 6.00"N)
Upper Right (   6.3522222,  50.9450000) (  6d21' 8.00"E, 50d56'42.00"N)
Lower Right (   6.3522222,  50.6850000) (  6d21' 8.00"E, 50d41' 6.00"N)
Center      (   6.1683333,  50.8150000) (  6d10' 6.00"E, 50d48'54.00"N)

Source SRS should look like this:
GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["arc-second",0.000004848136811095,AUTHORITY["EPSG","9104"]],AUTHORITY["EPSG","4326"]]

We know for sure that the data is in DS and not DD

Now the real problem comes when we want to warp this dataset to a different SRS. Since GDAL already reads a wrong SRS on the input data, all transformations end up wrong. We tried to override this by forcing a source coordinate systems on the data using gdal_translate, but to no avail:

gdal_translate -of vrt -a_srs /users/ds57201a/user1/uka/ukatiler/config/projections/ds_4326.prj aachen.vrt test.vrt

ds_4326.prj wkt:
GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["arc-second",0.000004848136811095,AUTHORITY["EPSG","9104"]],AUTHORITY["EPSG","4326"]]

gdalwarp --debug on -wm 500 -wo SKIP_NOSOURCE=YES -wo NUM_THREADS=ALL_CPUS -wo OPTIMIZE_SIZE=TRUE -tr 1.0 1.0 -of vrt -s_srs ds_4326.prj -t_srs ds_4258.prj vrt_with_target_srs.vrt

That command generated this result:

GDAL: GDALOpen(/users/ds57201a/user3/rbda/de/uka30qs/dlu/grd/617bc9d0-52ce-11e3-9de2-f1b38c98f22b/uka_parameter/vrt_with_source_srs.vrt, this=44f20) succeeds as VRT.
OGRCT: PROJ >= 4.8.0 features enabled
OGRCT: Wrap source at 6.16833.
OGRCT: Source: +proj=longlat +ellps=WGS84 +no_defs
OGRCT: Target: +proj=longlat +ellps=GRS80 +no_defs
OGRCT: Wrap target at 6.16833.
OGRCT: Source: +proj=longlat +ellps=GRS80 +no_defs
OGRCT: Target: +proj=longlat +ellps=WGS84 +no_defs
OGRCT: Wrap source at 6.16833.
OGRCT: Source: +proj=longlat +ellps=WGS84 +no_defs
OGRCT: Target: +proj=longlat +ellps=GRS80 +no_defs
OGRCT: Wrap target at 6.16833.
OGRCT: Source: +proj=longlat +ellps=GRS80 +no_defs
OGRCT: Target: +proj=longlat +ellps=WGS84 +no_defs
Creating output file that is 0P x 0L.
ERROR 1: Attempt to create 0x0 dataset is illegal,sizes must be larger than zero.

Gdalinfo on the warped dataset(when I don't get the above error) shows this:
Corner Coordinates:
- Upper Left  (   30301.000,  187402.000) (Invalid angle,Invalid angle)
Lower Left  (   30301.000,  178241.000) (Invalid angle,Invalid angle)
Upper Right (   39346.000,  187402.000) (Invalid angle,Invalid angle)
Lower Right (   39346.000,  178241.000) (Invalid angle,Invalid angle)
Center      (   34823.500,  182821.500) (Invalid angle,Invalid angle)

Is there anyone who has experienced the same thing or  knows how to deal with data in decimal seconds within GDAL?


Thanks and regards
Moses

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20131122/7eb490eb/attachment.html>