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

[Anonymous User]

Leo,

When you say "stuck on the the sum raster process" what do you mean? If it just looks unresponsive that is probably because you are running this analysis on 37 rasters, so this could take a long time depending on the resolution of your rasters.

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?

I do not think this is the case...From the tool" rel="nofollow" target="_blank">http://resources.arcgis.com/en/help/main/10.1/index.html#//009z... help documents:

[in_raster_or_constant,...] A list of input rasters for which a statistic will be calculated for each cell within the Analysis window.

This code:

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

Should be correct because it is passing in a list of rasters for the sum analysis.

This part of your code is not correct:

# you do not want to have this in a loop, remove the for loop
for raster in outputList: # take this out
    sumRaster = CellStatistics(outputList, "SUM", "DATA")
    sumRaster.save(r'E:\Test_awd\test.tif')

because you are running the the cell statistics on the entire list of rasters for every iteration of your list (which would be 37 times if there are still 37 rasters in your list). You are also overwriting the test.tif every time as well.

Also, judging from your error:

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.

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

The line highlighted in red needs to be changed to:

arcpy.CheckOutExtension("Spatial")

Because you were checking the license in before running the tool, therefore you were now allowing access to spatial analyst tools. The full code I supplied in my last post should do the trick. I would suggest trying that and if it is hanging up on the sum statistics part, try running that as a separate script.

[Leo_KrisPalao]

Hi Caleb,

Thanks for your great help. The problem with the code is the loop in the cell statistics portion and the arcpy.CheckInExtension("Spatial). I changed the code as you suggested and it worked.

I am sharing the final code so other users can use it:

Again, thanks for your help.
-Leo

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

raster_ws = r'C:\raster\folderpath'
shp_ws = r'E:\shapefile\folderpath'

# output workspace
out_ws = r'C:\outputpath'

# 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 = [] # container
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) # aopend all the output to the container

# compute the sum of all rasters with ignore no data option
sumRaster = CellStatistics(outputList, "SUM", "DATA")
sumRaster.save(r'C\outputpath\ofsumrasters.tif') # remove .tif if GRID format is desired
print "Created Sum Raster"

[Leo_KrisPalao]

Hi Caleb,

Sorry again for this request. I am again needed for some assistance in completing the python code that I need.

In the code that we completed before, I am trying to reiterate the resulting raster (one raster only) to multiple shapefiles (provinces).

I added some code below, but I am having an error, i.e.,

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\clip_raster_with_multifeatures_v2.py", line 46, in <module>
    arcpy.Clip_management(sumRaster, "", out_prov, fc, "", "clippingGeometry", "NO_MAINTAIN_EXTENT")
  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\management.py", line 12803, in Clip
    raise e
ExecuteError: Failed to execute. Parameters are not valid.
ERROR 000445: Extension is invalid for the output raster format.
Failed to execute (Clip).

I think the error is something to do with the filename.

Below is the code the modified code that I am trying to run:

import arcpy, os, sys
from arcpy import env
from arcpy.sa import *
arcpy.env.overwriteOutput = True
arcpy.CheckOutExtension("Spatial")
arcpy.env.workspace = r'E:\output_workspace'
arcpy.env.extent = "MAXOF"

raster_ws = r'E:\raster_workspace'
shp_ws = r'E:\shape_workspace'

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

# 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 = [] # container
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) # append all the output to the container

# compute the sum of all rasters with ignore no data option
sumRaster = CellStatistics(outputList, "SUM", "DATA")
sumRaster.save(r'E:\output_workspace\output.tif')
print "Created Sum Raster"

# clip sumRaster with shapefiles
province = r'E:\workspace_of_shapefiles_per_province'

# Process: Iterate Feature Classes
arcpy.env.workspace = province
fcList = [os.path.join(province, p) for p in arcpy.ListFeatureClasses()] # Get list
for fc in fcList:
    out_prov = os.path.join(out_ws, os.path.basename(fc))
    arcpy.Clip_management(sumRaster, "", out_prov, fc, "", "clippingGeometry", "NO_MAINTAIN_EXTENT")

Thanks very much,
-Leo

[Leo_KrisPalao]

Hi Caleb,

I am almost corrected the code. It seems to work properly. But my problem is the naming of the file. I am having trouble removing the .shp extension in my features list. For instance my filename input is out_prov + "clip.tif". When it saved the file it carries the .shp extension name of my feature list.

This is the output file name that I get: province.shpclip.tif

What I like to have is provinceclip.tif

Here is the last 10 lines of the code that I am using:

# clip sumRaster with shapefiles
province = r'E:\path\province_shapefiles'

# Process: Iterate Feature Classes
arcpy.env.workspace = province
fcList = [os.path.join(province, p) for p in arcpy.ListFeatureClasses()] # Get list
for fc in fcList: # clip the sumRaster with all the shapefiles in the folder
    out_prov = os.path.join(out_ws, os.path.basename(fc))
    arcpy.Clip_management(sumRaster, "", out_prov + "clip.tif", fc, "", "clippingGeometry", "NO_MAINTAIN_EXTENT")
    print "scoring maps clipped by province"

Thanks,
-Leo

[Anonymous User]

Leo,

Try changing the line in red:

# clip sumRaster with shapefiles
province = r'E:\path\province_shapefiles'

# Process: Iterate Feature Classes
arcpy.env.workspace = province
fcList = [os.path.join(province, p) for p in arcpy.ListFeatureClasses()] # Get list
for fc in fcList: # clip the sumRaster with all the shapefiles in the folder
    out_prov = os.path.join(out_ws, os.path.splitext(os.path.basename(fc))[0])
    arcpy.Clip_management(sumRaster, "", out_prov + "clip.tif", fc, "", "clippingGeometry", "NO_MAINTAIN_EXTENT")
    print "scoring maps clipped by province"

Also, you should mark this thread as answered to let others know your issue has been resolved. If you do not know how, see this [url=http://resources.arcgis.com/en/help/forums-mvp/]page[/url].

[Leo_KrisPalao]

Hi Caleb,

Sorry for this late reply. Why is it in the thread there is no option to mark the reply as best (check mark). I can only vote your reply but not mark it as best to signify that my query has been resolved.

Thanks,
-Leo