This is my first attempt at using arcpy: I have two CSVs which have columns:
CSV1: join|id_a|lat_a|lon_a
CSV2: join|id_b|lat_b|lon_b
Each has got 100 rows of data. I would like to write a script that gives me the drive-time from row[join] in csv1 to row[join] in csv2.
First I convert the CSVs into layer files:
env.workspace = ".../Documents/ArcGIS/Test"
for x in ["a", "b"]:
try:
in_Table = "test_%s.csv" % x
x_coords = "lon_%s" % x
y_coords = "lat_%s" % x
out_Layer = "test_%s_layer" % x
saved_Layer = "test_%s.lyr" % x
spRef = "Coordinate Systems/Geographic Coordinate Systems/World/WGS 1984.prj"
arcpy.MakeXYEventLayer_management(in_Table, x_coords, y_coords, out_Layer, spRef)
print arcpy.GetCount_management(out_Layer)
arcpy.SaveToLayerFile_management(out_Layer, saved_Layer)
except:
# If an error occurred print the message to the screen
print "Something went wrong"
print arcpy.GetMessages()
Then I try to do the route layer:
#http://desktop.arcgis.com/en/arcmap/10.3/tools/network-analyst-toolbox/make-route-layer.htm
import arcpy
import time
from arcpy import env
try:
# Check out the Network Analyst extension license
arcpy.CheckOutExtension("Network")
# Set environment settings
env.workspace = ".../ArcGIS/Default.gdb"
env.overwriteOutput = True
# Set local variables
inNetworkDataset = ".../ArcGIS/osm.gdb/france_osm/france_osm_nd"
outNALayerName = "BestRoute"
impedanceAttribute = "DriveTime"
outLayerFile = ".../ArcGIS/Test/outRoutes"
# Locations(should have common field: join)
route_a = ".../ArcGIS/Test/test_a.lyr"
route_b = ".../ArcGIS/Test/test_b.lyr"
# My original excel had columns: join|id_a|lat_a|lon_a
# Create a new Route layer
outNALayer = arcpy.na.MakeRouteLayer(inNetworkDataset, outNALayerName, impedanceAttribute)
# Get the layer object from the result object
outNALayer = outNALayer.getOutput(0)
# Get the names of all the sub-layers within the route layer.
subLayerNames = arcpy.na.GetNAClassNames(outNALayer)
# Stores the layer names that we will use later
stopsLayerName = subLayerNames["Stops"]
routesLayerName = subLayerNames["Routes"]
fieldMappings = arcpy.na.NAClassFieldMappings(outNALayer, stopsLayerName)
fieldMappings["RouteName"].mappedFieldName = "join"
# Add my FROM location
arcpy.na.AddLocations(outNALayer, stopsLayerName, route_a,
fieldMappings, "",
exclude_restricted_elements="EXCLUDE")
# Add my TO location
arcpy.na.AddLocations(outNALayer, stopsLayerName, route_b,
fieldMappings, "",
exclude_restricted_elements="EXCLUDE")
# Solve
print "Solving..."
stime = time.time()
arcpy.na.Solve(outNALayer, "SKIP")
print "Finished in %.0f seconds" % (time.time() - stime)
# Save the solved route layer as a layer file on disk with relative paths
RoutesSubLayer = arcpy.mapping.ListLayers(outNALayer, routesLayerName)[0]
arcpy.management.CopyFeatures(RoutesSubLayer, outLayerFile)
print "Script completed successfully"
except Exception as e:
import traceback, sys
tb = sys.exc_info()[2]
print "An error occured on line %i" % tb.tb_lineno
print str(e)
However, this is extremely slow (and only for 100 locations?) taking more than 10 minutes. Is this because the OSM network database is broken? Or just unoptimised.
Edit: 100 rows ended up taking 428 seconds!!
It works but super slow. I hope my network dataset is not corrupt as it took a week to get it from OSM into ArcCatalog
And my network dataset looks like this:
Can I check it's health somehow?
So this solves the speed issue but leaves me even more confused. Instead of crossing all my stores and imposing a cutoff of 50km with scipy.spatial.cKDTree to generate around 20mill O-D pairings within 50km -> proxy for being within 30 minutes, which I then wanted to solve with the route layer (as above) I decided to just try running an OD cost matrix with a cut-off of 30 minutes.
This went through all 36,000 origins (and = destinations) and generated 5.7 million drive-times in just 55 minutes!
So (36253 origins * 36253 destinations), 30min cut-off, ND = France OSM:
stime = time.time()
arcpy.CheckOutExtension("Network")
env.workspace = wpath
env.overwriteOutput = True
# Set local variables
inNetworkDataset = "france_osm/france_osm_nd"
outNALayerName = "ODCost_30min_36k"
impedanceAttribute = "DriveTime"
# Search tolerance: limits the distance from the location point
# to network lines for network analyst to consider them valid
searchTolerance = "1000 Meters"
inBoth = "test_b.lyr"
outLayerFile = outNALayerName + ".lyr"
#Create a new OD Cost matrix layer
print "Creating OD layer"
outNALayer = arcpy.na.MakeODCostMatrixLayer(inNetworkDataset, outNALayerName,
impedanceAttribute, 30, "",
output_path_shape="NO_LINES")
#Get the layer object from the result object
outNALayer = outNALayer.getOutput(0)
#Get the names of all the sublayers within the OD cost matrix layer.
subLayerNames = arcpy.na.GetNAClassNames(outNALayer)
#Stores the layer names that we will use later
originsLayerName = subLayerNames["Origins"]
destinationsLayerName = subLayerNames["Destinations"]
#Load the warehouse locations as origins using a default field mappings and
#a search tolerance of 1000 Meters.
print "Loading Origin"
arcpy.na.AddLocations(outNALayer, originsLayerName, inBoth, "",
searchTolerance)
#Load the store locations as destinations
fieldMappings = arcpy.na.NAClassFieldMappings(outNALayer, destinationsLayerName)
fieldMappings["Name"].mappedFieldName = "codeinsee"
print "Loading Destination"
arcpy.na.AddLocations(outNALayer, destinationsLayerName, inBoth,
fieldMappings, searchTolerance)
#Solve the OD cost matrix layer
print "Solving..."
arcpy.na.Solve(outNALayer)
print "Finished in %.0f" % (time.time() - stime)
#Save the solved OD cost matrix layer as a layer file on disk with relative
#paths
arcpy.management.SaveToLayerFile(outNALayer,outLayerFile,"RELATIVE")
print "Script completed successfully
How is possible to go from 4 seconds for a distance calculation to nearly 2k a second? Is it because the raw C code vectorises the calculation so 1 origin to 1 destination roughly takes the same as 1 origin to 10,000 destination? Is the OD cost matrix just faster than a route layer?
Maybe since no cut-off was imposed for the route layer it was harmed by small breaks in the network (as it would have to solve a maze through all of France), whereas here it just cuts off after 30minutes. Also I specified not to generate lines which I forgot with the route layer (I guess that speeds things up too).
Can anyone give me some average times for OD cost matrix and route layer so that I can orientate myself (e.g. seconds, or minutes, or hours ...)?
arcpy.CopyFeatures_management(saved_Layer,out_shape)and then using the .shp file, however it hasn't really made any difference – mptevsion Feb 12 '16 at 17:08