Raster durch mehrere Datensätze oder Polygone schneiden?

Ich möchte ein DEM mit einem Gitter aus Polygonen ausschneiden. Es ist wahrscheinlich einfacher, mehrere Polygone in einer Formdatei zu verwenden, aber ich habe dies nicht geschafft, daher versuche ich, eine for-Schleife zu verwenden, damit ich jedes Dataset in einer GDB durchlaufen kann (jedes enthält nur ein Polygon).

Hier ist mein Code (im Python-Fenster).

#creating a workspace and a list of feature classes
arcpy.env.workspace = "C:/data/lidar/lidar.gdb"
fcs = arcpy.ListFeatureClasses()

#looping through each feature class and creating a raster based on the extent of   
#feature class

for fc in fcs:
    arcpy.Clip_management("perth", "#", "C:/data/lidar", fc, "", "ClippingGeometry")

Mein Code wird jedoch nicht ausgeführt, er sitzt nur da und wartet auf etwas anderes ... aber was? Ich kann es für einen Clip zum Laufen bringen, aber nicht mit der Schleife.

Ich bin sicher, ich sollte etwas anderes für die Ausgabe tun, um jedes neue Raster nach Feature-Class oder so zu benennen ... aber ich weiß auch nicht, wie. Bitte lassen Sie mich wissen, ob ich weitere Informationen hinzufügen sollte.

Eine Sache, die mir auffällt, ist, dass Ihr dritter Parameter eine fest codierte Ausgabe ist (C: / data / lidar). Die Art und Weise, wie es jetzt geschrieben wird, durchläuft jede Ihrer Funktionen und überschreibt die Ausgabe jedes Mal. Da Sie jedoch möglicherweise das automatische Überschreiben von Dateien nicht zugelassen haben, kann dies möglicherweise zum Auflegen führen. Versuchen Sie, für jede Iteration einen eindeutigen Ausgabenamen zu erstellen:

#creating a workspace and a list of feature classes
arcpy.env.workspace = "C:/data/lidar/lidar.gdb"
fcs = arcpy.ListFeatureClasses()

#looping through each feature class and creating a raster based on the extent of   
#feature class

i=0
for fc in fcs:
    outputPath = "C:/data/lidar" + str(i)
    i+=1
    arcpy.Clip_management("perth", "#", outputPath, fc, "", "ClippingGeometry")

Ich bin mir auch nicht sicher, ob Sie beabsichtigt haben, die Ausgaben im Ordner C: / data mit dem Namen lidar abzulegen. Beachten Sie, dass der dritte Parameter im Clip der vollständige Pfad Ihres Ausgabe-Rasters ist und nicht ein Ordner, in dem er abgelegt wird Wenn Sie in Ihrem Ausgabepfadnamen keine Erweiterung angeben und diese in einem Standardordner ablegen, handelt es sich um ein Raster. Daher versucht Ihr Programm derzeit, ein neues Raster-Dataset mit dem Namen "lidar" in C: zu erstellen. / Datenordner.

Für zukünftige Suchende: Hier ist eine modifizierte Version des USGS-Raster-Split-Tool- Skripts, für das nichts über der ArcGIS Basic-Lizenzstufe (ArcView) erforderlich ist:

"""
Raster Split Tool 6/16/2011
ArcGIS 10 Script Tool
Python 2.6.5

Contact:
Douglas A. Olsen
Geographer
Upper Midwest Environmental Sciences Center
U.S. Geological Survey
2630 Fanta Reed Road
La Crosse, Wisconsin 54603
Phone: 608.781.6333

#####
modified 2015-04-15 by Jeremiah Poling ([email protected])
Now using only tools available at the ArcGIS Basic license level 
#####

"""
import arcpy
from arcpy import env
import os

