Rotación de imágenes hiperspectrales ENVI con GDAL

Question

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.

Answer 1

La rotación se puede cambiar desde la API de python.

Los parámetros de tamaño de cuadrícula, posición y rotación se pueden acceder con el método GetGeoTransform () y se pueden cambiar con el método SetGeoTransform. La geotransformación es una tupla de seis elementos:

 (x_origin,
cos(rotation) * x_pixel_size,
-sin(rotation) * x_pixel_size,
y_origin,
sin(rotation) * y_pixel_size,
cos(rotation) * y_pixel_size)
 

la rotación está en radianes, no grados.

Así que para su ejemplo, el código podría ser:

 import math
from osgeo import gdal

# make sure to open in update mode
dataset = gdal.Open(image_name, gdal.GA_Update)

gt = dataset.GetGeoTransform()

# convert -66 degrees to radians
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])

# should return 0
dataset.SetGeoTransform(new_geotransform)

dataset = None