[gdal-dev] keep metadata tags when creating multiband file
Thornton, Michele M.
thorntonmm at ornl.gov
Wed Oct 23 08:27:06 PDT 2013
Thank you Even,
Very useful information.
Kind regards,
Michele
-----Original Message-----
From: Even Rouault [mailto:even.rouault at mines-paris.org]
Sent: Wednesday, October 23, 2013 10:41 AM
To: gdal-dev at lists.osgeo.org
Cc: Thornton, Michele M.
Subject: Re: [gdal-dev] keep metadata tags when creating multiband file
Michele,
>
> I have a number of individual Landsat bands that I'm converting to
> GeoTiff and would then like to merge to a multiband file using either
> gdal_merge.py or another gdal utility (gdalbuildvrt, gdal_translate).
> In the individually derived geotif files, I am able to internally
> "label" the corresponding Landsat band by creating a metadata tag.
> For example, in using the following gdal_translate command from a file
> that is Landsat data for Band 1 called landsat.b1.dat:
>
> gdal_translate -of GTiff -strict -stats -mo "Band = LandsatBand ${f: -5:1}"
> -co "COMPRESS=LZW" -co "PROFILE=GeoTIFF" $f $f.tif
>
> This creates a metadata tag in the tif file:
> -------------------------------------------------
> Metadata:
> AREA_OR_POINT=Area
> Band =LandsatBand 1
>
> And there is a separate .aux.xml for each of the tif files:
> -------------------------------------------------------------------
> <PAMDataset>
> <Metadata>
> <MDI key="Band ">LandsatBand 1</MDI>
> </Metadata>
> <PAMRasterBand band="1">
> <Metadata>
> <MDI key="STATISTICS_MAXIMUM">255</MDI>
> <MDI key="STATISTICS_MEAN">40.826976974681</MDI>
> <MDI key="STATISTICS_MINIMUM">0</MDI>
> <MDI key="STATISTICS_STDDEV">31.471508896152</MDI>
> </Metadata>
> </PAMRasterBand>
>
Instead of setting the band name at the dataset level, you should set it at the band level instead. There's currently no way to do that with gdal_translate. You have to edit the .aux.xml file manually :
For example let's consider band1.tif with band1.tif.aux.xml :
<PAMDataset>
<PAMRasterBand band="1">
<Metadata>
<MDI key="Band">LandsatBand 1</MDI>
</Metadata>
</PAMRasterBand>
</PAMDataset>
and band2.tif with band2.tif.aux.xml :
<PAMDataset>
<PAMRasterBand band="1">
<Metadata>
<MDI key="Band">LandsatBand 2</MDI>
</Metadata>
</PAMRasterBand>
</PAMDataset>
and then you would want to run "gdalbuildvrt -separate merged.vrt band1.tif band2.tif". Unfortunately I see that it also loses band metadata. I don't think there's a fundamental reason to explain that they are lost, except that it hasn't been coded yet.
You can still put then again by editing the merged.vrt and adding the
<Metadata> <MDI key="Band">LandsatBand X</MDI> </Metadata> at
appropriate place.
See below :
<VRTDataset rasterXSize="20" rasterYSize="20">
<SRS>PROJCS["NAD27 / UTM zone
11N",GEOGCS["NAD27",DATUM["North_American_Datum_1927",SPHEROID["Clarke
1866",6378206.4,294.9786982139006,AUTHORITY["EPSG","7008"]],AUTHORITY["EPSG","6267"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4267"]],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"]],AUTHORITY["EPSG","26711"]]</SRS>
<GeoTransform> 4.4072000000000000e+05, 6.0000000000000000e+01, 0.0000000000000000e+00, 3.7513200000000000e+06, 0.0000000000000000e+00, -6.0000000000000000e+01</GeoTransform>
<VRTRasterBand dataType="Byte" band="1">
<Metadata>
<MDI key="Band">LandsatBand 1</MDI>
</Metadata>
<SimpleSource>
<SourceFilename relativeToVRT="1">byte1.tif</SourceFilename>
<SourceBand>1</SourceBand>
<SourceProperties RasterXSize="20" RasterYSize="10" DataType="Byte"
BlockXSize="20" BlockYSize="10" />
<SrcRect xOff="0" yOff="0" xSize="20" ySize="10" />
<DstRect xOff="0" yOff="0" xSize="20" ySize="10" />
</SimpleSource>
</VRTRasterBand>
<VRTRasterBand dataType="Byte" band="2">
<Metadata>
<MDI key="Band">LandsatBand 2</MDI>
</Metadata>
<SimpleSource>
<SourceFilename relativeToVRT="1">byte2.tif</SourceFilename>
<SourceBand>1</SourceBand>
<SourceProperties RasterXSize="20" RasterYSize="20" DataType="Byte"
BlockXSize="20" BlockYSize="20" />
<SrcRect xOff="0" yOff="0" xSize="20" ySize="20" />
<DstRect xOff="0" yOff="0" xSize="20" ySize="20" />
</SimpleSource>
</VRTRasterBand>
</VRTDataset>
And then
$ gdal_translate merged.vrt merged.tif
$ gdalinfo merged.tif
Driver: VRT/Virtual Raster
Files: merged.tif
Size is 20, 20
Coordinate System is:
PROJCS["NAD27 / UTM zone 11N",
GEOGCS["NAD27",
DATUM["North_American_Datum_1927",
SPHEROID["Clarke 1866",6378206.4,294.9786982139006,
AUTHORITY["EPSG","7008"]],
AUTHORITY["EPSG","6267"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4267"]],
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"]],
AUTHORITY["EPSG","26711"]]
Origin = (440720.000000000000000,3751320.000000000000000)
Pixel Size = (60.000000000000000,-60.000000000000000)
Corner Coordinates:
Upper Left ( 440720.000, 3751320.000) (117d38'28.21"W, 33d54' 8.47"N) Lower Left ( 440720.000, 3750120.000) (117d38'27.92"W, 33d53'29.51"N) Upper Right ( 441920.000, 3751320.000) (117d37'41.48"W, 33d54' 8.71"N) Lower Right ( 441920.000, 3750120.000) (117d37'41.20"W, 33d53'29.75"N)
Center ( 441320.000, 3750720.000) (117d38' 4.70"W, 33d53'49.11"N)
Band 1 Block=20x20 Type=Byte, ColorInterp=Gray
Metadata:
Band=LandsatBand 1
Band 2 Block=20x20 Type=Byte, ColorInterp=Undefined
Metadata:
Band=LandsatBand 2
Best rergards,
Even
--
Geospatial professional services
http://even.rouault.free.fr/services.html
More information about the gdal-dev
mailing list