Batch Solving Many Routes

[ShaharLevenson]

Hello all,

I periodically need to calculate the shortest walking distance for over 1.8 million routes (2 points per route).
With such a large amount of routes I cannot use the desktop ArcGIS - loading that many locations isn't realistic.

What would be the best practice? Creating a python script? Using ArcObjects? Something else?

I have been doing the calculations in smaller pieces (100K, more or less) in version 9.x but now, in version 10, I want to improve the process and do it in one "chunk".

I have access to a strong workstation (20 GB RAM, 24 cores). How can I utilize all the processing potential?

I have no need for hierarchy, restrictions or preserving stop order and all I want is the distance - not the actual routes layer.

Any thoughts?

Thanks,
Shahar.

[NaAn]

Hi,

There are several ways to do the batch solving, for example, use ArcObject or Python. A simple way is to create a route layer template which includes all your analysis settings, i.e. no hierarchy, uncheck all the restrictions. You can add a pair of origin and destination points into the route layer, solve and get the cost, and delete stops and the route from the layer afterwards. Repeat the procedure over all OD pairs, you will get the cost for all.

Here is a simple example with python assuming you have a text input file with contents like:
ID,X1,Y1,X2,Y2
1,-117.195533,34.057058,-117.197204,34.055463
2,-117.17157,46.736774,-117.261162,47.657173

Open ArcMap, add a route layer MyRouteLayer, open python window and run the following code:

import arcpy, os, sys
arcpy.env.overwriteOutput = True

try:
    costDic1 = {}
    f = open(r"C:\temp\OriginDestinations.txt", 'r+')
    layer1 = "MyRouteLayer" #an open in-memory layer
     
    for line in f:
        twoPs = line.split(",")
 
        # if the line contains number
        if twoPs[0].isdigit():
            OPoint = arcpy.Point()
            DPoint = arcpy.Point()
            pointGeometryList = []
            
            OPoint.X = twoPs[1]
            OPoint.Y = twoPs[2]
            DPoint.X = twoPs[3]
            DPoint.Y = twoPs[4]
            pointGeometry = arcpy.PointGeometry(OPoint)
            pointGeometryList.append(pointGeometry)
            pointGeometry = arcpy.PointGeometry(DPoint)
            pointGeometryList.append(pointGeometry)    
     myPoints = "Mypoints"
     arcpy.CopyFeatures_management(pointGeometryList, myPoints) 
            arcpy.na.AddLocations(layer1,"Stops",myPoints)
            arcpy.Solve_na(layer1)
            
            # get the cost info from the layer
            rows = arcpy.SearchCursor(lyrName + "\\Routes")
            for row in rows:
                costDic1[twoPs[0]] = row.Total_TravelTime
                
            # delete stops and route from the layer
            arcpy.DeleteRows_management(lyrName + "\\Stops")
            arcpy.DeleteRows_management(lyrName + "\\Routes")
    f.close 
    
    # Output the result to a txt file    
    fOut = open(r"C:\temp\Out.txt", 'w')
    for costKey in costDic1.keys():
        line = str(costKey) + "," + str(costDic1[costKey])
        fOut.writelines(line)
    fOut.close

except:
    print "Exception"

You can modify the code to open the template layer on disk and run it outside ArcMap. Hope this helps!

Anna

[MichaelVorster]

Hi Anna 

I'm hoping you can help. I am new to python coding and have tried to create a script to solve and display multiple routes between OD pairs simultaneously. In my particular case the origins are also the destinations. The script used is as shown below.

# Name: MakeRouteLayer_MultiRouteWorkflow.py
# Description: Calculate the home-work commutes for cyclists and save
# the output to a feature class
# Requirements: Network Analyst Extension

#Import system modules
import arcpy
from arcpy import env

try:
      #Check out the Network Analyst extension license
      arcpy.CheckOutExtension("Network")

      #Set environment settings
      env.overwriteOutput = True

      #Set local variables
      inNetworkDataset = "C:/Users/michael.vorster/Google Drive/Masters/Thesis/Model Info/NMBM Roads/5km       Zone/roads_in_5km_classification/Commuter_Perspective_ND.nd"
      inStops_Home = "C:/Users/michael.vorster/Google Drive/Masters/Thesis/Model Info/Cadastral/ArcGIS/Allotment       Centroids/Allotment_Centroids.shp"
      inStops_Work = "C:/Users/michael.vorster/Google Drive/Masters/Thesis/Model Info/Cadastral/ArcGIS/Allotment       Centroids/Allotment_Centroids.shp"
      outNALayerName = "CycleRoutes"
      outRoutesFC = "C:/Users/michael.vorster/Google Drive/Masters/Thesis/Model Info/Route Data/outRoutes"
      impedanceAttribute = "Commuter"

      #Create a new Route layer. Optimize on Commuter, but compute the
      #distance travelled by accumulating the Length attribute.
      outRouteResultObject = arcpy.na.MakeRouteLayer(inNetworkDataset, outNALayerName, impedanceAttribute,       "FIND_BEST_ORDER", "PRESERVE_BOTH", "", ["Commuter", "Length"], "NO_UTURNS", "", "NO_HIERARCHY",       "", "TRUE_LINES_WITH_MEASURES","")

      #Get the layer object from the result object. The route layer can now be
      #referenced using the layer object.
      outNALayer = outRouteResultObject.getOutput(0)

      #Get the names of all the sublayers within the route layer.
      subLayerNames = arcpy.na.GetNAClassNames(outNALayer)
      #Store the layer names that we will use later
      stopsLayerName = subLayerNames["Stops"]
      routesLayerName = subLayerNames["Routes"]

      #Before loading the commuters' home and work locations as route stops, set
      #up field mapping. Map the "ALLOTMENT" field from the input data to
      #the RouteName property in the Stops sublayer, which ensures that each
      #unique ALLOTMENT will be placed in a separate route. Matching
      #ALLOTMENT from inStops_Home and inStops_Work will end up in the same
      #route.
      fieldMappings = arcpy.na.NAClassFieldMappings(outNALayer, stopsLayerName)
      fieldMappings["RouteName"].mappedFieldName = "ALLOTMENT"

      #Add the commuters' home and work locations as Stops. The same field mapping
      #works for both input feature classes because they both have a field called
      #"ALLOTMENT"
      arcpy.na.AddLocations(outNALayer, stopsLayerName, inStops_Home, fieldMappings, "300 meters", "FID",       [["Commuter_Perspective_ND", "SHAPE"], ["Commuter_Perspective_ND_Junctions", "NONE"]], "", "", "SNAP", "", "",       "")
      arcpy.na.AddLocations(outNALayer, stopsLayerName, inStops_Work, fieldMappings, "300 meters", "FID",       [["Commuter_Perspective_ND", "SHAPE"], ["Commuter_Perspective_ND_Junctions", "NONE"]], "", "", "SNAP", "", "",       "")

      #Solve the route layer.
      arcpy.na.Solve(outNALayer)

      # Get the output Routes sublayer and save it to a feature class
      RoutesSubLayer = arcpy.mapping.ListLayers(outNALayer, routesLayerName)[0]
      arcpy.management.CopyFeatures(RoutesSubLayer, outRoutesFC)

      print "Script completed successfully"

