[GRASS-dev] [GRASS GIS] #1988: towgs84 datum transform options not maintained with shapefile export

GRASS GIS GRASS GIS developer
,
1

#1988: towgs84 datum transform options not maintained with shapefile export
-------------------------------+--------------------------------------------
Reporter: voncasec | Owner: grass-dev@…
     Type: defect | Status: new
Priority: normal | Milestone: 6.4.3
Component: Default | Version: unspecified
Keywords: g.proj, v.out.ogr | Platform: Linux
      Cpu: x86-64 |
-------------------------------+--------------------------------------------
Perhaps this issue is related to the GDAL ticket 4345
(http://trac.osgeo.org/gdal/ticket/4345).

Including it here in case GRASS handles the data differently.

In general, if I am using a projection defined with towgs84 datum
transform parameters and I want to export the layer to a shapefile, the
+towgs84 datum transform parameters are not included.

Often times, I will just manually modify the .prj file after the fact to
include the TOWGS84 flag to ensure that there is no shift when viewing the
locations in a GIS that allow for on the fly reprojections.

Example of g.proj -e

GRASS 6.4.3svn (Bayaguana_NAD27_Caribbean):/usr/local/grass-6.4.3svn/etc >
g.proj -e
PROJCS["UTM_Zone_19_Northern_Hemisphere",
     GEOGCS["GCS_clark66",
         DATUM["D_North_American_1927",
             SPHEROID["Clarke_1866",6378206.4,294.9786982]],
         PRIMEM["Greenwich",0],
         UNIT["Degree",0.017453292519943295]],
     PROJECTION["Transverse_Mercator"],
     PARAMETER["latitude_of_origin",0],
     PARAMETER["central_meridian",-69],
     PARAMETER["scale_factor",0.9996],
     PARAMETER["false_easting",500000],
     PARAMETER["false_northing",0],
     UNIT["Meter",1]]

Desired output of g.proj -e + .prj file from v.out.ogr

GRASS 6.4.3svn (Bayaguana_NAD27_Caribbean):/usr/local/grass-6.4.3svn/etc >
g.proj -e
PROJCS["UTM_Zone_19_Northern_Hemisphere",
     GEOGCS["GCS_clark66",
         DATUM["D_North_American_1927",
             SPHEROID["Clarke_1866",6378206.4,294.9786982],
             TOWGS84[-3.0,142.0,183.0,0,0,0,0]],
         PRIMEM["Greenwich",0],
         UNIT["Degree",0.017453292519943295]],
     PROJECTION["Transverse_Mercator"],
     PARAMETER["latitude_of_origin",0],
     PARAMETER["central_meridian",-69],
     PARAMETER["scale_factor",0.9996],
     PARAMETER["false_easting",500000],
     PARAMETER["false_northing",0],
     UNIT["Meter",1]]

--
Ticket URL: <https://trac.osgeo.org/grass/ticket/1988&gt;
GRASS GIS <http://grass.osgeo.org>

2

#1988: towgs84 datum transform options not maintained with shapefile export
--------------------------------+-------------------------------------------
Reporter: voncasec | Owner: grass-dev@…
     Type: defect | Status: new
Priority: normal | Milestone: 6.4.3
Component: Projections/Datums | Version: unspecified
Keywords: g.proj, v.out.ogr | Platform: Linux
      Cpu: x86-64 |
--------------------------------+-------------------------------------------
Changes (by hamish):

  * component: Default => Projections/Datums

Comment:

Yes, it's a known issue, and it needs to be solved in GDAL/OGR first. If I
recall correctly, WKT doesn't have a defined way to include datum
transforms, ISTR some time ago (older than that ticket) Frank suggested it
could be packed into something like a comment term, but then that would
only be known to other osgeo-stack software which knew to look for it.

Also, the proj4 epsg file includes explicit +towgs84 terms for some SRSs,
but others not (several options for the user to choose from under the
umbrella of the generalized code); the one referred to in the gdal ticket
(epsg:4632) is one of the ones with explicit +towgs84 terms in the epsg
def'n.

?,
Hamish

--
Ticket URL: <https://trac.osgeo.org/grass/ticket/1988#comment:1&gt;
GRASS GIS <http://grass.osgeo.org>

3

Thanks Hamish,

I assumed it needed to be looked at with GDAL/OGR first when exporting to a shapefile, I just wasn't sure if this was also the case when using g.proj -e.

On 13-06-03 04:19 PM, GRASS GIS wrote:

#1988: towgs84 datum transform options not maintained with shapefile export
--------------------------------+-------------------------------------------
  Reporter: voncasec | Owner: grass-dev@…
      Type: defect | Status: new
  Priority: normal | Milestone: 6.4.3
Component: Projections/Datums | Version: unspecified
  Keywords: g.proj, v.out.ogr | Platform: Linux
       Cpu: x86-64 |
--------------------------------+-------------------------------------------
Changes (by hamish):

   * component: Default => Projections/Datums

Comment:

  Yes, it's a known issue, and it needs to be solved in GDAL/OGR first. If I
  recall correctly, WKT doesn't have a defined way to include datum
  transforms, ISTR some time ago (older than that ticket) Frank suggested it
  could be packed into something like a comment term, but then that would
  only be known to other osgeo-stack software which knew to look for it.

  Also, the proj4 epsg file includes explicit +towgs84 terms for some SRSs,
  but others not (several options for the user to choose from under the
  umbrella of the generalized code); the one referred to in the gdal ticket
  (epsg:4632) is one of the ones with explicit +towgs84 terms in the epsg
  def'n.

  ?,
  Hamish

--

4

#1988: towgs84 datum transform options not maintained with shapefile export
--------------------------------+-------------------------------------------
Reporter: voncasec | Owner: grass-dev@…
     Type: defect | Status: new
Priority: normal | Milestone: 6.4.3
Component: Projections/Datums | Version: unspecified
Keywords: g.proj, v.out.ogr | Platform: Linux
      Cpu: x86-64 |
--------------------------------+-------------------------------------------

Comment(by hamish):

Casey wrote: (grass-dev ML)
> Thanks Hamish,
>
> I assumed it needed to be looked at with GDAL/OGR first when exporting
> to a shapefile, I just wasn't sure if this was also the case when using
> g.proj -e.

Hi,

I did a little digging. You get the TOWGS84 in the WKT if you use 'g.proj
-w' and your PERMANENT/PROJ_INFO file has a towgs84 line, but not 'g.proj
-we' (-e modifies -w, and triggers it because it implies it). Grid datum
transform files are not handled by any variant of Well Known Text afaik (a
slight correction to my last comment), except perhaps by using a custom
osgeo-family "EXTENSION" keyword, which I'm not sure has been commonly
defined yet.

So the good news is that for generic Shapefile export you could leave off
the v.out.ogr/g.prok '-e' flag to use generic WKT, which supports TOWGS84.
For export to ESRI products it would be lost though, as ESRI apparently
(according to the OGR code) doesn't support TOWGS84 statements in the WKT.

specifically 'g.proj -we' calls GPJ_grass_to_wkt() which if the ESRI-
compatibility flag is given calls OGR's OSRMorphToESRI() (aka
morphToESRI()).

