Estoy tratando de grabar un shapefile en un ráster utilizando GDAL del RasterizeLayer. Me pre-crear un área de interés de la trama desde un shapefile, dado un determinado tamaño de píxel. Esta ZONA sirve como una base para todos los siguientes rasterizations (el mismo número de collumns y filas, la misma proyección y geotransform).
El problema se produce, sin embargo, cuando voy a grabar las formas de su propia trama, basada en el mismo tamaño de píxel y proyecciones. El enlace de más abajo (no tienen suficiente rep de subir la imagen), se muestra el original shapefile en color beige y el rosa oscuro donde RasterizeLayer ha quemado datos. La luz de color rosa es nodata los valores de la rosa oscuro datos raster. El gris es el AOI basa en que el shapefile de grabación se completó.
Dada la extensión de los shapefile de polígonos, yo esperaría a ver quemar los valores en la parte inferior dos esquinas, así como los dos píxeles debajo de los datos que muestra. Obviamente, sin embargo, este no es el caso.
Imagen para el Problema de - Terminado Trama Quemaduras
La siguiente es el código que he utilizado para generar estas. Todas las figuras fueron creadas usando QGIS, y fueron creados en la misma proyección. (Cabe señalar que la disposición en cuadrícula en la imagen que se muestra fue sólo para dar una idea de el tamaño de píxel de que estaba utilizando.)
from osgeo import ogr
from osgeo import gdal
aoi_uri = 'AOI_Raster.tif'
aoi_raster = gdal.Open(aoi_uri)
def new_raster_from_base(base, outputURI, format, nodata, datatype):
cols = base.RasterXSize
rows = base.RasterYSize
projection = base.GetProjection()
geotransform = base.GetGeoTransform()
bands = base.RasterCount
driver = gdal.GetDriverByName(format)
new_raster = driver.Create(str(outputURI), cols, rows, bands, datatype)
new_raster.SetProjection(projection)
new_raster.SetGeoTransform(geotransform)
for i in range(bands):
new_raster.GetRasterBand(i + 1).SetNoDataValue(nodata)
new_raster.GetRasterBand(i + 1).Fill(nodata)
return new_raster
shape_uri = 'activity_3.shp'
shape_datasource = ogr.Open(shape_uri)
shape_layer = shape_datasource.GetLayer()
raster_out = 'new_raster.tif'
raster_dataset = new_raster_from_base(aoi_raster, raster_out, 'GTiff',
-1, gdal.GDT_Int32)
band = raster_dataset.GetRasterBand(1)
nodata = band.GetNoDataValue()
band.Fill(nodata)
gdal.RasterizeLayer(raster_dataset, [1], shape_layer, burn_values=[1])
Mi pregunta por tanto es: Es esto un error en GDAL? O es RasterizeLayer la grabación de datos, basado en algo más que sólo la presencia o ausencia de un polígono dentro de una determinada área de píxeles?
Edit: Los archivos que yo estaba usando se puede encontrar aquí.
