Estoy desarrollando una aplicación con la GDAL biblioteca para combinar ráster georreferenciado las imágenes hiperespectrales en el MEDIO del formato de la NASA Aire Visual/Infrared Imaging Spectrometer (AVIRIS) con otros conjuntos de datos geográficos. Al importar las imágenes en una aplicación de SIG como QGIS, me parece que es girado incorrectamente:
El flightline debe ir de suroeste a noreste, pero es muestra de sur a norte. El map info
de campo de la comisión ENVI archivo de encabezado contiene un rotation
parámetro de -66° que supongo que es relativo a la parte superior izquierda de píxeles:
$ grep "map info" ang20150422t163638_corr_v1e_img_zero.hdr
map info = { UTM , 1.000 , 1.000 , 736600.089 , 4078126.750 , 2.7000000000e+00 , 2.7000000000e+00 , 12 , North , WGS-84 , units=Meters , rotation=-66.00000000 }
Sin embargo, este campo es no documentados y al parecer es ignorado por GDAL. Aquí está la actual proyección de la imagen:
$ gdalsrsinfo ang20150422t163638_corr_v1e_img_zero.img
PROJ.4 : '+proj=utm +zone=12 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs '
OGC WKT :
PROJCS["UTM Zone 12, 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],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0,
AUTHORITY["EPSG","8901"]],
UNIT["degree",0.0174532925199433,
AUTHORITY["EPSG","9108"]],
AUTHORITY["EPSG","4326"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["latitude_of_origin",0],
PARAMETER["central_meridian",-111],
PARAMETER["scale_factor",0.9996],
PARAMETER["false_easting",500000],
PARAMETER["false_northing",0],
UNIT["Meter",1]]
Lo que me gustaría hacer es girar la imagen utilizando el GDAL utilidades o, preferentemente, de la GDAL API de Python por lo que está orientado con el norte arriba sin cambiar nada más acerca de la imagen.
He leído que las rotaciones se puede lograr con una combinación de gdalwarp
y gdal_translate
, ya sea mediante la aplicación de un geotransform o el cambio de la +towgs84
parámetros.
Sin embargo, yo soy un SIG novato y buscar orientación específica sobre qué parámetros necesarios para lograr esto.
Usted puede encontrar una banda única versión de la imagen y su encabezado aquí.
Primer Intento
Basado en lpdudleysugerencia, he hecho una copia de la imagen:
$ cp ang20150422t163638_corr_v1e_img_zero.hdr ang20150422t163638_corr_v1e_img_zero_rot.hdr
$ cp ang20150422t163638_corr_v1e_img_zero.img ang20150422t163638_corr_v1e_img_zero_rot.img
Y entonces se aplica el geotransform en el interprete de Python:
>>> image_name = "ang20150422t163638_corr_v1e_img_zero_rot.img"
>>> dataset = gdal.Open(image_name, gdal.GA_Update)
>>> gt = dataset.GetGeoTransform()
>>> rotation = (math.pi/180) * -66
>>> new_geotransform = (gt[0],
math.cos(rotation) * gt[1],
-math.sin(rotation) * gt[1],
gt[3],
math.sin(rotation) * gt[5],
math.cos(rotation) * gt[5])
... ... ... ... ... >>> dataset.SetGeoTransform(new_geotransform)
0
>>> dataset = None
Sin embargo, la imagen se escala en lugar de girar:
Error
Rotación de la comisión ENVI parece ser un error en GDAL:
- https://lists.osgeo.org/pipermail/gdal-dev/2013-January/035146.html
- https://trac.osgeo.org/gdal/ticket/1778
El error parece provenir de la comisión ENVI conductor. Si encuentro una solución que voy a compartir aquí.
Revisión
He publicado una solicitud de extracción y compartido en el GDAL lista de correo utilizando los valores sugeridos por lpdudley. Ahora, cuando me aplican gdalwarp
a la imagen:
$ gdalwarp -of envi ang20150422t163638_corr_v1e_img_zero.img ang20150422t163638_corr_v1e_img_zero_rot.img
Creating output file that is 4048P x 2461L.
Processing input file ang20150422t163638_corr_v1e_img_zero.img.
0...10...20...30...40...50...60...70...80...90...100 - done.
La rotación será detectado correctamente:
Esperemos que las futuras versiones de la biblioteca de la manija de la girado ENVI archivos correctamente.