morphToESRI() has this comment as it strips away any TOWGS84:
{{{
/* -------------------------------------------------------------------- */
/* Strip all CT parameters (AXIS, AUTHORITY, TOWGS84, etc). */
/* -------------------------------------------------------------------- */
}}}

OGR's StripCTParms() is called in to do that:

{{{
  * This method will remove all components of the coordinate system
  * that are specific to the OGC CT Specification. That is it will attempt
  * to strip it down to being compatible with the Simple Features 1.0
  * specification.
...
     poCurrent->StripNodes( "AUTHORITY" );
     poCurrent->StripNodes( "TOWGS84" );
     poCurrent->StripNodes( "AXIS" );
     poCurrent->StripNodes( "EXTENSION" );

     return OGRERR_NONE;
}
}}}

see grass's lib/proj/convert.c and GDAL/OGR's ogr/ogr_srs_esri.cpp +
ogrspatialreference.cpp.

so in the case of 'g.proj -w -e' it behaves as expected, TOWGS84 is
removed on purpose.

hope it helps,
Hamish

--
Ticket URL: <https://trac.osgeo.org/grass/ticket/1988#comment:2&gt;
GRASS GIS <http://grass.osgeo.org>

5

That does help, thanks.

Casey

On 13-06-03 09:34 PM, GRASS GIS wrote:

#1988: towgs84 datum transform options not maintained with shapefile export
--------------------------------+-------------------------------------------
  Reporter: voncasec | Owner: grass-dev@…
      Type: defect | Status: new
  Priority: normal | Milestone: 6.4.3
Component: Projections/Datums | Version: unspecified
  Keywords: g.proj, v.out.ogr | Platform: Linux
       Cpu: x86-64 |
--------------------------------+-------------------------------------------

Comment(by hamish):

  Casey wrote: (grass-dev ML)
  > Thanks Hamish,
  >
  > I assumed it needed to be looked at with GDAL/OGR first when exporting
  > to a shapefile, I just wasn't sure if this was also the case when using
  > g.proj -e.

  Hi,

  I did a little digging. You get the TOWGS84 in the WKT if you use 'g.proj
  -w' and your PERMANENT/PROJ_INFO file has a towgs84 line, but not 'g.proj
  -we' (-e modifies -w, and triggers it because it implies it). Grid datum
  transform files are not handled by any variant of Well Known Text afaik (a
  slight correction to my last comment), except perhaps by using a custom
  osgeo-family "EXTENSION" keyword, which I'm not sure has been commonly
  defined yet.

  So the good news is that for generic Shapefile export you could leave off
  the v.out.ogr/g.prok '-e' flag to use generic WKT, which supports TOWGS84.
  For export to ESRI products it would be lost though, as ESRI apparently
  (according to the OGR code) doesn't support TOWGS84 statements in the WKT.

  specifically 'g.proj -we' calls GPJ_grass_to_wkt() which if the ESRI-
  compatibility flag is given calls OGR's OSRMorphToESRI() (aka
  morphToESRI()).

  morphToESRI() has this comment as it strips away any TOWGS84:
  {{{
  /* -------------------------------------------------------------------- */
  /* Strip all CT parameters (AXIS, AUTHORITY, TOWGS84, etc). */
  /* -------------------------------------------------------------------- */
  }}}

  OGR's StripCTParms() is called in to do that:

  {{{
   * This method will remove all components of the coordinate system
   * that are specific to the OGC CT Specification. That is it will attempt
   * to strip it down to being compatible with the Simple Features 1.0
   * specification.
  ...
      poCurrent->StripNodes( "AUTHORITY" );
      poCurrent->StripNodes( "TOWGS84" );
      poCurrent->StripNodes( "AXIS" );
      poCurrent->StripNodes( "EXTENSION" );

      return OGRERR_NONE;
  }
  }}}

  see grass's lib/proj/convert.c and GDAL/OGR's ogr/ogr_srs_esri.cpp +
  ogrspatialreference.cpp.

  so in the case of 'g.proj -w -e' it behaves as expected, TOWGS84 is
  removed on purpose.

  hope it helps,
  Hamish

--