I want to batch extract raster(.tif) covered by a vector file(.shp). how to operate in python？ Could you give an example?
closed as off-topic by Vince, John Powell, PolyGeo♦ Dec 15 '16 at 4:04
This question appears to be off-topic. The users who voted to close gave this specific reason:
- "Questions seeking help to debug/write/improve code must include the desired behavior, a specific problem or error and the shortest code necessary to reproduce it in the question itself. Providing a clear problem statement and evidence of a code attempt will help others to help you. See: How to create a Minimal, Complete, and Verifiable example." – John Powell, PolyGeo
If you are using Linux/Mac, you could try the following. It makes use of glob (filename matching) to search for '.tif' files in the current directory. gdalwarp is used to subset the rasters.
import glob import subprocess import os ImageList = sorted(glob.glob('*.tif')) Shapefile = 'Name_of_Vector.shp' TargetRes = '10 10' # Defines the target resolution # Create output directory OutDir = './Clipped_Rasters/' if not os.path.exists(OutDir): os.makedirs(OutDir) for Image in ImageList: print('Processing ' + Image) OutImage = OutDir + Image.replace('.tif', '_Clipped.tif') # Defines Output Image # Clip image subprocess.call('gdalwarp -of GTiff -co "COMPRESS=DEFLATE" -tr '+TargetRes+ ' -cutline ' +Shapefile+ ' -crop_to_cutline ' +Image+ ' ' +OutImage, shell=True) # Build image overviews subprocess.call('gdaladdo --config COMPRESS_OVERVIEW DEFLATE ' +OutImage+ ' 2 4 8 16 32 64', shell=True) print('Done.' + '\n') print('All images processed.')
There are other useful options for gdalwarp. For example, nearest neighbour resampling is used by default, but you can change this using the '-r' option. You can also specify the nodata value using the '-dstnodata' command.