Skip to content

Instantly share code, notes, and snippets.

@rajanski
Created November 7, 2015 18:32
Show Gist options
  • Save rajanski/2e7d3df6654ce6d96afb to your computer and use it in GitHub Desktop.
Save rajanski/2e7d3df6654ce6d96afb to your computer and use it in GitHub Desktop.
OSM to directed igraph graph with node coordinates
"""
Read graphs in Open Street Maps osm format
Based on osm.py from brianw's osmgeocode
http://github.com/brianw/osmgeocode, which is based on osm.py from
comes from Graphserver:
http://github.com/bmander/graphserver/tree/master and is copyright (c)
2007, Brandon Martin-Anderson under the BSD License
"""
import xml.sax
import copy
import igraph as ig
from decimal import *
#highway_cat = 'motorway|trunk|primary|secondary|tertiary|road|residential|service|motorway_link|trunk_link|primary_link|secondary_link|teriary_link'
def download_osm(left,bottom,right,top,highway_cat):
"""
Downloads OSM street (only highway-tagged) Data using a BBOX,
plus a specification of highway tag values to use
Parameters
----------
left,bottom,right,top : BBOX of left,bottom,right,top coordinates in WGS84
highway_cat : a list of highway tag strings
Returns
----------
stream object with osm xml data
"""
#Return a filehandle to the downloaded data."""
from urllib import urlopen
#fp = urlopen( "http://api.openstreetmap.org/api/0.6/map?bbox=%f,%f,%f,%f"%(left,bottom,right,top) )
#fp = urlopen( "http://www.overpass-api.de/api/xapi?way[highway=*][bbox=%f,%f,%f,%f]"%(left,bottom,right,top) )
hw_cat_string = ''
for i in highway_cat:
hw_cat_string += i+'|'
hw_cat_string = hw_cat_string[:-1]
print "trying to download osm data from "+str(left),str(bottom),str(right),str(top)+" with highways of categories "+hw_cat_string
try:
print "downloading osm data from "+str(left),str(bottom),str(right),str(top)+" with highways of categories"+hw_cat_string
fp = urlopen( "http://www.overpass-api.de/api/xapi?way[highway=%s][bbox=%f,%f,%f,%f]"%(hw_cat_string,left,bottom,right,top) )
#slooww only ways,and in ways only "highways" (i.e. roads)
#fp = urlopen( "http://open.mapquestapi.com/xapi/api/0.6/way[highway=*][bbox=%f,%f,%f,%f]"%(left,bottom,right,top) )
return fp
except Exception,ex:
print "osm data download unsuccessful: ", ex
def read_osm(filename_or_stream, highway_cat = None,only_roads = True):
"""Read graph in OSM format from file specified by name or by stream object.
Filter by highway tags
Parameters
----------
filename_or_stream : filename or stream object
highway_cat : list of highway tag strings
only_roads : restriction to roads highways only
Returns
-------
G : Graph
Examples
--------
>>> G=nx.read_osm(nx.download_osm(-122.33,47.60,-122.31,47.61))
>>> plot([G.node[n]['data'].lat for n in G], [G.node[n]['data'].lon for n in G], ',')
"""
osm = OSM(filename_or_stream)
G = ig.Graph(directed = True)
getcontext().prec = 12
for w in osm.ways.itervalues():
if only_roads and 'highway' not in w.tags:
continue
#filter by highway_cat list
if w.tags['highway'] in highway_cat or highway_cat == None:
#Igraph_version
nodeslist = []
for n in w.nds:
node = osm.nodes[n]
try:
print G.vs.find(name = node.id)
except:
G.add_vertex(name = node.id, lat = Decimal(node.lat), lon = Decimal(node.lon))
nodeslist.append(node.id)
z1 = nodeslist[:-1]
z2 = nodeslist[1:]
edges = zip(z1, z2)
#print edges
for e in edges:
#print w.tags['highway']
G.add_edge(e[0],e[1], highway = w.tags['highway'] )
if ('oneway' not in w.tags) or (w.tags['oneway'] != 'yes'):
edges = zip(reversed(z2),reversed(z1))
for e in edges:
G.add_edge(e[0],e[1], highway = w.tags['highway'] )
print G.is_directed()
return G
class Node:
def __init__(self, id, lon, lat):
self.id = id
self.lon = lon
self.lat = lat
self.tags = {}
class Way:
def __init__(self, id, osm):
self.osm = osm
self.id = id
self.nds = []
self.tags = {}
def split(self, dividers):
# slice the node-array using this nifty recursive function
def slice_array(ar, dividers):
for i in range(1,len(ar)-1):
if dividers[ar[i]]>1:
#print "slice at %s"%ar[i]
left = ar[:i+1]
right = ar[i:]
rightsliced = slice_array(right, dividers)
return [left]+rightsliced
return [ar]
slices = slice_array(self.nds, dividers)
# create a way object for each node-array slice
ret = []
i=0
for slice in slices:
littleway = copy.copy( self )
littleway.id += "-%d"%i
littleway.nds = slice
ret.append( littleway )
i += 1
return ret
class OSM:
def __init__(self, filename_or_stream):
""" File can be either a filename or stream/file object."""
nodes = {}
ways = {}
superself = self
class OSMHandler(xml.sax.ContentHandler):
@classmethod
def setDocumentLocator(self,loc):
pass
@classmethod
def startDocument(self):
pass
@classmethod
def endDocument(self):
pass
@classmethod
def startElement(self, name, attrs):
if name=='node':
self.currElem = Node(attrs['id'], (attrs['lon']), (attrs['lat']))
elif name=='way':
self.currElem = Way(attrs['id'], superself)
elif name=='tag':
self.currElem.tags[attrs['k']] = attrs['v']
elif name=='nd':
self.currElem.nds.append( attrs['ref'] )
@classmethod
def endElement(self,name):
if name=='node':
nodes[self.currElem.id] = self.currElem
elif name=='way':
ways[self.currElem.id] = self.currElem
@classmethod
def characters(self, chars):
pass
xml.sax.parse(filename_or_stream, OSMHandler)
self.nodes = nodes
self.ways = ways
#"""
#count times each node is used
node_histogram = dict.fromkeys( self.nodes.keys(), 0 )
for way in self.ways.values():
if len(way.nds) < 2: #if a way has only one node, delete it out of the osm collection
del self.ways[way.id]
else:
for node in way.nds:
node_histogram[node] += 1
#use that histogram to split all ways, replacing the member set of ways
new_ways = {}
for id, way in self.ways.iteritems():
split_ways = way.split(node_histogram)
for split_way in split_ways:
new_ways[split_way.id] = split_way
self.ways = new_ways
#"""
@rajanski
Copy link
Author

rajanski commented Nov 7, 2015

Admittedly, the igraph version is not so neat as the networkx version. This is due to that there is no add_path function in igraph, and wiritng one (which I looked into) has to happen in the C part of igraph, which means I can#t do it right now.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment