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")
Aangemelde leden kunnen berichten plaatsen, updates volgen en meer. Nieuw hier? Registreer een gratis account.
Find useful guides, FAQs, and documents to help you navigate and make the most of Esri Community.
