This sample notebook shows how different types of network analysis operations can be performed using the ArcGIS API for Python.

Creating a connection to an ArcGIS Online organization¶

First, you would need to establish a connection to your organization which could be an ArcGIS Online organization or an ArcGIS Enterprise. To be able to run the code in this sample notebook, you would need to provide credentials of a user within an ArcGIS Online organization.

In [1]:

from copy import deepcopy
from datetime import datetime
from IPython.display import HTML
import json
import pandas as pd
from arcgis.gis import GIS
import arcgis.network as network
import arcgis.geocoding as geocoding

user_name = 'arcgis_python'
password = 'P@ssword123'
my_gis = GIS('https://www.arcgis.com', user_name, password)

Creating a network analysis layer¶

Thereafter, you would need to create a NetworkLayer object which is used for operations such as finding the closest facility, the best route between multiple stops, or service area around a facility. To create a NetworkLayer object which you would use later for network analysis, you would need to provide the URL to the network analysis service. However, since all ArcGIS Online organizations already have access to those routing services, you can access this URL through the GIS object. If you have your own ArcGIS Server based map service with network analysis capability enabled, you would need to provide the REST end point URL for this service.

Creating a `RouteLayer` layer¶

In [2]:

route_service_url = my_gis.properties.helperServices.route.url
route_service_url

Out[2]:

'https://route.arcgis.com/arcgis/rest/services/World/Route/NAServer/Route_World'

In [3]:

route_service = network.RouteLayer(route_service_url, gis=my_gis)
route_service

Out[3]:

<RouteLayer url:"https://route.arcgis.com/arcgis/rest/services/World/Route/NAServer/Route_World">

ArcGIS Online Premium Routing services have worldwide coverage which means that you don't have to worry about specifying the analysis region as the service will automatically calculate the region name based on the location of the input points.

In [4]:

my_gis.content.search('title: ArcGIS Online Directions and Routing Services Coverage type:Web Map owner:esri',
                      item_type='Web Map', outside_org=True)[0]

Out[4]:

ArcGIS Online Directions and Routing Services Coverage
This map shows which countries are partially or fully covered by ArcGIS Online's directions and routing services.

Web Map by esri
Last Modified: June 22, 2016
0 comments, 74,157 views

ArcGIS Online also provides access to the traffic map service which you can use to visualize real time traffic speeds and incidents such as accidents, construction sites, or street closures. Traffic visualization displays information about how travel speeds on specific road segments change over time.

In [5]:

my_gis.content.search('title: World Traffic Map type:Web Map owner:esri',
                      item_type='Web Map', outside_org=True)[0]

Out[5]:

World Traffic Map
This map features near real-time traffic information for different regions in the world. The data is updated every 5 minutes. This map requires an ArcGIS Online organizational subscription.

Web Map by esri
Last Modified: June 29, 2017
5 comments, 104,124 views

Calculating the total distance and travel time between two stops using the coordinates pairs, addresses and JSON data¶

To solve the route between multiple stops, you would need to provide the coordinates for the input stops. The set of stops is loaded as network locations during analysis. Stops can be specified using a simple comma / semi-colon based syntax or as a JSON structure.

Depending on what information you need to get returned from the solve() method, you can access the returned JSON object and retrieve the necessary details. Since in this particular case we are not interested in getting the route directions or the route shape, we will set those parameters to False to speed up the process of result generation.

To learn more about the parameters of the RouteLayer.solve() method, please refer to the ArcGIS REST API Help page Route service with synchronous execution.

In this example, we have specified the stops coordinates in the latitude/longitude format. The result object returned from the service is a Python dictionary which makes it very easy to access the needed information about the solved route.

In [7]:

route_layer = network.RouteLayer(route_service_url, gis=my_gis)
result = route_layer.solve(stops='''18.068598,59.329268; 18.068598,59.429268''',
                           return_directions=False, return_routes=True, 
                           output_lines='esriNAOutputLineNone',
                           return_barriers=False, return_polygon_barriers=False, 
                           return_polyline_barriers=False)

travel_time = result['routes']['features'][0]['attributes']['Total_TravelTime']
print("Total travel time is {0:.2f} min".format(travel_time))

