Ich verwende ArcGIS 10.1 und möchte ein neues Raster erstellen, das auf zwei bereits vorhandenen Rastern basiert. Das RasterToNumPyArray hat ein gutes Beispiel, das ich anpassen möchte.
import arcpy
import numpy
myArray = arcpy.RasterToNumPyArray('C:/data/inRaster')
myArraySum = myArray.sum(1)
myArraySum.shape = (myArray.shape[0],1)
myArrayPerc = (myArray * 1.0)/ myArraySum
newRaster = arcpy.NumPyArrayToRaster(myArrayPerc)
newRaster.save("C:/output/fgdb.gdb/PercentRaster")
Problem ist, dass es den Raumbezug und auch die Zellengröße entfernt. Ich dachte, es muss arcpy.env sein, aber wie stelle ich sie basierend auf dem Eingabe-Raster ein? Ich kann es nicht herausfinden.
Nach Lukes Antwort ist dies meine vorläufige Lösung.
Beide von Lukes Lösung eingestellten Raumbezug, Ausdehnung und Zellgröße korrekt. Bei der ersten Methode wurden die Daten im Array nicht korrekt übertragen, und das Ausgabe-Raster ist überall mit Knotendaten gefüllt. Seine zweite Methode funktioniert meistens, aber wo ich große Region der Knoten habe, füllt es mit blocky Nullen und 255s. Dies kann damit zusammenhängen, wie ich mit Knotenzellen umgegangen bin, und ich bin mir nicht ganz sicher, wie ich es gemacht habe (sollte aber ein anderes Q sein). Ich habe Bilder von dem beigefügt, wovon ich spreche.
#Setting the raster properties directly
import arcpy
import numpy
inRaster0='C:/workspace/test0.tif'
inRaster1='C:/workspace/test1.tif'
outRaster='C:/workspace/test2.tif'
dsc=arcpy.Describe(inRaster0)
sr=dsc.SpatialReference
ext=dsc.Extent
ll=arcpy.Point(ext.XMin,ext.YMin)
# sorry that i modify calculation from my original Q.
# This is what I really wanted to do, taking two uint8 rasters, calculate
# the ratio, express the results as percentage and then save it as uint8 raster.
tmp = [ np.ma.masked_greater(arcpy.RasterToNumPyArray(_), 100) for _ in inRaster0, inRaster1]
tmp = [ np.ma.masked_array(_, dtype=np.float32) for _ in tmp]
tmp = ((tmp[1] ) / tmp[0] ) * 100
tmp = np.ma.array(tmp, dtype=np.uint8)
# i actually am not sure how to properly carry the nodata back to raster...
# but that's another Q
tmp = np.ma.filled(tmp, 255)
# without this, nodata cell may be filled with zero or 255?
arcpy.env.outCoordinateSystem = sr
newRaster = arcpy.NumPyArrayToRaster(myArrayPerc,ll,dsc.meanCellWidth,dsc.meanCellHeight)
newRaster.save(outRaster)
Bild zeigt Ergebnisse. In beiden Fällen sind die Knotenzellen gelb dargestellt.
Lukes zweite Methode
Meine vorläufige Methode
quelle
Ich hatte einige Probleme, ArcGIS dazu zu bringen, NoData-Werte mit den hier gezeigten Beispielen korrekt zu verarbeiten. Ich habe das Beispiel aus dem Blog reomtesensing.io (das den hier gezeigten Lösungen mehr oder weniger ähnlich ist) erweitert, um besser mit NoData umgehen zu können.
Anscheinend gefällt ArcGIS (10.1) der Wert -3.40282347e + 38 als NoData. Also konvertiere ich zwischen numpy NaN und -3.40282347e + 38 hin und her. Der Code ist hier:
quelle