# Get Split Shapfile Name
splitShape = arcpy.GetParameterAsText(0)
# Get Field Name
splitField = arcpy.GetParameterAsText(1)
# Get Output Directory
outDirectory = arcpy.GetParameterAsText(2)
# Get Raster to Split
splitRaster = arcpy.GetParameterAsText(3)
# Get type of output tif, img, or GRID
rasterType = arcpy.GetParameterAsText(4)
# Set Workspace environment
env.workspace = outDirectory

# Loop through the rows in the clipping shapefile
cursor = arcpy.SearchCursor(splitShape)
for row in cursor:
    resultName = row.getValue(splitField)

    # Create feature layer of current clipping polygon
    whereClause = '"' + splitField + '" = ' + "'" + resultName + "'"
    arcpy.MakeFeatureLayer_management(splitShape, 'currentMask', whereClause)

    # Replace spaces with underscores
    resultName = resultName.replace(' ','_')

    # Remove .shp suffix
    if resultName[-4:] == '.shp':
        resultName = resultName[:-4]

    arcpy.AddMessage("Processing: " + resultName)

    if rasterType == 'img':
        resultName = resultName + "." + rasterType
    elif rasterType == 'tif':
        resultName = resultName + "." + rasterType
    else:
        print ('No extension')

    # Save the clipped raster
    arcpy.Clip_management(
        in_raster = splitRaster,
        rectangle = "#",
        out_raster = resultName,
        in_template_dataset = 'currentMask',
        nodata_value="255",
        clipping_geometry="ClippingGeometry",
        maintain_clipping_extent="NO_MAINTAIN_EXTENT"
        )

    if arcpy.Exists('currentMask'):
        arcpy.Delete_management('currentMask')

Ein paar Ideen:

1. Probieren Sie das Raster Split Tool aus, das im USGS als Skript- Tool erhältlich ist (siehe beigefügten Quellcode).
2. Verwenden Sie zum einfachen Ausschneiden / Kacheln von Raster das integrierte ArcGIS 10.1-Tool namens Raster teilen. Sie können die Raster entweder nach der Anzahl der Kacheln oder nach der Größe der Kacheln aufteilen.

Quellcode für das USGS Raster Split-Tool:

"""
Raster Split Tool 6/16/2011
ArcGIS 10 Script Tool
Python 2.6.5

Contact:
Douglas A. Olsen
Geographer
Upper Midwest Environmental Sciences Center
U.S. Geological Survey
2630 Fanta Reed Road
La Crosse, Wisconsin 54603
Phone: 608.781.6333
"""
import arcpy
from arcpy import env
import os
from arcpy.sa import *

# Get Shapefile Name
inShape = arcpy.GetParameterAsText(0)
# Get Split Shapfile Name
splitShape = arcpy.GetParameterAsText(1)
# Get Field Name
splitField = arcpy.GetParameterAsText(2)
# Get Output Directory
outDirectory = arcpy.GetParameterAsText(3)
# Get Raster to Split
splitRaster = arcpy.GetParameterAsText(4)
# Get type of output tif, img, or GRID
rasterType = arcpy.GetParameterAsText(5)
# Set Workspace environment
env.workspace = outDirectory

# Run Split command on Input Shapefile
arcpy.Split_analysis(inShape, splitShape, splitField, outDirectory)

# Loop through a list of feature classes in the workspace
for fc in arcpy.ListFeatureClasses():

    # Execute ExtractByMask
    rfc = ExtractByMask(splitRaster, fc)

    # Clean up shp file from work directory
    arcpy.Delete_management(fc, "")    

    arcpy.AddMessage("Processing: " + fc)

    # Replace spaces with underscores
    fc = fc.replace(' ','_')

    # Remove .shp suffix
    fc = fc[:-4]

    # Trim to 13 chars
    fc = fc[0:13]

    if rasterType == 'img':
        fc = fc + "." + rasterType
    elif rasterType == 'tif':
        fc = fc + "." + rasterType
    else:
        print ('No extension')

    # Save the output
    rfc.save(fc)