Total travel time is 18.91 min

ArcGIS API for Python also provides access to the geocoding functionality available within your ArcGIS Online Organization or ArcGIS Enterprise. So, if your input data doesn't come as pairs of coordinates, you could geocode the input addresses and then pass the coordinates to the solve() method:

In [9]:

stop1_address = '1200 E Colton Ave, Redlands, CA 92374'
stop2_address = '11711 Sand Canyon Rd, Yucaipa, CA 92399'

stop1_geocoded = geocoding.geocode(stop1_address)
stop2_geocoded = geocoding.geocode(stop2_address)

stops = '{0},{1}; {2},{3}'.format(stop1_geocoded[0]['attributes']['X'],
                                  stop1_geocoded[0]['attributes']['Y'],
                                  stop2_geocoded[0]['attributes']['X'],
                                  stop2_geocoded[0]['attributes']['Y'])

route_layer = network.RouteLayer(route_service_url, gis=my_gis)
result = route_layer.solve(stops=stops, return_directions=False, return_routes=True, 
                           output_lines='esriNAOutputLineNone', return_barriers=False, 
                           return_polygon_barriers=False, return_polyline_barriers=False)

travel_time = result['routes']['features'][0]['attributes']['Total_TravelTime']
print("Total travel time is {0:.2f} min".format(travel_time))

Total travel time is 8.05 min

In other cases, your input stops coordinates may be stored using another spatial reference or you may need to pass in additional information about your stops such as names or wait times. In this example, we are passing the input stops in the JSON format with stops' coordinates defined using the Web Mercator coordinate system. This example would be helpful if you store your stops in another system and can programmatically read the stops and parse JSON.

In [10]:

route_layer = network.RouteLayer(route_service_url, gis=my_gis)
input_stops = json.loads('''{"features": [
{"geometry": {"x": -13108894.078499999, "y": 4034551.2082, 
"spatialReference": {"wkid":3857, "latestWkid":3857}},
"attributes": {"Name": "945 Azalea Dr, Pomona, California"}}, 
{"geometry": {"x": -13105731.496399999, "y": 4038487.7151999995, 
"spatialReference": {"wkid":3857, "latestWkid":3857}},
"attributes": {"Name":"1321 Loranne Ave, Pomona, California"}}]}''')

result = route_layer.solve(stops=input_stops, return_directions=False, return_routes=True, 
                           output_lines='esriNAOutputLineNone', return_barriers=False, 
                           return_polygon_barriers=False, return_polyline_barriers=False)

total_distance = result['routes']['features'][0]['attributes']['Total_Kilometers']
print("Total distance is {0:.2f} km".format(total_distance))

Total distance is 5.89 km

Getting the directions between two stops applying the style¶

The solve() operation can also optionally return the directions which you would see when using the Directions widget in the Web App Builder, in the ArcGIS Online map viewer, or when building your own custom web app using the ArcGIS API for JavaScript. The directions object can be loaded into an HTML object to be printed, but it can also be loaded into a pandas data frame and presented as a table for reference purposes.

In this example, we have specified the stops coordinates in the latitude/longitude format.

In [11]:

stops = '''-3.203062,55.906437; -3.190080,55.935570'''

As the ArcGIS Online routing service is traffic aware, one can specify the start time of the departure to get the most accurate travel time. In this example, we will pass as the start time the current time (UTC).

In [13]:

import datetime

start_time = int(datetime.datetime.utcnow().timestamp() * 1000)

route_layer = network.RouteLayer(route_service_url, gis=my_gis)
result = route_layer.solve(stops=stops, 
                           directions_language='en-US', return_routes=False,
                           return_stops=False, return_directions=True,
                           directions_length_units='esriNAUKilometers',
                           return_barriers=False, return_polygon_barriers=False,
                           return_polyline_barriers=False, start_time=start_time,
                           start_time_is_utc=True)

Creating a `pandas` data frame from returned directions¶

After the result object is returned, we would create a table formatted directions panel which would include the information on when the route started, the duration and distance of every segment in the route as well as the direction text.

In [14]:

records = []
travel_time, time_counter = 0, 0
distance, distance_counter = 0, 0

