Batch Raster Clipping

[Leo_KrisPalao]

Hi ArcGIS Users,

Good Day!

I have been searching for answers in the web but haven't find any plausible answers yet. What I want to do is automate the clipping of my raster files using a set of shapefiles. Basically, what I want to do is:

raster 1 clip with shapefile 1
raster 2 clip with shapefile 2
raster 3 clip with shapefile 3
.
.
.
.
raster 37 clip with shapefile 37

My raster files and shapefiles are in different directory.

My initial code is like the one below which I got from exporting the model builder to a python script.

How can I automate my clipping process in such a way that it will loop to all my raster data and match each raster data with a specified shapefile. I am currently using ArcGIS 10.2. Can you help me build the code using the initial code that I have below?

Thanks for any help.

-Leo

# Import arcpy module
import arcpy


# Local variables:
Shapefile = r'E:\folder\folder\folder\folder\shapefile1.shp'
Raster_Dataset = r'E:\folder\folder\folder\raster1.tif'
Output_Raster_Dataset = r'E:\folder\test.tif'

# Process: Clip
arcpy.Clip_management(Raster_Dataset, "", Output_Raster_Dataset, Shapefile, "", "ClippingGeometry", "NO_MAINTAIN_EXTENT")

[DuncanHornby]

This can all be done in model builder no need to create a python script. You need to read the help file about iterators and inline substitution. Read those pages to understand the model below which will do what you want.

[ATTACH=CONFIG]33864[/ATTACH]

[CharlottePeters95540]

Even though your suggestion was posted 10 years ago, this worked great for me! I was pulling my hair out trying to get a sub-model with an iterator to work within a model that also had an iterator, but your suggestion was so much easier to implement. Brilliant workaround.

[Leo_KrisPalao]

If I use the model that you suggested, can I input multiple raster layers (*.tif)? and when I run the model, it will automatically match the first shapefile in my folders with my first raster?

I have tried something similar before, but what it does is the first shapefile in my folder reiterates through all of my rasters in my Raster List ( i have used raster layer and enable "List of values" option) like:

SHP 1 clip with raster 1
SHP 1 clip with raster 2
SHP 1 clip with raster 3
.
.
.
SHP 1 clip with raster 37

what I would love to do is to match each SHP with my 1st raster in my list, then 2nd shapefile with my 2nd raster, and so on.

But I will try what I can do in the model builder as you suggested.

Thanks,
-Leo

[DuncanHornby]

This is why I parse the output of the iterator to get the shapefile name so c:\temp\fred.shp will parse to "fred". I then use inline substitution to grab the equivalent raster from c:\raster\fred.tif.

[Anonymous User]

If your file names really are named like "raster1" and "shapefile1" and there is a 1:1 relationship, this code should work:

import arcpy, os

raster_ws = r'E:\path_to_raster_fold'
shp_ws = r'E:\path_to_shp_fold'

# output workspace
out_ws = r'E:\path_to_output_fold'

# get lists
arcpy.env.workspace = raster_ws
rasters = [os.path.join(raster_ws, r) for r in arcpy.ListRasters()]

arcpy.env.workspace = shp_ws
shapes = [os.path.join(shp_ws, s) for s in arcpy.ListFeatureClasses()]

# get dictionary {raster1 : shapefile1.shp, ...}
clip_dict = dict(zip(rasters, shapes))

# clip all reasters
for raster, shapefile in clip_dict.iteritems():
    out_raster = os.path.join(out_ws, os.path.basename(raster))
    arcpy.Clip_management(raster, "", out_raster, shapefile, "", "ClippingGeometry", "NO_MAINTAIN_EXTENT")
    print 'Clipped {0}'.format(os.path.basename(raster))

[Leo_KrisPalao]

Hi Caleb,

Thanks for your reply and sharing the code. I would like to modify the code to incorporate calculation of cell statistics after the clipping process. Would you suggest incorporating it in the code that you gave me, or it will be much simpler to write another code just for calculating cell statistics.

Thanks,
-Leo

[Leo_KrisPalao]

Hi Caleb and other ArcGIS Python users,

I modified the code a bit to incorporate computation of cell statistics in the code. Basically what I have done so far is to create a container >>> outputlist = [] that will store the output of the clipping. Then used this list to compute for cell statistics. I think I got the code correctly until Line 28 >>> outputList.append(out_raster) but I have a problem in computing the cell statistics. The error says:

Traceback (most recent call last):
  File "C:\Python27\ArcGIS10.2\Lib\site-packages\pythonwin\pywin\framework\scriptutils.py", line 326, in RunScript
    exec codeObject in __main__.__dict__
  File "E:\Python_scripts\Clip\rasterClip_with_cellstats.py", line 29, in <module>
    for raster in outputList():
TypeError: 'list' object is not callable

The code with a little modification is below:

import arcpy, os, sys
from arcpy import env
from arcpy.sa import *
arcpy.env.overwriteOutput = True

raster_ws = r'E:\Test_awd\Criteria_SP_2mm_v2'
shp_ws = r'E:\Test_awd\Season1'

# output workspace
out_ws = r'E:\Test_awd'

# get lists
arcpy.env.workspace = raster_ws
rasters = [os.path.join(raster_ws, r) for r in arcpy.ListRasters()]

arcpy.env.workspace = shp_ws
shapes = [os.path.join(shp_ws, s) for s in arcpy.ListFeatureClasses()]
         
# get dictionary {raster1 : shapefile1.shp, ...}
clip_dict = dict(zip(rasters, shapes))

