Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

another newbie question from me.

I would like to clip a DEM using a grid of polygons. It's probably easier to use multiple polygons in one shape file, but I haven't managed this so I'm trying to use a for loop so I can loop through each dataset in a gdb (each contains only one polygon).

Here's my code (doing it in the python window).

#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")

My code doesn't execute however, it just sits there, waiting for something else... but what? I can get it to work for one clip, but not with the loop.

I'm sure I should be doing something else for the output, to name each new raster by feature class or something... but again, don't know how. Please let me know if I should add any more info. I'm new at this scripting forum thing. Ta.

share|improve this question
2  
Perhaps test first whether your loop is working by commenting out the Clip and placing a print or AddMessage to see if it spits out the name of each feature class OK. –  PolyGeo Oct 23 '12 at 7:46
    
Could you explain why you need to do this clipping? Unless the output is needed for communication with other software, usually there are far more efficient methods to analyze raster data polygon-by-polygon, so perhaps the best answer to your question is to suggest a different approach altogether. –  whuber Oct 23 '12 at 15:38
    
Apologies for posting and running - thanks everyone for your help and I will be back to follow up on this, as soon as I can. –  Rosie Nov 13 '12 at 6:16
    
wuber, I wanted to do this so I could slice our DEMs into manageable chunks to then create contours and then stitch back together. Probably not the most efficient way of doing it. –  Rosie Dec 4 '12 at 6:30
1  
What software are using? ArcGIS, QGIS, etc? –  Chad Cooper Feb 4 '13 at 13:33
add comment

2 Answers

One thing I notice is that your third parameter is a hard coded output (C:/data/lidar). The way it's written now will loop through each of your features and overwrite the output every time, but since you may not have allowed the automatic overwriting of files, this could potentially be the hang-up. Try creating a unique output name for each iteration:

#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")

Also, I'm not certain that you intended to place the outputs in the C:/data folder named lidar... note that the third parameter in clip is the full path of your output raster, not a folder that it will be placed in. If you don't specify an extension in your output path name and place them in a standard folder, it will be a grid, so right now your program is attempting to create a new grid dataset named 'lidar' in the C:/data folder.

share|improve this answer
    
Many thanks mbenedetti, it was the lack of a unique name for each new clip that was one of my problems. I also had something strange going on with it being a gdb - when I tried this using shapefiles in a folder instead, it worked fine. So I probably didn't have right file path (not used to working with gdb). –  Rosie Dec 4 '12 at 6:32
add comment

A few ideas:

  1. Try the Raster Split Tool, available as a script tool from the USGS (see attached source code).
  2. For simple raster clipping/tiling, use the built-in ArcGIS 10.1 tool called Split Raster. You can split the rasters by either the number of tiles or by the size of tile.

enter image description here

Source code for 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)
share|improve this answer
    
Is there a file size limit for the Raster Split Tool? My raster is 20.8 GB! I keep getting "failed to execute". All files have the same projection. Running ArcMap 10.1 –  SoilSciGuy Aug 16 '13 at 19:51
    
You can always run it in 64 bit background geoprocessing mode: resources.arcgis.com/en/help/main/10.1/index.html#//… –  Aaron Aug 17 '13 at 2:08
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.