I was trying to get this script done yesterday, but puzzled why half of it doesn't work. Because ESRI is doing site maintenance, I couldn'tget back in to ask this question. Client still waiting.
I was trying to get a script to walk down a tree of LIDAR data, find any rasters and feature classes, and clip them to an area of interest. My test data tree just has some tifs in one folder, some shapefiles in another to make things simpler. And I need to make this run in ArcGIS 10.0, so da.walk and da.Cursor are not available.
When my script hit the first tif file, the AOI polygon didn't overlap it at all and I got an error message that it didn't overlap. (Failed to create raster dataset. The reason could be: The clip feature is outside the raster extent.
Failed to execute (Clip).) I was going to put off solving the problem of testing for overlap until I could get the file finding and clipping figured out. So, following an example at http://gis.stackexchange.com/questions/13571/arcpy-compare-extents-of-rasters, I added the commented out section, lines 26 to 49, to check for overlap of extent for the tifs. However, the script hits line 48, and says name 'Subset' is not defined. Also none of the print statements in the section are executed, so it looks like none of the extent tests worked. What went wrong? If I just run the script down to the first print statement, all the tifs are found and printed.
In the next section, an if statement looks for just shapefiles, and successfully clips them and writes the output to a new output tree. No check for overlap is done for the shapefiles but I probably should to eliminate empty output shapefile and to save execution time.
Hey, how do you insert code in this new forum?
#testclip.py
#test trying to walk a directory tree, find rasters and shapefiles, and clip to a area of Interest (AOI)
# and make it run in 10.0 so without da.walk.
# Paul Huffman, Yakama Fisheries, 7/7/2014
# Import arcpy module
import arcpy
import os
arcpy.CheckOutExtension("spatial")
# Set Geoprocessing environments
arcpy.env.scratchWorkspace = "c:\\avdata\\PythonTest\\ClipLidar\\Scratch"
arcpy.env.workspace = "c:\\avdata\\PythonTest\\ClipLidar\TestData"
arcpy.env.overwiteOutput = True
clipshape = "c:\\avdata\\PythonTest\\ClipLidar\\AOIPoly.shp"
InputFolder = "c:\\avdata\\PythonTest\\ClipLidar\\TestData\\Input"
OutputFolder = "c:\\avdata\\PythonTest\\ClipLidar\\TestData\\Output"
#Find AOI extent
AOI = arcpy.mapping.Layer(clipshape)
AOIextent=AOI.getExtent()
#Find the input feature classes
for root, dirs, files in os.walk(InputFolder):
for name in files:
if os.path.splitext(name)[1] == ".tif":
print os.path.join(root, name), '\n'
FileObj = arcpy.mapping.Layer(os.path.join(root, name))
FileExtent = FileObj.getExtent()
print "Processing: " + name
print FileExtent
## Extent_Overlap = str(FileExtent.overlaps(AOIextent))
## Extent_Within = str(FileExtent.within(AOIextent))
## Extent_Touches = str(FileExtent.touches(AOIextent))
## Extent_Crosses = str(FileExtent.crosses(AOIextent))
## if Extent_Overlap == 'True':
## Subset = 1
## print "Extent_Overlaps"
## elif Extent_Within == 'True':
## Subset = 1
## print "Extent_Within"
## elif Extent_Touches == 'True':
## Subset = 1
## print "Extent_Touches"
## elif Extent_Crosses == 'True':
## Subset = 1
## print "Extent_Crosses"
## #print(os.path.join(root, name))
## #arcpy.Clip_management(name,clipshape,(os.path.join(OutputFolder,name)))
## if Subset == 1:
## arcpy.Clip_management((os.path.join(root, name)),"#",(os.path.join(OutputFolder,name)),clipshape,"#","ClippingGeometry")
## if os.path.splitext(name)[1] == ".shp":
## print(os.path.join(root, name))
## #arcpy.Clip_management((os.path.join(root, name)),"#",(os.path.join(OutputFolder,name)),clipshape,"#","ClippingGeometry")
## arcpy.Clip_analysis((os.path.join(root, name)),clipshape,(os.path.join(OutputFolder,name)))
End of code.
I also found the python glob method. It also can find all the tifs by searching through the tree like:
#Find the input feature classes
print (glob.glob("c:\\avdata\\PythonTest\\ClipLidar\\TestData\\Input\\*\\*.tif"),'\n')
##for root, dirs, files in os.walk(InputFolder):
## for name in files:
## if os.path.splitext(name)[1] == ".tif":
## print os.path.join(root, name), '\n'
where the glob statement returns the same list as the next four commented out statements. However, I eventually want to parse out the folders in which the files were found so I can reconstruct the input tree structure in the output tree. Don't know I can do that with the glob returns. Actually I don't know how to to that yet with the walk returns either.