for i in result['directions'][0]['features']:
    time_of_day = datetime.datetime.fromtimestamp(i['attributes']['arriveTimeUTC'] / 1000).strftime('%H:%M:%S')
    time_counter = i['attributes']['time']
    distance_counter = i['attributes']['length']
    travel_time += time_counter
    distance += distance_counter
    records.append( (time_of_day, i['attributes']['text'], 
                     round(travel_time, 2), round(distance, 2))  )

In [15]:

pd.set_option('display.max_colwidth', 100)
df = pd.DataFrame.from_records(records, index=[i for i in range(1, len(records) + 1)], 
                               columns=['Time of day', 'Direction text', 
                                        'Duration (min)', 'Distance (km)'])
HTML(df.to_html(index=False))

Out[15]:

Time of day	Direction text	Duration (min)	Distance (km)
13:38:42	Start at Location 1	0.00	0.00
13:38:42	Go northwest on Buckstone Place toward Buckstone Road	0.59	0.09
13:39:17	Turn left on Buckstone Road	0.77	0.18
13:39:28	Turn right on Buckstone Gardens	1.59	0.44
13:40:17	Turn right on Buckstone Terrace (A702)	1.97	0.60
13:40:40	Make sharp right on Braid Road	4.22	1.98
13:42:55	Turn right on Hermitage Drive	4.82	2.35
13:43:31	Turn left on Braid Avenue	5.62	2.79
13:44:19	Turn right on Cluny Gardens	6.63	3.48
13:45:19	Turn left on Oswald Road	7.72	3.98
13:46:25	Turn left on Kilgraston Road	8.85	4.55
13:47:33	Turn right on Beaufort Road	9.24	4.77
13:47:56	Continue on Grange Road	9.29	4.80
13:47:59	Finish at Location 2, on the left	9.29	4.80

Creating HTML from returned directions¶

If you would need to pass the returned directions to another system which expects to receive an HTML object, you can apply custom formating to the directions object. This is done by using the IPython built-in HTML object and the display() method.

In [16]:

td_data = ['<td align="left">{:s}</td> \
           <td align="left">{:s}</td> \
           <td align="left">{:.2f}</td> \
           <td align="left">{:.2f}</td>'.format(*values) 
           for values in records]
tr_data = ['<tr>{}</tr>'.format(i) for i in td_data]

display(HTML('<table> \
             <tr> <th> Time of day </th> \
             <th>Instruction</th> \
             <th>Time</th> \
             <th>Distance</th> </tr> {0} </table>'.format(''.join(tr_data))))

Time of day	Instruction	Time	Distance
13:38:42	Start at Location 1	0.00	0.00
13:38:42	Go northwest on Buckstone Place toward Buckstone Road	0.59	0.09
13:39:17	Turn left on Buckstone Road	0.77	0.18
13:39:28	Turn right on Buckstone Gardens	1.59	0.44
13:40:17	Turn right on Buckstone Terrace (A702)	1.97	0.60
13:40:40	Make sharp right on Braid Road	4.22	1.98
13:42:55	Turn right on Hermitage Drive	4.82	2.35
13:43:31	Turn left on Braid Avenue	5.62	2.79
13:44:19	Turn right on Cluny Gardens	6.63	3.48
13:45:19	Turn left on Oswald Road	7.72	3.98
13:46:25	Turn left on Kilgraston Road	8.85	4.55
13:47:33	Turn right on Beaufort Road	9.24	4.77
13:47:56	Continue on Grange Road	9.29	4.80
13:47:59	Finish at Location 2, on the left	9.29	4.80

Drawing the result route on a web map as a layer¶

The solve() method can also return the geometry of the route between the stops. The geometry can be added as a graphics to your web map and drawn interactively in the notebook. The solved route as well as the input stops can optionally have a pop-up window configured to include the key information.

In this example, we have again specified the stops coordinates in the latitude/longitude format.

In [17]:

stops = '''-3.203062,55.906437; -3.190080,55.935570'''

route_layer = network.RouteLayer(route_service_url, gis=my_gis)
result = route_layer.solve(stops=stops, 
                           directions_language='en-US', return_routes=True,
                           return_stops=True, return_directions=False,
                           return_barriers=False, return_polygon_barriers=False,
                           return_polyline_barriers=False)

