[gdal-dev] assigning vertical datum to image files

Even Rouault even.rouault at spatialys.com
Tue Jun 9 06:32:47 PDT 2015


Selon "Newcomb, Doug" <doug_newcomb at fws.gov>:

> Thanks!
>
> That should work in this case to build the virtual raster.  I'm still
> curious about assigning the vertical datum to the rasters.

Doug,

I'm a bit surprised that the HFA driver reports a SRS with a VERTCS node as a
child of the PROJCS node. I believe this is an incorrect WKT definition. There
should rather be a COMPD_CS root note with 2 children: a PROJCS node for the
horizontal CRS and a VERTCS node for the vertical CRS. Perhaps this SRS string
is defined in the .aux.xml sidecar files that are apparently present ?

Regarding your question, there are a few EPSG codes that correspond to
aCOMPD_CS, for example 6349. Otherwise, GDAL allows to build a custom COMPD_CS
with the syntax EPSG:XXXX+YYYY where XXXX is the EPSG code of the horizontal SRS
and YYYY the EPSG code of the vertical SRS. In your case that should be
EPSG:6543+5703 I think

So you can do things like "gdal_translate in.tif out.tif -a_srs EPSG:6543+5703"

If you read back out.tif, you will only get the horizontal SRS by default. This
was for backward compatibility reason I think. But if you specify "--config
GTIFF_REPORT_COMPD_CS YES" you will get the COMPD_CS.

Even

