Thank you Even,
Very useful information.
Kind regards,
Michele
-Original Message-
From: Even Rouault [mailto:even.roua...@mines-paris.org]
Sent: Wednesday, October 23, 2013 10:41 AM
To: gdal-dev@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_MAXIMUM255/MDI
MDI key=STATISTICS_MEAN40.826976974681/MDI
MDI key=STATISTICS_MINIMUM0/MDI
MDI key=STATISTICS_STDDEV31.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=BandLandsatBand 1/MDI
/Metadata
/PAMRasterBand
/PAMDataset
and band2.tif with band2.tif.aux.xml :
PAMDataset
PAMRasterBand band=1
Metadata
MDI key=BandLandsatBand 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=BandLandsatBand X/MDI /Metadata at
appropriate place.
See below :
VRTDataset rasterXSize=20 rasterYSize=20
SRSPROJCS[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,50],PARAMETER[false_northing,0],UNIT[metre,1,AUTHORITY[EPSG,9001]],AUTHORITY[EPSG,26711]]/SRS
GeoTransform 4.4072e+05, 6.e+01,
0.e+00, 3.75132000e+06, 0.e+00,
-6.e+01/GeoTransform
VRTRasterBand dataType=Byte band=1
Metadata
MDI key=BandLandsatBand 1/MDI
/Metadata
SimpleSource
SourceFilename relativeToVRT=1byte1.tif/SourceFilename
SourceBand1/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=BandLandsatBand 2/MDI
/Metadata
SimpleSource
SourceFilename relativeToVRT=1byte2.tif/SourceFilename
SourceBand1/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,50],
PARAMETER[false_northing,0],
UNIT[metre,1,
AUTHORITY[EPSG,9001]],
AUTHORITY[EPSG,26711]]
Origin = (440720.000,3751320.000)
Pixel Size
