Get extent of Georeferenced Rasters in Python and output to polygon shapefile

المشرف العام

Administrator
طاقم الإدارة
I would like to create a shape file containing the extents of each of the rasters in a directory. Is it possible to capture the extent of a raster using Python?

I have tried

extent1=arcgisscripting.Raster.extent('stg1_05.jpg') Runtime error : 'getset_descriptor' object is not callable

and I can't seem to find any help on the module.

I also tried

arcpy.RasterToPolygon_conversion(inRaster, outPolygons, "SIMPLIFY", field)

Runtime error : ERROR 010247: The size of the output shapefile exceeds the 2 GB limit. ERROR 010067: Error in executing grid expression. Failed to execute (RasterToPolygon).

Anyway this will be a polygon of the whole raster file and not just the extents -even if this is generated, I guess we could then run a merge/dissolve but the files created are to big.

Another option I was thinking of was to convert the raster to layer

import arcpy, glob, os, sys from arcpy import env, mapping path = os.getcwd() env.workspace = path RasterType = 'TIF' FileList = glob.glob(path + "\*." + RasterType) print 'Reading files from ' + path os.chdir(path) #print FileList x=0 z=1005 File=FileList[x] LayerWorking=arcpy.mapping.Layer(File) print File LayerExtent=LayerWorking.getExtent() XMAX = LayerExtent.XMax XMIN = LayerExtent.XMin YMAX = LayerExtent.YMax YMIN = LayerExtent.YMin pnt1 = arcpy.Point(XMIN, YMIN) pnt2 = arcpy.Point(XMIN, YMAX) pnt3 = arcpy.Point(XMAX, YMAX) pnt4 = arcpy.Point(XMAX, YMIN) array = arcpy.Array() array.add(pnt1) array.add(pnt2) array.add(pnt3) array.add(pnt4) array.add(pnt1) polygon = arcpy.Polygon(array) ShapeFile = path + "\Polygon_Extent" + "_" + str(z) + ".shp" print ShapeFile print arcpy.GetMessages() arcpy.CreateFeatureclass_management(path, ShapeFile, "POLYGON") print arcpy.GetMessages() arcpy.CopyFeatures_management(polygon, ShapeFile) Gives following output

Reading files from C:\Python26\script_tests\rectify C:\Python26\script_tests\rectify\Stage1_01a.TIF C:\Python26\script_tests\rectify\Polygon_Extent_1005.shp Executing: CopyFeatures in_memory\fB4DC2172_7D02_44B9_B55C_9E71427AE96E C:\Python26\script_tests\rectify\Polygon_Extent_1004.shp # 0 0 0 Start Time: Thu Jul 14 16:11:38 2011 Failed to execute. Parameters are not valid. ERROR 000725: Output Feature Class: Dataset C:\Python26\script_tests\rectify\Polygon_Extent_1004.shp already exists. Failed to execute (CopyFeatures). Failed at Thu Jul 14 16:11:38 2011 (Elapsed Time: 0.00 seconds) Traceback (most recent call last): File "C:\Python26\script_tests\rectify\temp.py", line 36, in arcpy.CreateFeatureclass_management(path, ShapeFile, "POLYGON") File "C:\Program Files\ArcGIS\Desktop10.0\arcpy\arcpy\management.py", line 1539, in CreateFeatureclass raise e ExecuteError: ERROR 999999: Error executing function. Failed to execute (CreateFeatureclass).

using IDLE...Which is wierd as when you run it from the Python window in ArcMap it works fine when the Create/Copy feature commands are run individually.

ShapeFile = "Polygon_Extent" + "_" + str(z) + ".shp" arcpy.CreateFeatureclass_management(path, ShapeFile, "POLYGON")

Result 'C:\Python26\script_tests\rectify\Polygon_Extent_1005.shp'>

arcpy.CopyFeatures_management(polygon, ShapeFile)

Result 'C:\Python26\script_tests\rectify\Polygon_Extent_1005.shp'>

this seems to be a very convoluted process...

