Trails and Graph Theory 11: Add Connections

We have built tools to find a fairly long circuit in the Gila National Forest, without repeating any trail segments. We can make a few improvements.

  • Add forest roads, or even a cross-country route if necessary, to connect the Aldo Leopold Wilderness trail system to the rest of the Gila.
  • Add forest roads, paved roads, and BLM trails to get near resupply towns: Silver City, Pie Town, Alma/Glenwood, Gila Hot Springs (Doc Campbell’s), and perhaps others.
  • Include trails in the Gila Cliff Dwellings. Since this is a national monument, its trails were not added when we imported Gila National Forest trails.
  • When trails enter short sections of non-USFS land, such as private or state trust, add connections for those areas, since they would not be included in the original Gila import command.
  • Perhaps add part of Bursum Road, NM 159. It is one of my favorite road walks of the Grand Enchantment Trail, and connects several trails together.
  • Add the Continental Divide Wilderness Study Area, south of Pie Town, that contains a good span of the CDT and connects two northern protrusions of the Gila National Forest.
  • Add several brief road-walks that are part of the Continental Divide Trail and Grand Enchantment Trail, that would not be imported if we only add foot-paths.

The OSMnx library does not allow us a direct way to import GPX files of specific routes, since we need to use the node values already defined in OpenStreetMaps. As far as I know, OSMnx does not support importing specific trails or roads by name. We can import a rectangular region into a graph, and then merge it with our existing Gila trail graph.

As an example, Alma and Glenwood are near Mineral Creek Trail (quite an interesting route with a lot of history and geology) but this trail would normally be excluded from our loop because it dead-ends at a trailhead. By connecting Mineral Creek Trail with short road-walks to Alma (good cafe and small convenience store) and Glenwood (post office and maybe convenience store) and then to Little Whitewater Creek Trail (or Whitewater Creek Trail with the Catwalk), we can include these trails in our loop.

To add these roads, which might include other roads we do not care about, we can get away with a single rectangle.

# custom filter imports roads AND foot-paths
cf_add_routes = '["highway"~"path | footway | 
    motorway | motorway_link | trunk | primary | secondary | 
    tertiary | unclassified | residential | living_street | 
    residential | service | track"]'

n,s,e,w = 33.42443, 33.26350, -108.78571, -108.92392
GT = ox.graph_from_bbox(n,s,e,w,simplify=False,truncate_by_edge=True,
    retain_all=True,custom_filter=cf_add_routes)
# merge with existing graph G
G2 = nx.compose(GT,G2)

Whenever we add a connection, we add code for a test to make sure any roads and trails that we add properly connect to each other.

description = 'glenwood mineral creek to little whitewater'
Point1 = (33.42070, -108.82503)
Point2 = (33.32292, -108.81121)

node1 = ox.nearest_nodes(T,Point1[1], Point1[0])
node2 = ox.nearest_nodes(T,Point2[1], Point2[0])
try:
    path_length = nx.shortest_path_length(T,node1,node2,weight='length')
except:
    print('no path found: ', description)
else:
    print(description,' distance(miles): ', str(round(path_length*meters_to_mile,1)))
glenwood mineral creek to little whitewater  distance(miles):  14.7

For Silver City, the CDT connects to the city going northwest on a long road walk, and there is a network of trails about 5 miles to the east at Arenas Valley. Corre Caminos bus line goes from Silver City to Arenas Valley a few times each weekday, so we will include these connecting roads to pass through Silver City.

Several short roads connect between trails in the Aldo Leopold Wilderness and Gila Wilderness. When importing only trails into NetworkX, these paths would be interpreted as disconnected. Fixing this might include more of the Aldo. (However, the Aldo is blocked on the south and east, and to the north has few connected trails going out to the rest of the Gila, aside from the CDT, so the effect may be limited.)

We will not detail all the other route additions here, but they are documented in the source code, available for download below. Pie Town did not seem to have a loop opportunity, so that trail town will be added as an optional supply route.

An option was added to our draw() command to display all the connection rectangles we added. This was essential for debugging.

Finding and adding all the brief road connections for the CDT was very tedious and took many iterations, and I am still not finished. When one is making this many additions by hand, it is usually a sign that the problem should be solved another way.

Might it be possible to write a script that finds these gaps? And speaking of added connections, might it be possible to write a program to find forest roads that connect the ends of trails, and span less than 0.5 miles in length or so, and automagically add those road segments to our network?

That will be the subject of the next blog post.

Download the source code, available here. (For illustrative purposes only, because I am sure there is a better solution.)

Related Posts:

Trails and Graph Theory: The Series

Landing page for the entire series of posts on graph theory and trails.

Trails and Graph Theory 10: Databook

