sr = arcpy.SpatialReference() sr.factoryCode = 32759 # case sensitive sr.create() # essential step that is non-intuitive! print sr.name # to check
import arcpy from arcpy import env env.overwriteOutput = True prjFile = r"C:\Program Files (x86)\ArcGIS\Desktop10.0\Coordinate Systems\Geographic Coordinate Systems\World\WGS 1984.prj" spatRef = arcpy.SpatialReference(prjFile) # list each raster and extract coordinates env.workspace = r"C:\temp\python\rasters" lstRasters = arcpy.ListRasters("*", "TIF") for ras in lstRasters: rasName = ras.split('.')[0] lat_lon = rasName.split('_')[1] if 'N' in lat_lon: lat = int(lat_lon[1:3]) else: lat = int(lat_lon[1:3]) * -1 if 'W' in lat_lon: lon = int(lat_lon[4:]) * -1 else: lon = int(lat_lon[4:]) # Create a point based on the coordinate values point = arcpy.Point(lon, lat) UTMfc = r"C:\TEMP\Python\UTM.shp" # create an empty point feature class with a field called 'Name' env.workspace = r"C:\temp\python\test.gdb" arcpy.CreateFeatureclass_management(env.workspace, "LatLon", "POINT", "", "", "", spatRef) arcpy.AddField_management("LatLon", "Name", "TEXT") # Populate the feature class with the coordinate value and the raster name rows = arcpy.InsertCursor("LatLon") row = rows.newRow() row.Name = rasName row.Shape = point rows.insertRow(row) del row, rows # Perform an intersect to find which zone the raster falls in arcpy.Intersect_analysis(["LatLon", UTMfc], "Intersect", "", "", "POINT") rows = arcpy.SearchCursor("Intersect") for row in rows: print str(row.Name) + " is in zone UTM" + str(row.Zone) arcpy.Delete_management("LatLon") arcpy.Delete_management("Intersect")
Angemeldete Mitglieder können Beiträge verfassen, Updates folgen und mehr. Neu hier? Registriere ein kostenloses Konto.
Find useful guides, FAQs, and documents to help you navigate and make the most of Esri Community.
