ArcPy: Clipping Rasters by Shapefile
Clipping rasters using shapefiles is a fundamental GIS operation that allows you to extract portions of raster datasets based on vector boundaries. Using ArcPy, Python’s interface to ArcGIS, you can automate this process efficiently for single rasters or batch operations across multiple datasets.
What is Raster Clipping?
Raster clipping extracts a subset of a raster dataset using the boundary of a vector feature (shapefile, feature class, or geometry). The output raster contains only the pixels that fall within the specified boundary, with areas outside the boundary either removed or set to NoData values.
Prerequisites
- ArcGIS Desktop (ArcMap) or ArcGIS Pro with Spatial Analyst extension
- Python with ArcPy library
- Input raster dataset
- Clipping shapefile or feature class
- Appropriate licensing for Spatial Analyst tools
Basic Syntax
The primary function for clipping rasters in ArcPy is arcpy.management.Clip:
python
arcpy.management.Clip(in_raster, rectangle, out_raster,
in_template_dataset, nodata_value,
clipping_geometry, maintain_clipping_extent)Step-by-Step Implementation
1. Import Required Libraries
import arcpy
from arcpy.sa import *
import os
# Check out Spatial Analyst extension
arcpy.CheckOutExtension("Spatial")
2. Set Environment Settings
# Set workspace
arcpy.env.workspace = r"C:\GIS_Data"
arcpy.env.overwriteOutput = True
# Set coordinate system (optional)
arcpy.env.outputCoordinateSystem = arcpy.SpatialReference(4326) # WGS84
3. Basic Clipping Operation
def clip_raster_by_shapefile(input_raster, clip_shapefile, output_raster):
"""
Clip a raster using a shapefile boundary
Parameters:
input_raster (str): Path to input raster
clip_shapefile (str): Path to clipping shapefile
output_raster (str): Path for output clipped raster
"""
try:
# Perform the clip operation
arcpy.management.Clip(
in_raster=input_raster,
rectangle="#", # Use shapefile extent instead of manual rectangle
out_raster=output_raster,
in_template_dataset=clip_shapefile,
nodata_value="", # Use default NoData value
clipping_geometry="ClippingGeometry", # Clip to exact shapefile boundary
maintain_clipping_extent="NO_MAINTAIN_EXTENT"
)
print(f"Successfully clipped {input_raster} to {output_raster}")
except arcpy.ExecuteError:
print(f"ArcPy Error: {arcpy.GetMessages(2)}")
except Exception as e:
print(f"Python Error: {str(e)}")
# Example usage
input_raster = "elevation.tif"
clip_boundary = "study_area.shp"
output_raster = "elevation_clipped.tif"
clip_raster_by_shapefile(input_raster, clip_boundary, output_raster)
Advanced Examples
Batch Clipping Multiple Rasters
def batch_clip_rasters(raster_folder, clip_shapefile, output_folder):
"""
Clip multiple rasters using the same shapefile
"""
# Set workspace to raster folder
arcpy.env.workspace = raster_folder
# Get list of raster files
raster_list = arcpy.ListRasters("*.tif") # Adjust extension as needed
for raster in raster_list:
try:
# Create output filename
output_name = f"clipped_{raster}"
output_path = os.path.join(output_folder, output_name)
# Clip raster
arcpy.management.Clip(
in_raster=raster,
rectangle="#",
out_raster=output_path,
in_template_dataset=clip_shapefile,
clipping_geometry="ClippingGeometry"
)
print(f"Processed: {raster}")
except Exception as e:
print(f"Error processing {raster}: {str(e)}")
Clipping with Multiple Polygons
def clip_by_multiple_features(input_raster, shapefile, output_folder, id_field):
"""
Clip raster by each feature in a shapefile separately
"""
# Create search cursor to iterate through features
with arcpy.da.SearchCursor(shapefile, [id_field, "SHAPE@"]) as cursor:
for row in cursor:
feature_id = row[0]
geometry = row[1]
try:
# Create temporary feature class for single feature
temp_fc = f"temp_feature_{feature_id}"
arcpy.management.CreateFeatureclass(
out_path=arcpy.env.workspace,
out_name=temp_fc,
geometry_type="POLYGON",
spatial_reference=geometry.spatialReference
)
# Insert the geometry
with arcpy.da.InsertCursor(temp_fc, ["SHAPE@"]) as insert_cursor:
insert_cursor.insertRow([geometry])
# Clip raster
output_raster = os.path.join(output_folder, f"clipped_{feature_id}.tif")
arcpy.management.Clip(
in_raster=input_raster,
rectangle="#",
out_raster=output_raster,
in_template_dataset=temp_fc,
clipping_geometry="ClippingGeometry"
)
# Clean up temporary feature class
arcpy.management.Delete(temp_fc)
print(f"Clipped raster for feature {feature_id}")
except Exception as e:
print(f"Error processing feature {feature_id}: {str(e)}")
Parameter Options Explained
Clipping Geometry Options
- “NONE”: Uses rectangular extent only
- “ClippingGeometry”: Clips to exact shape boundary (recommended)
Maintain Clipping Extent Options
- “MAINTAIN_EXTENT”: Output raster has same extent as input
- “NO_MAINTAIN_EXTENT”: Output raster extent matches clipping boundary
NoData Handling
# Set custom NoData value
arcpy.management.Clip(
in_raster=input_raster,
out_raster=output_raster,
in_template_dataset=clip_shapefile,
nodata_value="-9999", # Custom NoData value
clipping_geometry="ClippingGeometry"
)
Best Practices
1. Coordinate System Alignment
# Check coordinate systems before clipping
input_sr = arcpy.Describe(input_raster).spatialReference
clip_sr = arcpy.Describe(clip_shapefile).spatialReference
if input_sr.name != clip_sr.name:
print("Warning: Coordinate systems don't match")
# Consider reprojecting one dataset
2. Memory and Performance Optimization
# Set processing extent to reduce memory usage
arcpy.env.extent = clip_shapefile
arcpy.env.snapRaster = input_raster
# Use compression for large outputs
arcpy.env.compression = "JPEG 90" # For imagery
arcpy.env.compression = "LZ77" # For other data types
3. Error Handling and Validation
def robust_clip_operation(input_raster, clip_shapefile, output_raster):
"""
Robust clipping with comprehensive error handling
"""
# Validate inputs
if not arcpy.Exists(input_raster):
raise ValueError(f"Input raster does not exist: {input_raster}")
if not arcpy.Exists(clip_shapefile):
raise ValueError(f"Clipping shapefile does not exist: {clip_shapefile}")
# Check geometry type
desc = arcpy.Describe(clip_shapefile)
if desc.shapeType.lower() != "polygon":
raise ValueError("Clipping feature must be polygon geometry")
# Perform clip with error handling
try:
arcpy.management.Clip(
in_raster=input_raster,
rectangle="#",
out_raster=output_raster,
in_template_dataset=clip_shapefile,
clipping_geometry="ClippingGeometry"
)
# Validate output
if arcpy.Exists(output_raster):
print(f"Clip operation successful: {output_raster}")
return True
else:
print("Clip operation failed - no output created")
return False
except arcpy.ExecuteError as e:
print(f"ArcPy execution error: {e}")
return False
Common Issues and Solutions
Issue 1: Memory Errors with Large Rasters
Solution: Process in tiles or use pyramid levels
# Build pyramids before clipping
arcpy.management.BuildPyramids(input_raster)
# Set tile size for processing
arcpy.env.tileSize = "1024 1024"
Issue 2: Coordinate System Mismatches
Solution: Reproject data to common coordinate system
# Project shapefile to match raster
projected_shapefile = "shapefile_projected.shp"
arcpy.management.Project(
in_dataset=clip_shapefile,
out_dataset=projected_shapefile,
out_coor_system=arcpy.Describe(input_raster).spatialReference
)
Issue 3: Empty Output Rasters
Solution: Check spatial overlap and geometry validity
# Check if shapefile overlaps with raster
raster_extent = arcpy.Describe(input_raster).extent
shapefile_extent = arcpy.Describe(clip_shapefile).extent
if not raster_extent.overlaps(shapefile_extent):
print("Warning: No spatial overlap between raster and shapefile")
Performance Considerations
- Use appropriate cell sizes: Match output cell size to analysis requirements
- Consider data types: Use appropriate bit depth (8-bit, 16-bit, 32-bit)
- Implement parallel processing: For batch operations, consider multiprocessing
- Use local storage: Process data on local drives rather than network locations
Clipping rasters by shapefiles using ArcPy is a powerful technique for spatial data processing. By following the examples and best practices outlined above, you can implement efficient, robust clipping operations for both individual rasters and batch processing workflows. Remember to always validate inputs, handle errors gracefully, and optimize for your specific use case and data volumes.