In [19]:

my_map = my_gis.map('Braidburn Valley', zoomlevel=13)
my_map

my_map

In the code below, we would reverse geocode the input stops to enrich the pop-ups of stops with the address information as well as apply custom formatting to the pandas data frame to draw it's contents within the pop-up windows both for the route and the stops.

In [20]:

my_map.clear_graphics()
stop_count = result['routes']['features'][0]['attributes']['StopCount']
travel_time = result['routes']['features'][0]['attributes']['Total_TravelTime']
distance = result['routes']['features'][0]['attributes']['Total_Kilometers']

data = [('Number of stops', stop_count),
        ('Total travel time', '{0:.2f} min'.format(travel_time)),
        ('Total travel distance', '{0:.2f} km'.format(distance))]

df = pd.DataFrame.from_records(data)
styles = [    
    dict(selector="td", props=[("padding", "2px")]),
    dict(selector='.row_heading, .blank', props=[('display', 'none;')]),
    dict(selector='.col_heading, .blank', props=[('display', 'none;')])]

symbol = {
    "type": "esriSLS",
    "style": "esriSLSSolid",
    "color": [128,0,128,90],
    "width": 4
}

popup_route = {"title": "Route", 
               "content": df.style.set_table_styles(styles).render()}
popup_stop = {"title": "Stop {}", 
              "content": df.style.set_table_styles(styles).render()}

my_map.draw(result['routes']['features'][0]['geometry'], popup_route, symbol)

for stop in result['stops']['features']:
    address = geocoding.reverse_geocode(stop['geometry'])['address']['Match_addr']
    my_map.draw(stop['geometry'], 
                {"title": "Stop {}".format(stop['attributes']['Sequence']), 
                 "content": address})
my_map.zoom = 12

Click the stops markers as well as the route polyline to explore their attributes. Because the content of the pop-ups may include HTML source code, it is also possible to have the pop-up windows include other resources such as tables and images.

Solving multiple routes¶

In this example, we will solve multiple routes by reading a .csv file which contains the start/end coordinates for each route. The .csv file looks as below:

Start_X	Start_Y	End_X	End_Y	RouteName
34.032885	-117.875613	33.950965	-118.007080	Route11
34.1032324	-118.0830204	33.961414	-118.215706	Route22
34.015697	-118.359059	33.890443	-118.336239	Route33
34.197766	-118.639790	34.251779	-118.500158	Route44
34.254439	-118.771523	34.188270	-119.032007	Route55

We will read the .csv file into a pandas data frame which we will later enrich with the travel time and distance.

In [ ]:

df = pd.read_csv(r'C:\Data\routes.csv')
df['TravelTime'], df['Distance'] = None, None
df

In [18]:

route_layer = network.RouteLayer(route_service_url, gis=my_gis)

for idx, values in df.iterrows():
    stops = '{0},{1}; {2},{3}'.format(values['Start_Y'], values['Start_X'], 
                                      values['End_Y'], values['End_X'])
    result = route_layer.solve(stops=stops)['routes']['features'][0]['attributes']
    df.loc[idx, 'TravelTime'] = result['Total_TravelTime']
    df.loc[idx, 'Distance'] = result['Total_Kilometers']

df

Out[18]:

	Start_X	Start_Y	End_X	End_Y	RouteName	TravelTime	Distance
0	34.032885	-117.875613	33.950965	-118.007080	Route11	21.7414	21.0819
1	34.103232	-118.083020	33.961414	-118.215706	Route22	31.0056	26.1765
2	34.015697	-118.359059	33.890443	-118.336239	Route33	21.1058	19.8738
3	34.197766	-118.639790	34.251779	-118.500158	Route44	23.5477	20.1751
4	34.254439	-118.771523	34.188270	-119.032007	Route55	30.632	36.9144

Solving the traveling salesperson problem (TSP)¶

In this example, we will solve the traveling salesperson problem (TSP). For this scenario, we will have one single input address which will be used both as the start and the end point. The orders' locations will be loaded from a feature layer hosted within an ArcGIS Online organization. However, the orders location could also come from a .csv file, an ArcGIS Server feature service or from a disk-based data source such a shapefile or a geodatabase feature class.

