[gdal-dev] Geocentric Coordinate System Support

[Ben Discoe]

Frank et al.,

I am testing the Geocentric support; in particular, the ability to transform
from Geocentric to a SRS with a vertical datum.

Since Geocentric was committed in March, and GDAL 1.8.1 was released in
July, i thought that 1.8.1 might contain the Geocentric support.  But, it
does not seem to be there.  The ogr_fromepsg.cpp in 1.8.1 dates from last
December.  So, for testing, i have built GDAL from today's SVN.

Basic geocentric support seems to work fine:

    OGRSpatialReference srs1, srs2;
    srs1.SetWellKnownGeogCS("WGS84");
    srs2.importFromEPSG(4978);	// geocentric (ECEF)
    OGRCoordinateTransformation *oct =
OGRCreateCoordinateTransformation(&srs1, &srs2);

    double x, y, z;
    x = -122.272778;
    y = 37.871667;
    z = 0;
    oct->Transform(1, &x, &y, &z);

This gives (-2691744, -4262401, 3894209) which is a reasonable-looking
geocentric value.  However, the first question is:

1. How is the input z (i.e. elevation) being interpreted?  The input is only
a 2D coordinate system, with no definition of what Z is.  Is it assuming
that z is an ellipsoidal height, i.e. in this case relative to the WGS84
ellipsoid?

Ellipsoidal elevation is not very useful.  To be useful, we must be able to
specify the vertical coordinates of the input, e.g. in order to interpret is
as orthometric height, i.e. relative to a geoid.

I know that the underlying support for Vertical Datums in PROJ.4 works,
using an SRS like:
	+proj=longlat +datum=WGS84 +no_defs +geoidgrids=g2009conus.gtx

However, how to do this from OGR is not clear.  It might involve
OGRSpatialReference::SetVertCS, but i could not find any example anywhere
(even googling) of how to call SetVertCS.  I tried to guess, by working
backwards:

    srs.importFromProj4("+proj=longlat +datum=WGS84 +no_defs
+geoidgrids=g2009conus.gtx");
    srs.exportToPrettyWkt(...);

Which produces:
  COMPD_CS["WGS 84 + Unnamed Vertical Datum",
    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"]],
    VERT_CS["Unnamed",
        VERT_DATUM["Unnamed",2005,
            EXTENSION["PROJ4_GRIDS","g2009conus.gtx"]],
        AXIS["Up",UP]]]

2. How to specify the vertical with OGR (using SetVertCS, or any other way)
such that it will result in PROJ.4 using the necessary GTX file?

Thanks,
Ben

> -----Original Message-----
> From: gdal-dev-bounces at lists.osgeo.org [mailto:gdal-dev-
> bounces at lists.osgeo.org] On Behalf Of Frank Warmerdam
> Sent: Wednesday, March 09, 2011 8:41 AM
> To: gdal-dev
> Subject: [gdal-dev] Geocentric Coordinate System Support
> 
> Folks,
> 
> Yesterday I introduced Geocentric coordinate system support in
> GDAL/OGR.  This includes OGRSpatialReference support for geocentric
> coordinate systems, the ability to looking them up from EPSG, and
> the ability to translate them to/from PROJ.4 format so they can be
> used for transformation.  This is an example of a geocentric coordinate
> system looked up from EPSG.  (Note the root is "GEOCCS")
> 
> GEOCCS["WGS 84 (geocentric)",
>      DATUM["WGS_1984",
>          SPHEROID["WGS 84",6378137,298.257223563,
>              AUTHORITY["EPSG","7030"]],
>          AUTHORITY["EPSG","6326"]],
>      PRIMEM["Greenwich",0,
>          AUTHORITY["EPSG","8901"]],
>      UNIT["metre",1,
>          AUTHORITY["EPSG","9001"]],
>      AXIS["Geocentric X",OTHER],
>      AXIS["Geocentric Y",OTHER],
>      AXIS["Geocentric Z",NORTH],
>      AUTHORITY["EPSG","4328"]]
> 
> I think there will be relatively little need for this coordinate system,
> but I would appreciate feedback from anyone using it.
> 
>    http://trac.osgeo.org/gdal/changeset/21916
> 
> Best regards,
> --
>
---------------------------------------+------------------------------------
-
> -
> I set the clouds in motion - turn up   | Frank Warmerdam,
warmerdam at pobox.com
> light and sound - activate the windows | http://pobox.com/~warmerdam
> and watch the world go round - Rush    | Geospatial Programmer for Rent
> 
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/gdal-dev