# clip all rasters
outputList = []
for raster, shapefile in clip_dict.iteritems():
    out_raster = os.path.join(out_ws, os.path.basename(raster))
    arcpy.Clip_management(raster, "", out_raster, shapefile, "", "ClippingGeometry", "NO_MAINTAIN_EXTENT")
    print 'Clipped {0}'.format(os.path.basename(raster))
    outputList.append(out_raster)
for raster in outputList():
    sumRaster = CellStatistics(["outputList"], "SUM", "DATA")
    sumRaster.save(r'E:\Test_awd')

Can you help me solve the error that I am getting.

Thanks,
-Leo

[Anonymous User]

Leo,

You were close...You do not need the () after outputList. Try modifying the last portion of the code:

import arcpy, os, sys
from arcpy import env
from arcpy.sa import *
arcpy.env.overwriteOutput = True
arcpy.CheckOutExtension('Spatial')

raster_ws = r'E:\Test_awd\Criteria_SP_2mm_v2'
shp_ws = r'E:\Test_awd\Season1'

# output workspace
out_ws = r'E:\Test_awd'

# get lists
arcpy.env.workspace = raster_ws
rasters = [os.path.join(raster_ws, r) for r in arcpy.ListRasters()]

arcpy.env.workspace = shp_ws
shapes = [os.path.join(shp_ws, s) for s in arcpy.ListFeatureClasses()]
         
# get dictionary {raster1 : shapefile1.shp, ...}
clip_dict = dict(zip(rasters, shapes))

# clip all rasters
outputList = []
for raster, shapefile in clip_dict.iteritems():
    out_raster = os.path.join(out_ws, os.path.basename(raster))
    arcpy.Clip_management(raster, "", out_raster, shapefile, "", "ClippingGeometry", "NO_MAINTAIN_EXTENT")
    print 'Clipped {0}'.format(os.path.basename(raster))
    outputList.append(out_raster)

# output list is already a list, so no need to wrap it in [ ]
sumRaster = CellStatistics(outputList, "SUM", "DATA")
sumRaster.save(r'E:\Test_awd')
print "Created Sum Raster"
arcpy.CheckInExtension('Spatial')

[Leo_KrisPalao]

Hi Caleb,

I tried revising the code this morning and run it again, but the process is stuck on sumRaster = CellStatistics(ouputList, "SUM", "DATA"). Is this an indication that we cannot use list of raster values in Cellstatistics function?

Below is the revised code:

import arcpy, os, sys
from arcpy import env
from arcpy.sa import *
arcpy.env.overwriteOutput = True
arcpy.CheckInExtension("Spatial")

raster_ws = r'E:\Test_awd\Criteria_SP_2mm_v2'
shp_ws = r'E:\Test_awd\Season1'

# output workspace
out_ws = r'E:\Test_awd'

# get lists
arcpy.env.workspace = raster_ws
rasters = [os.path.join(raster_ws, r) for r in arcpy.ListRasters()]

arcpy.env.workspace = shp_ws
shapes = [os.path.join(shp_ws, s) for s in arcpy.ListFeatureClasses()]
         
# get dictionary {raster1 : shapefile1.shp, ...}
clip_dict = dict(zip(rasters, shapes))

# clip all rasters
outputList = []
for raster, shapefile in clip_dict.iteritems():
    out_raster = os.path.join(out_ws, os.path.basename(raster))
    arcpy.Clip_management(raster, "", out_raster, shapefile, "", "ClippingGeometry", "NO_MAINTAIN_EXTENT")
    print 'Clipped {0}'.format(os.path.basename(raster))
    outputList.append(out_raster)

for raster in outputList:
    sumRaster = CellStatistics(outputList, "SUM", "DATA")
    sumRaster.save(r'E:\Test_awd\test.tif')
print "Created Sum Raster"

This is the code that I am getting:

Traceback (most recent call last):
  File "C:\Python27\ArcGIS10.2\Lib\site-packages\pythonwin\pywin\framework\scriptutils.py", line 326, in RunScript
    exec codeObject in __main__.__dict__
  File "E:\Python_scripts\Clip\rasterClip_with_cellstats.py", line 32, in <module>
    sumRaster = CellStatistics(outputList, "SUM", "DATA")
  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\sa\Functions.py", line 2931, in CellStatistics
    ignore_nodata)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\sa\Utils.py", line 47, in swapper
    result = wrapper(*args, **kwargs)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\sa\Functions.py", line 2927, in Wrapper
    [function] + Utils.flattenLists(in_rasters_or_constants))
RuntimeError: ERROR 000824: The tool is not licensed.

I tried to use cell statistics in a separate code, and it works. Hence, my question on using list of rasters as input in Cellstatistics. The code looks like this:

import arcpy
from arcpy import env
from arcpy.sa import *

env.workspace = r'E:\Test_awd\output'

Raster01 = "RE_suitable_dek001.tif"
Raster02 = "RE_suitable_dek002.tif"
Raster03 = "RE_suitable_dek003.tif"
Raster04 = "RE_suitable_dek004.tif"
Raster05 = "RE_suitable_dek005.tif"

arcpy.CheckOutExtension("Spatial")

outCellStatistics = CellStatistics([Raster01, Raster02, Raster03, Raster04, Raster05], "SUM", "DATA")
outCellStatistics.save(r'E:\Test_awd\output.tif')

I think the code is working okay, but get stuck with the line >>> sumRaster = CellStatistics(outputList, "SUM", "DATA")

Thanks,
-Leo