except Exception as e:
      # If an error occurred, print line number and error message
      import traceback, sys
      tb = sys.exc_info()[2]
      print "An error occured on line %i" % tb.tb_lineno
      print str(e)

The script runs through the python window and the route layer is created.  The stops are added to the map but there are no routes displayed. There is also no route information in the attribute table. I suspect this has to do with the fact that the origins and destinations are the same points. I understand that the route info from point "A" to point "A" is zero (0) but there should be info for "A to B" and "A to C", "A to D"....etc. 

Any help with the above would be greatly appreciated. 

Thanks, 

Mike

[JaySandhu]

What are you trying to achieve? Are you trying to get a distance matrix for all points that you have? In that case the route solver is the wrong approach. It will only solve the routes for stop pairs loaded with the same ROUTENAME property, which I think in your case is stop A to stop A and then then Stop B to Stop B, etc. It is not going to solve Stop A to Stop B, unless you load that pair as well. So if you want all pairs, use the OD Cost Matrix and load your locations as origins and destinations and solve. It will solve and give you a complete matrix of travel times. But it only computes a straight line route shape (for speed). And in case you need the actual route geometry of the shortest paths, then use the Closest Facility solver.

Jay Sandhu

[MichaelVorster]

Hi Jay 

Thank you for the response.

I'm actually looking for a combination of the OD Cost Matrix and Closest Facility solvers. The OD Cost Matrix gives me everything I need except the actual routes. The resulting attribute table is perfect and if I could get the actual routes then I would be a very happy man.

I have tried the Closest Facility solver as well but experienced the same problem as with the python code in my original post (i.e. no routes were generated, just points). I suspect this is because the origins and destinations are the same points and are all located on the same layer. I could go through the process of moving each of the points onto a new layer and then solving the routes between a single origin on one layer and multiple destinations on another but this is extremely time consuming, especially if there are lots of points.

Is there any way to use the OD Cost Matrix coding but to change to the output shape from a STRAIGHT_LINE to a TRUE_LINE_WITH_MEASURES?

Regards, 

Mike

[JaySandhu]

Mike,

CF is like OD. One difference is that by default it is set to find ONE closest. You need to change it to find the number of locations you have. It will give you all the routes. Look at the TargetFacilityCount parameter when you are doing the make CF layer.

Closest facility analysis—Help | ArcGIS Desktop 

Jay Sandhu

[MichaelVorster]

Hi Jay 

You are an absolute legend. Thanks for the help! I just ran the CF solver and it seems to be working just fine.

Thanks again.

Regards, 

Mike

[MichaelVorster]

Hi Jay 

I'm hoping you can help me with something else. I am now trying to calculate the straight line or Euclidean distances between the same points. My reasoning for doing this is so that I can get the ratio of straight line distance to route distance for each of the OD pairs.

To get the Euclidean distances between all OD pairs I was going to use the, "Generate Near Table" tool but for some reason it is not available in my ArcToolbox. Are you aware of another method I can use to get the straight line distances? 

Kind regards, 

Mike

[JaySandhu]

Within ArcGIS, the Generate Near Table tool is the best way to compute Euclidean distances.  The other way would be to script it up in python. You can project the point data to the right projection and then run the Add XY Coordinates tool. Then with the X,Y coordinates, you will have to compute the straight line distance and write it out.

Jay Sandhu

[FrancisSenyah]

Hi Jay, I wonder if you can assist with my use case...I'm actually trying to use the closest facility solver to perform a large calculation, (2.5 million origins to 1 destination initially).  I'm running in ArcGIS Pro and using a Network Dataset stored in SQL server.  After 5 hours or so it fails with the empty geometry error message below which I presume is related to issues with the SDE Road Network being used.  Can you advise how to get around this as the network was built for the whole country and inevitably has some invalid geometries.  Perhaps a batch approach is best but I am unsure how to go about this.  Any advice would be appreciated.  Regards, Francis.