Long distance hikers, especially old-school hikers before hiking apps, have been accustomed to rely upon “databooks“, a brief text summary of the important intersections, water sources, features, and distances of a long hike. Here is a portion of a page from the Appalachian Trail 1987 databook, a slim pocket-size saddle-stapled booklet.

And this is an image from a portion of the databook for the Northern New Mexico Loop (by Blisterfree, the creator of the Grand Enchantment Trail).

Not all backpackers use a databook during their hike, even pre-app, but it is still useful for planning prior to the trip. Some hikers do find an advantage in reading a text summary on trail, rather than squinting at a topo map representation on their smartphone screen. Hike your own hike; your mileage may vary; use tools and do not let them use you.

In our previous post, we were able to find a long loop route in the Gila National Forest (that does not repeat any segment), in the form of a numeric list of nodes. Our graph imported via OSMnx from OpenStreetMap includes extra information as node attributes and edge attributes, that we can extract in a Pythony way to make our own databook of sorts. Here is an example of the text information in one edge of our graph.

{'osmid': 243511168, 'ref': '165', 'name': 'Pine Flat Trail #165', 'highway': 'path', 'oneway': False, 'reversed': [False, True], 'length': 4374.8279999999995, 'from': [2509057410, 2509594627, 2509057291, 2509594252, 2509057424, 2509594385, 2509594643, 2509594263, 2509594657, 2509594529, 2509057314, 2509594279, 2509594664, 2509594284, 2509594543, 2509594673, 2509594294, 2509057334, 2509594680, 2509594556, 2509594176, 2509594306, 2509594691, 2509594435, 2509057350, 2509594568,
<truncated>

From that information we can generate a summary from our Euler graph.

  MILE   ELEVATION  TRAIL                         INTERSECTION                      
...
  105.0  6509       Hells Hole Trail #268         West Fork Trail #151          
  105.5  6529       Trail #785                    Hells Hole Trail #268         
  106.0  6483       West Fork Trail #151          Trail #785                    
  109.7  6834       Turkey Creek Trail #155       West Fork Trail #151          
  111.8  7851       Mogollon Creek Trail #153     Turkey Creek Trail #155     
...

The Python code is fairly short (excluding the altitude part, which is not included in our OSM data). The column format is enabled by the library columnar, exceedingly easy to use.

def report_path(Q, txt='',freedom_units=False,lat_long=False):
    meters_to_mile = 0.0006213712
    meters_to_feet = 3.28084
    meters_to_km = 0.001
    
    if nx.is_eulerian(Q):
        total_length = 0
        data = []
        for u,v,k in nx.eulerian_circuit(Q, keys=True):
             data_line = []
             if freedom_units:
                 data_line.append(str(round(total_length*meters_to_mile,1)))
             else:
                 data_line.append(str(round(total_length*meters_to_km,1)))
                 
             d = Q.get_edge_data(u,v)[k]
             length=0                 
             if 'length' in d:
                 length = d['length'] #meters
                 total_length += length
             node_attr = Q.nodes(data=True)[u]
             long = node_attr['x']
             lat = node_attr['y']
             if lat_long:
                 data_line.append(lat)
                 data_line_append(long)
             elevation = get_elevation(lat,long) #in meters
             elevation *= meters_to_feet #freedom units
             if freedom_units:
                 data_line.append(str(round(elevation*meters_to_feet)))
             else:
                 data_line.append(str(round(elevation)))

             trail_name=''
             if 'name' in d:
                 trail_name = d['name']
                 if not isinstance(trail_name, str):
                     trail_name = ','.join(trail_name) #in case is a list
             data_line.append(trail_name)                     
             
             crossing = set()
             crossing_txt = ''
             for u,neightbor,k,d in Q.edges(u,data=True, keys=True):
                 if 'name' in d:
                     name = d['name']
                     if not isinstance(name, str): #FIX with a foreach
                         name = ','.join(name) #in case is a list
                     if name != trail_name:
                         crossing.add(name)
             if len(crossing)&gt;0:
                 crossing_txt = ','.join(crossing)
             data_line.append(crossing_txt)
                              
             data.append(data_line)                 
        header = []
        if freedom_units:
            header.append('MILE')
        else:
            header.append('KM')
        if lat_long:
            header.append('LAT')
            header.append('LONG')
        header.append('ELEVATION')
        header.append('TRAIL')
        header.append('INTERSECTION')
        table = columnar(data, header, no_borders=True)
        print(table)     
    else:
        print('The graph is not eulerian')

One alternative, that only occurred to me after I had written this, is to output the route as CSV (comma separated values), import into a spreadsheet, and make any changes to presentation and format within the spreadsheet, and eventually export as PDF.

We still do not have water sources, towns, fence or road crossings, etc. I will go off and have a long think on how to include those features automagically.

Download source code here.

Related Posts: