[gdal-dev] Handling data raster data in decimal seconds
Even Rouault
even.rouault at mines-paris.org
Fri Nov 22 10:38:14 PST 2013
Le vendredi 22 novembre 2013 13:12:32, Moses.Gone at t-systems.com a écrit :
> 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.0000
> 04848136811095,AUTHORITY["EPSG","9104"]],AUTHORITY["EPSG","4326"]]
>
> We know for sure that the data is in DS and not DD
I can see in the AIGRID driver that there's a special case to deal with "Units
DS" when reading in prj.adf.
I'm not sure when your source SRS comes from ?
It would be interesting that you create a ticket with the dataset attached and
the expect corner coordinates.
>
> 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.0000
> 04848136811095,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
--
Geospatial professional services
http://even.rouault.free.fr/services.html
More information about the gdal-dev
mailing list