[gdal-dev] Setting vertical units in a raster
Even Rouault
even.rouault at spatialys.com
Thu Oct 29 13:03:20 PDT 2015
Jack,
> Hey,
> So I am exporting a raster and would like to set the VerticalUnitsGeoKey to
> feet. I have tried to use the SetUnitType( const char * ) member function
> of GDALRasterBand without success.
SetUnitType() will write a GDAL metadata in the TIFFTAG_GDAL_METADATA XML tag
> I have looked at the internal code the
> driver uses to set this flag, but would like to avoid having to reopen the
> file to do this after I have finished the export. Is there any standard
> wayto set GeoTIFF tags?
See the following Python script for an example :
from osgeo import gdal,osr
sr = osr.SpatialReference()
sr.SetFromUserInput("""COMPD_CS["unknown",
PROJCS["WGS 84 / UTM zone 11N",
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["degree",0.0174532925199433,
AUTHORITY["EPSG","9122"]],
AUTHORITY["EPSG","4326"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["latitude_of_origin",0],
PARAMETER["central_meridian",-117],
PARAMETER["scale_factor",0.9996],
PARAMETER["false_easting",500000],
PARAMETER["false_northing",0],
UNIT["metre",1,
AUTHORITY["EPSG","9001"]],
AXIS["Easting",EAST],
AXIS["Northing",NORTH],
AUTHORITY["EPSG","32611"]],
VERT_CS["Unnamed",
VERT_DATUM["North American Vertical Datum 1988",2005,
AUTHORITY["EPSG","5103"]],
UNIT["foot",0.3048,AUTHORITY["EPSG","9002"]],
AXIS["Up",UP]]]""")
ds = gdal.GetDriverByName('GTiff').Create('/vsimem/tiff_srs_compd_cs.tif',1,1)
ds.SetProjection(sr.ExportToWkt())
ds = None
# Check back
gdal.SetConfigOption('GTIFF_REPORT_COMPD_CS','YES')
ds = gdal.Open('/vsimem/tiff_srs_compd_cs.tif')
wkt = ds.GetProjectionRef()
gdal.SetConfigOption('GTIFF_REPORT_COMPD_CS',None)
print(wkt)
Even
--
Spatialys - Geospatial professional services
http://www.spatialys.com
More information about the gdal-dev
mailing list