[Gdal-dev] reading linear units from GeoTIFF broken?

[Vladimir Slepnev]

Hello.

I'm really sorry about bugging you again...
but linear units don't seem to be read correctly from the GeoTIFF format.

Here's what I do.
take a PROJ.4 string:

+proj=utm +zone=13 +datum=WGS84 +units=ft

convert it to WKT via an OGRSpatialReference:

PROJCS["UTM Zone 13, Northern Hemisphere",GEOGCS["WGS
84",DATUM["WGS_1984",SPHEROID["WGS
84",6378137,298.257223563,AUTHORITY["EPSG",7030]],TOWGS84[0,0,0,0,0,0,0],AUT
HORITY["EPSG",6326]],PRIMEM["Greenwich",0,AUTHORITY["EPSG",8901]],UNIT["degr
ee",0.0174532925199433,AUTHORITY["EPSG",9108]],AXIS["Lat","NORTH"],AXIS["Lon
g","EAST"],AUTHORITY["EPSG",4326]],PROJECTION["Transverse_Mercator"],PARAMET
ER["latitude_of_origin",0],PARAMETER["central_meridian",-105],PARAMETER["sca
le_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northi
ng",0],UNIT["Foot (International)",0.3048]]

write it into a newly created GeoTIFF file through GDAL (SetProjection()),
here's a piece of listgeo info for the newly written file:

   Keyed_Information:
      GTModelTypeGeoKey (Short,1): ModelTypeProjected
      GTRasterTypeGeoKey (Short,1): RasterPixelIsArea
      GTCitationGeoKey (Ascii,33): "UTM Zone 13, Northern Hemisphere"
      GeographicTypeGeoKey (Short,1): GCS_WGS_84
      ProjectedCSTypeGeoKey (Short,1): User-Defined
      ProjectionGeoKey (Short,1): Proj_UTM_zone_13N
      ProjLinearUnitsGeoKey (Short,1): Linear_Foot
      End_Of_Keys.
   End_Of_Geotiff.

read the same file through GDAL, get the WKT projection string
(GetProjectionRef()):

PROJCS["unnamed",GEOGCS["WGS
84",DATUM["WGS_1984",SPHEROID[,6378137,298.2572235630016,AUTHORITY["EPSG",70
30]],AUTHORITY["EPSG",6326]],PRIMEM["Greenwich",0],UNIT["degree",0.017453292
5199433],AUTHORITY["EPSG",4326]],PROJECTION["Transverse_Mercator"],PARAMETER
["latitude_of_origin",0],PARAMETER["central_meridian",-105],PARAMETER["scale
_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing
",0],UNIT["unknown",1]]

convert it to PROJ.4 via an OGRSpatialReference:

+proj=utm +zone=13 +ellps=WGS84 +datum=WGS84 +units=m

That is, the units got written by GDAL, but not read back correctly.
Do you have any idea what the matter is?
(I am using today's version of GDAL from CVS.)

Once again, sorry for the inconvenience...
our correspondence seems to be a major part of gdal-dev by now :-)
And thanks a lot for your help with all this stuff.

Vladimir Slepnev
Programmer, RDC ScanEx