To be able to solve the TSP, you should specify True for the findBestSequence parameter. The True parameter value causes the service to switch from solving a shortest-path problem to a traveling salesperson problem.

Note Keep in mind that solving a TSP is computer-intensive operation and incurs additional service credits per route.

In [22]:

sample_cities = my_gis.content.search('title:"USA Major Cities" type:Feature Service owner:esri', 
                                      outside_org=True)[0]
sample_cities

Out[22]:

USA Major Cities
This layer presents the locations of cities within United States with populations of 10,000 or greater (based on Census 2000 figures), all state capitals, and the national capital.

Feature Layer Collection by esri
Last Modified: February 07, 2017
2 comments, 1,383,774 views

In [23]:

office_address = '380 New York St, Redlands, CA 92373'

In [24]:

fl = sample_cities.layers[0]
cities_to_visit = fl.query(where="ST = 'CA' AND POP2007 > 300000", 
                           out_fields='NAME', out_sr=4326)

office_coords = geocoding.geocode(office_address)[0]['attributes']
office_feature = {'geometry': {'x': office_coords['X'], 
                               'y': office_coords['Y']}, 
                  'attributes': {'Name': office_address}
                 }

cities_features = deepcopy(cities_to_visit.features)
cities_features.insert(0, office_feature)
cities_features.append(office_feature)

For this particular scenario, we are only interested in finding out the best sequence of visiting the stops; we are not interested in the route directions. However, as we want to draw the lines connecting the order points for visualization purposes only, we need to specify an appropriate value for the outputLines parameter.

In [25]:

route_layer = network.RouteLayer(route_service_url, gis=my_gis)
result = route_layer.solve(stops={'features': cities_features}, preserve_first_stop=True, 
                           preserve_last_stop=True, find_best_sequence=True, return_directions=False, 
                           return_stops=True, return_barriers=False, return_polygon_barriers=False, 
                           return_polyline_barriers=False, return_routes=True, 
                           output_lines='esriNAOutputLineStraight')

In [26]:

my_map = my_gis.map('Santa Maria, CA', zoomlevel=6)
my_map

my_map

In [27]:

my_map.clear_graphics()
my_map.draw(result['routes']['features'][0]['geometry'])

sequence_symbol = {
     "type": "esriTS",
     "color": [255,255,255,255],
     "backgroundColor": [0,0,0,0],
     "borderLineSize": 2,
     "borderLineColor": [255,0,255,255],
     "haloSize": 2,
     "haloColor": [0,0,0,255],
     "verticalAlignment": "bottom",
     "horizontalAlignment": "left",
     "rightToLeft": False,
     "angle": 0,
     "xoffset": 0,
     "yoffset": 0,
     "kerning": True,
     "font": {
      "family": "Segoe UI",
      "size": 9,
      "style": "normal",
      "weight": "bold",
      "decoration": "none"
     },
    "text": ""
}

for i in result['stops']['features'][0:-1]:
    sequence_symbol['text'] = i['attributes']['Sequence']
    my_map.draw(i['geometry'], symbol=sequence_symbol)

sequence_symbol['xoffset'] = 10
sequence_symbol['yoffset'] = 0
sequence_symbol['text'] = ' | {}'.format(len(result['stops']['features']))
my_map.draw(result['stops']['features'][-1]['geometry'], symbol=sequence_symbol)

The solved sequence of stops can also be presented in a table form for reference purposes.

In [28]:

records = []
for i in sorted(result['stops']['features'], key=lambda x: (x['attributes']['Sequence'])):
    records.append((i['attributes']['Name'], i['attributes']['Sequence'] ))
pd.DataFrame.from_records(records, columns=['Location', 'Sequence'])

Out[28]:

	Location	Sequence
0	380 New York St, Redlands, CA 92373	1
1	San Jose	2
2	San Francisco	3
3	Oakland	4
4	Sacramento	5
5	Fresno	6
6	Los Angeles	7
7	Long Beach	8
8	Anaheim	9
9	Santa Ana	10
10	San Diego	11
11	380 New York St, Redlands, CA 92373	12

Feedback on this topic?

ArcGIS API for Python Guide