>
> Doug
>
> On Tue, Jun 9, 2015 at 1:48 AM, Hermann Peifer <peifer at gmx.eu> wrote:
>
> >
> > About your initial problem:
> > > but gdalbuildvrt complained, informing me that
> > > they were not in the same projection.
> >
> > What about using the below gdalbuildvrt option?
> >
> > Hermann
> >
> > -allow_projection_difference:
> > (starting with GDAL 1.7.0) When this option is specified, the utility will
> > accept to make a VRT even if the input datasets have not the same
> > projection. Note: this does not mean that they will be reprojected. Their
> > projection will just be ignored.
> >
> > Source: http://www.gdal.org/gdalbuildvrt.html
> >
> >
> >
> > On 2015-06-08 21:01, Newcomb, Doug wrote:
> >
> >> Hi Folks,
> >> I have a directory of 783 .img format Lidar - based  DEMs in the same
> >> projection, North Carolina State Plane  Feet (2011) , NAD 83 , NVD88.  I
> >> was going to use gdalbuildvrt to create a single virtual image for the
> >> area, but gdalbuildvrt complained, informing me that they were not in
> >> the same projection.
> >>
> >> A couple of quick bash scripts/commands
> >>
> >> for x in *.img; do gdalinfo $x >>img_info.txt;done
> >>
> >> and
> >> grep PROJCS img_info.txt|sort|uniq -c
> >>
> >> gives me the following:
> >>
> >>      437 PROJCS["NAD_1983_2011_StatePlane_North_Carolina_FIPS_3200_Ft_US",
> >>      346 PROJCS["NAD_1983_StatePlane_North_Carolina_FIPS_3200_Feet_2011",
> >>
> >> gdalinfo gives the following for each type of file ( origin section
> >> clipped out) :
> >>
> >> Driver: HFA/Erdas Imagine Images (.img)
> >> Files: D05_37_20878102_20141231.img
> >>         D05_37_20878102_20141231.img.aux.xml
> >> Size is 1000, 1000
> >> Coordinate System is:
> >> PROJCS["NAD_1983_2011_StatePlane_North_Carolina_FIPS_3200_Ft_US",
> >>      GEOGCS["GCS_NAD_1983_2011",
> >>          DATUM["NAD_1983_2011",
> >>              SPHEROID["GRS_1980",6378137.0,298.257222101]],
> >>          PRIMEM["Greenwich",0.0],
> >>          UNIT["Degree",0.017453292519943295]],
> >>      PROJECTION["Lambert_Conformal_Conic_2SP"],
> >>      PARAMETER["False_Easting",2000000.0],
> >>      PARAMETER["False_Northing",0.0],
> >>      PARAMETER["Central_Meridian",-79.0],
> >>      PARAMETER["Standard_Parallel_1",34.3333333333333],
> >>      PARAMETER["Standard_Parallel_2",36.1666666666667],
> >>      PARAMETER["Latitude_Of_Origin",33.75],
> >>      UNIT["Foot_US",0.30480060960121924],
> >>      VERTCS["NAVD_1988_Foot_US",
> >>          VDATUM["North_American_Vertical_Datum_1988"],
> >>          PARAMETER["Vertical_Shift",0.0],
> >>          PARAMETER["Direction",1.0],
> >>          UNIT["Foot_US",0.3048006096012192]]]
> >>
> >>
> >> Driver: HFA/Erdas Imagine Images (.img)
> >> Files: D05_37_20957301_20141231.img
> >>         D05_37_20957301_20141231.img.aux.xml
> >> Size is 1000, 1000
> >> Coordinate System is:
> >> PROJCS["NAD_1983_StatePlane_North_Carolina_FIPS_3200_Feet_2011",
> >>      GEOGCS["GCS_NAD_1983_2011",
> >>          DATUM["NAD_1983_2011",
> >>              SPHEROID["GRS_1980",6378137.0,298.257222101]],
> >>          PRIMEM["Greenwich",0.0],
> >>          UNIT["Degree",0.0174532925199433]],
> >>      PROJECTION["Lambert_Conformal_Conic_2SP"],
> >>      PARAMETER["False_Easting",2000000.002616666],
> >>      PARAMETER["False_Northing",0.0],
> >>      PARAMETER["Central_Meridian",-79.0],
> >>      PARAMETER["Standard_Parallel_1",34.33333333333334],
> >>      PARAMETER["Standard_Parallel_2",36.16666666666666],
> >>      PARAMETER["Latitude_Of_Origin",33.75],
> >>      UNIT["Foot_US",0.30480060960121924],
> >>      VERTCS["NAVD_1988_Foot_US",
> >>          VDATUM["North_American_Vertical_Datum_1988"],
> >>          PARAMETER["Vertical_Shift",0.0],
> >>          PARAMETER["Direction",1.0],
> >>          UNIT["Foot_US",0.3048006096012192]]]
> >>
> >>
> >>
> >> In theory, these should all be EPSG:6543, so just assigning the correct
> >>   horizontal projection/datum with the epsg code should make things
> >> usable. (i.e,
> >> gdal_translate -a_"srs epsg:6543" --config GDAL_CACHEMAX 512 -of GTiff
> >> -co COMPRESS=DEFLATE -co PREDICTOR=3 in.img out.tif )  ( Thank you for
> >> the EPSG:6543 projection support in GDAL 2.0!)
> >>
> >> However, this only assigns the horizontal infomation, how does one
> >> assign a vertical datum with a horizontal EPSG code?
> >>
> >> I see the NVD88 code in vertcs.csv , but how would I implement it in the
> >> above command?
> >>
> >> Using gdal 2.0 Beta2
> >>
> >>
> >> Doug
> >> --
> >> Doug Newcomb
> >> USFWS
> >> Raleigh, NC
> >> 919-856-4520 ext. 14 doug_newcomb at fws.gov <mailto:doug_newcomb at fws.gov>
> >>
> >>
>
---------------------------------------------------------------------------------------------------------
> >> The opinions I express are my own and are not representative of the
> >> official policy of the U.S.Fish and Wildlife Service or Dept. of the
> >> Interior.   Life is too short for undocumented, proprietary data formats.
> >>
> >>
> >> _______________________________________________
> >> gdal-dev mailing list
> >> gdal-dev at lists.osgeo.org
> >> http://lists.osgeo.org/mailman/listinfo/gdal-dev
> >>
> >>
> >
>
>
> --
> Doug Newcomb
> USFWS
> Raleigh, NC
> 919-856-4520 ext. 14 doug_newcomb at fws.gov
>
---------------------------------------------------------------------------------------------------------
> The opinions I express are my own and are not representative of the
> official policy of the U.S.Fish and Wildlife Service or Dept. of the
> Interior.   Life is too short for undocumented, proprietary data formats.
>


-- 
Spatialys - Geospatial professional services
http://www.spatialys.com


More information about the gdal-dev mailing list