topic Re: Raster ExtractByMask with python in Geoprocessing Questions

Raster ExtractByMask with python

bogdanpalade1 — Tue, 29 Jun 2010 09:39:24 GMT

Hello,

I am trying to clip a raster with a polygon shapefile.The shapefile contains several polygons representing administrative regions. The method is to temporarily export each polygon to a new shapefile and to extract the coresponding piece of raster.

The problem is that the output of

gp.ExtractByMask_sa(raster, shapefile, newRaster)

is always just one band of the original RGB raster. Do you know how to extract all the 3 bands?

The original raster is .img

Thank you

Re: Raster ExtractByMask with python

ChrisSnyder — Tue, 29 Jun 2010 16:30:08 GMT

Check out the Clip_managment tool (not the Clip_analysis one). Be sure to set the snap raster in the Envr. Settings, otherwise the pixel allignment may differ from the original.

Re: Raster ExtractByMask with python

bogdanpalade1 — Wed, 07 Jul 2010 20:05:48 GMT

Hello chris,

Sorry for the late reply. Clip_management asks for a rectangle to do the clip...instead I want to use a shapefile to do the clip with.

Re: Raster ExtractByMask with python

DanPatterson_Retired — Wed, 07 Jul 2010 20:23:36 GMT

http://help.arcgis.com/en/arcgisdesktop/10.0/help/index.html#//009z0000002n000000.htm
indicates that a grid stack is the default output, I am not using 10 so I don't know what other output options there, so check those options

Re: Raster ExtractByMask with python

bogdanpalade1 — Thu, 08 Jul 2010 11:34:40 GMT

Dan,

Thank you for the reply but this works only in ArcGIS 10 environment (by the way, congrats to esri team for that although I still don't understand why only 3 bands can be extracted and not all of them since I guess that not all the people are working with RGB or pan images). I'm using 9.1 and 9.3 at the time so...I have to deal with it for the time being.

In the end I think I will have to do it separatly band by band and run some kind of layer stack afterwards....It will be soooo....slow!

Thanks,
Bogdan

Re: Raster ExtractByMask with python

ChrisSnyder — Thu, 08 Jul 2010 15:47:36 GMT

In v9.3.1 ArcInfo level, the Clip_managment tool does indeed have an option for using a vector featureclass as a psudo-mask (see the "in_template_dataset" parameter). This is tool to use when you have a multi-band image, since the output can be any ESRI-supported format (not just GRID like the ExtractByMask tool). Check it out...

Re: Raster ExtractByMask with python

bogdanpalade1 — Fri, 09 Jul 2010 06:37:06 GMT

Thank you Chris,

I will check it out now...

Re: Raster ExtractByMask with python

bogdanpalade1 — Sat, 11 Dec 2021 10:26:12 GMT

Ok, so the clip_management allows me to use a feature class extent to clip the raster but it's still necessary to introduce the rectangle parameters (I thought that I could find something to substitute the rectangle by a feature class).
Well, the initial purpose was to cut some rasters with a polygon feature class (administrative boundaries). The approach was first to extract single polygons to new feature classes and to use these to cut the rasters located in a specific folder. When I thought all was great I ran into an error: after extracting the single polygons to new shapefiles the clipping process is interupted by the following message:
ERROR 000732: Output Extent: Dataset 1.shp does not exist or is not supported
Failed to execute (Clip).
I attached a print screen of the error message.

The script is the following:

import arcgisscripting, os, string
gp = arcgisscripting.create()
gp.AddToolbox("C:/Program Files/ArcGIS/ArcToolbox/Toolboxes/Data Management Tools.tbx")
gp.SetProduct("ArcInfo")


clip = "D:/hyperspectralImages/torrejon/torrejon_ALI/clip/clip.shp"   #the polygon feature class
rows = gp.SearchCursor(clip)
row = rows.Next()
while row:  # here I am extracting each polygon to an individual feature class
    FID = row.GetValue("FID")
    hoja = row.GetValue("hoja")
    Expression = '\"FID\" = ' + str(FID)
    outFeatures = "D://hyperspectralImages//torrejon//torrejon_ALI//outClips//"+str(hoja)+".shp" #the single fc's will be named after the field "hoja"
    print "Extracting sheet number "+str(hoja)                                                   #like 1.shp, 2.shp...etc
    gp.Select_analysis(clip,outFeatures,Expression)
    row=rows.Next()
print "Finished extracting individual sheets.Begining to clip..."

inputFolder= "D:/hyperspectralImages/torrejon/torrejon_ALI/rasters" #the folder where all the rasters to be clipped are located
gp.Workspace = inputFolder
rasters = gp.ListRasters()
rasters.Reset()
raster = rasters.Next()

clipFeatures = "D://hyperspectralImages//torrejon//torrejon_ALI//outClips//" #the folder where the polygon feature classes used for clipping are located

while raster: # the loop to clip each raster by each polygon feature class
    for feature in os.listdir(clipFeatures):
        if feature.endswith(".shp"):
            outRaster = "D://hyperspectralImages//torrejon//torrejon_ALI//output//"+feature[0:-3]+".tif"
            print "Clipping "+raster+" with sheet number "+feature[0:-3]
            gp.clip_management(raster,"446716,298473 4425996,437405 499153,348855 4492682,68626",outRaster,feature,"255","ClippingGeometry")
                
print "Done"

The 1.shp, 2.shp ..and so on do exist in the right folder...

Do you have any idea of what's happening?

Thanks,
Bogdan

Re: Raster ExtractByMask with python

DanPatterson_Retired — Mon, 19 Jul 2010 05:48:48 GMT

Python paths use a single forward slash /
or double back slash \\
not double forward slash //
as you have in several places

Re: Raster ExtractByMask with python

bogdanpalade1 — Wed, 21 Jul 2010 14:00:43 GMT

Dan,
Thank you for your advice. I changed the paths accordingly but the problem remains the same.

Re: Raster ExtractByMask with python

ChrisSnyder — Wed, 21 Jul 2010 16:00:13 GMT

Try naming your shapefile something else than a number (1.shp). Generally it is a bad idea to name a layer or field with a number as the 1st character.

Re: Raster ExtractByMask with python

bogdanpalade1 — Tue, 03 Aug 2010 08:05:34 GMT

Tried changing the shapefile name but still the same error....

Re: Raster ExtractByMask with python

ChrisSnyder — Tue, 03 Aug 2010 16:18:02 GMT

Not sure what's wrong... Take look at the sample code: http://forums.arcgis.com/threads/9512-Rectangle-for-Clip_management-tool-in-arcpy (my post Aug 2nd @10:25).

Are you sure that your raster and shapefile are using the same projection?

Re: Raster ExtractByMask with python

JamieWilliams — Fri, 05 Nov 2010 10:53:48 GMT

I had the same problem. I solved it by adding the full file path to the feature name.
eg. pathFeature = "R:\\folderContainingShpFile\\" + str(feature)

gp.clip_management(raster,"446716,298473 4425996,437405 499153,348855 4492682,68626",outRaster,pathFeature,"255","ClippingGeometry")