UPDATE: Thanks for the fix. I am trying to add the filename to the table and for some reason it inserts a new row into the table and doesn't accept "updateRow(row)"? what am I doing wrong? Also the files don't seem to retain the projection I assign them.

Reading files from C:\Python26\script_tests\rectify C:\Python26\script_tests\rectify\Stage1_01a.TIF Created: Polygon_Extent_1.shp Filled in C:\Python26\script_tests\rectify\Stage1_01a.TIF Traceback (most recent call last): File "C:\Python26\script_tests\rectify\ResterExtent_toSHP.py", line 45, in rows.updateRow(row) File "C:\Program Files\ArcGIS\Desktop10.0\arcpy\arcpy\arcobjects\arcobjects.py", line 102, in updateRow return convertArcObjectToPythonObject(self._arc_object.UpdateRow(*gp_fixargs(args))) RuntimeError: ERROR 999999: Error executing function.

SCRIPT all working now.

import arcpy, glob, os, sys, arcgisscripting from arcpy import env, mapping path = os.getcwd() env.workspace = path env.overwriteOutput = True RasterType = 'TIF' FileList = glob.glob(path + "\*." + RasterType) print 'Reading files from ' + path os.chdir(path) #print FileList geometry_type = "POLYGON" template = "C:\\Python26\\GDA_1994_MGA_Zone_55.shp" has_m = "DISABLED" has_z = "DISABLED" # Creating a spatial reference object spatial_reference = arcpy.SpatialReference("C:\\Python26\\GDA_1994_MGA_Zone_55.prj") x=0 z=x+1 for File in FileList: #File=FileList[x] RasterFile = arcgisscripting.Raster(File) RasterExtent = RasterFile.extent print File XMAX = RasterExtent.XMax XMIN = RasterExtent.XMin YMAX = RasterExtent.YMax YMIN = RasterExtent.YMin pnt1 = arcpy.Point(XMIN, YMIN) pnt2 = arcpy.Point(XMIN, YMAX) pnt3 = arcpy.Point(XMAX, YMAX) pnt4 = arcpy.Point(XMAX, YMIN) array = arcpy.Array() array.add(pnt1) array.add(pnt2) array.add(pnt3) array.add(pnt4) array.add(pnt1) polygon = arcpy.Polygon(array) arcpy.CreateFeatureclass_management(path, "Temp_Polygon_Extent" + "_" + str(z), geometry_type, template, has_m, has_z, spatial_reference) arcpy.CopyFeatures_management(polygon, "Temp_Polygon_Extent" + "_" + str(z)) ShapeFile = "Temp_Polygon_Extent" + "_" + str(z) + ".shp" print "Created: " + ShapeFile arcpy.AddField_management(ShapeFile,'FileName','TEXT') desc = arcpy.Describe(ShapeFile) print desc, ShapeFile #rows = arcpy.InsertCursor(ShapeFile, desc) rows = arcpy.UpdateCursor(ShapeFile) #row = rows.newRow() #row.FileName = str(File) #row.FileName = File print 'Filled in: ', str(File) #rows.insertRow(row) for row in rows: row.FileName = str(ShapeFile) rows.updateRow(row) x=x+1 z=z+1 #cleaning up arcpy.CreateFeatureclass_management(path, "Extent.shp", geometry_type, template, has_m, has_z, spatial_reference) list =[] lstFCs = arcpy.ListFeatureClasses("Temp_Polygon_Extent*") print 'Merging Polygon_Extents* to Extent.shp' for fc in lstFCs: print fc list.append(fc) arcpy.Merge_management(list, "Extent.shp") #print 'Deleting identical entries and temp files' #arcpy.DeleteIdentical_management("Extent.shp", ["SHAPE"]) print 'Created Extent.shp and exiting' for item in list: arcpy.Delete_management(item) del row, rows GDAL Ver.

import os, gdal gdaltindex Extent.shp *.tif
File "". line 1 SyntaxError: invalid syntax

BTW FTTools Python window is a pain as you can't copy/past code/errors from it to.



أكثر...
 
أعلى