En este hilo de gis.stackexchange.com se refiere la necesidad de rotar un ráster 90 grados en sentido antihorario debido a que los productos de 3B43 TRMM Multisatellite Precipitation monthly están supuestamente rotados 90 grados en sentido horario. Asumiendo que esto es cierto, con la ayuda de un algoritmo encontrado aquí, se desarrolló el codigo siguiente para rotar el array 90 grados en sentido antihorario y luego producir el raśter con los métodos de GDAL.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 | from osgeo import gdal, osr import numpy as np dataset = gdal.Open("/home/zeito/pyqgis_data/test_0_5_0_5_clip.tif") band = dataset.GetRasterBand(1) array = dataset.ReadAsArray() #rotate raster 90 degrees anticlockwise array_rot = (list(list(x) for x in zip(*matrix))[::-1]) geotransform = dataset.GetGeoTransform() wkt = dataset.GetProjection() # Create gtif file driver = gdal.GetDriverByName("GTiff") output_file = "/home/zeito/pyqgis_data/rotated_raster.tif" dst_ds = driver.Create(output_file, band.XSize, band.YSize, 1, gdal.GDT_Float32) array_rot = np.array(array_rot) #writting output raster dst_ds.GetRasterBand(1).WriteArray( array_rot ) #setting extension of output raster # top left x, w-e pixel resolution, rotation, top left y, rotation, n-s pixel resolution dst_ds.SetGeoTransform(geotransform) # setting spatial reference of output raster srs = osr.SpatialReference() srs.ImportFromWkt(wkt) dst_ds.SetProjection( srs.ExportToWkt() ) #Close output raster dataset dataset = None dst_ds = None |
Aunque se probó con muchos píxeles, para ilustrar la comprobación se escogió el de la posición (16,19) en el ráster original. Este corresponde al valor de 803; tal como se desprende de la imagen siguiente con la ayuda del plugin Value Tool.
Este mismo valor, una vez realizada la rotación de 90 grados en sentido antihorario, debería estar en la posición (16,0) del ráster rotado. Es lo que se corrobora en la imagen siguiente:
No hay comentarios:
Publicar un comentario