Ich versuche, ein Shapefile mit GDALs RasterizeLayer in ein Raster zu brennen. Ich erstelle ein Interessenbereich-Raster aus einem anderen Shapefile bei einer bestimmten Pixelgröße. Dieser AOI dient dann als Basis für alle folgenden Rasterisierungen (gleiche Anzahl von Spalten und Zeilen, gleiche Projektion und Geotransformation).
Das Problem tritt jedoch auf, wenn ich die Formen auf der Grundlage derselben Pixelgröße und Projektionen in ein eigenes Raster brenne. Der unten stehende Link (hat nicht genügend Wiederholungen, um das Bild zu veröffentlichen) zeigt das ursprüngliche Shapefile in Tan und das dunkle Rosa, in dem RasterizeLayer Daten gebrannt hat. Das Hellrosa ist ein Knotenwert für die dunkelrosa Rasterdaten. Das Grau ist der AOI, auf dessen Grundlage der Shapefile-Brennvorgang abgeschlossen wurde.
Angesichts der Ausmaße der Shapefile-Polygone würde ich erwarten, dass in den unteren beiden Ecken Brennwerte sowie die beiden Pixel unter den angezeigten Daten angezeigt werden. Offensichtlich ist dies jedoch nicht der Fall.
Wie folgt ist der Code, den ich verwendet habe, um diese zu generieren. Alle Formen wurden mit QGIS erstellt und in derselben Projektion. (Es sollte beachtet werden, dass das Gitter in dem gezeigten Bild nur eine Vorstellung von der Pixelgröße gab, die ich verwendet habe.)
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])
Handelt es sich um einen Fehler in GDAL, oder brennt RasterizeLayer Daten, die nicht nur auf dem Vorhandensein oder Fehlen eines Polygons in einem bestimmten Pixelbereich beruhen?
Die von mir verwendeten Dateien finden Sie hier .
quelle
Antworten:
Ich habe diese Woche mit GDALRasterizeLayers gespielt und habe eine ziemlich gute Vorstellung davon, was es tut. Standardmäßig wird ein Pixel gerastert, wenn sich die Pixelmitte innerhalb des Polygons befindet. Wenn sich nichts in der Mitte befindet, wird es nicht gerastert, auch wenn sich Teile eines Polygons innerhalb der Pixelgrenzen befinden. Verwenden Sie die Option "ALL_TOUCHED", damit das Rastern wie gewünscht funktioniert:
quelle
['ALL_TOUCHED=TRUE']
leider so, dass nur Polygonebenen fixiert werden. Meine Punkt-Shapefile-Ebenen sind immer noch super-wackelig und erscheinen einen Pixel über der Stelle, an der sie platziert sind.