Python – How to Convert Shapefile into Graph for Dijkstra’s Algorithm

dijkstragraphpythonshapefile

I am working in a simulation using NetLogo, the idea is to make the agents interact inside the roads of a simplified city, in order to do that I need to use some algorithms to find the fastest path between two points in the city (for example Dijkstra or A*), and to do that I need to work with a graph of the city (with nodes and vertices) but I don't know how to convert the shapefile into a graph.

My idea was to get all the intersections between two streets and start building the graph connecting those nodes. I have tried networkx for python but I had some problems installing the library gdal. Currently I'm using the s2g library (here) to do it, but for example for a shp file like this

example shapefile

I always get this:

INFO:root:Validating pair-wise line connections of raw shapefiles (total 3 lines)
100% (3 of 3) |##########################################################################################################| Elapsed Time: 0:00:00 Time: 0:00:00
Major components statistics:
    Total components: 0
    Component size: max 0, median 0.0, min 0, average 0.0
    Top comp. sizes: 0

I suppose that is an empty graph so I cannot work with it. Perhaps there is a better approach to my problem or an easy tool to convert shp to graphs. I mac user so I have no access to spatialite.

Best Answer

The sg module uses Fiona to read the shapefiles (see shapegraph.py ) and if you can use the module, therefore Fiona is installed.

If you cannot use nx_shp.py (because of osgeo) and you have problems with sg, you can use Fiona and Networkx to create a Networkx Graph. (GSE: How to calculate edge length in Networkx for example).

from shapely.geometry import shape
import fiona
geoms =[shape(feature['geometry']) for feature in fiona.open("stac_graphe.shp")]
import itertools
# create a Graph
import networkx as nx
G = nx.Graph()
for line in geoms:
   for seg_start, seg_end in itertools.izip(list(line.coords),list(line.coords)[1:]):
    G.add_edge(seg_start, seg_end)

Result

enter image description here

You can also create a Planar Graph

from shapely.ops import unary_union
res = unary_union(geoms)
G = nx.Graph()
for line in res:
    for seg_start, seg_end in itertools.izip(list(line.coords),list(line.coords)[1:]):
        G.add_edge(seg_start, seg_end)

enter image description here

New

There is a problem if some geometries are Multigeometries
For example:

 geoms =[shape(feature['geometry']) for feature in fiona.open("multiline.shp")]
 for line in geoms:
     print line
 MULTILINESTRING ((3 4, 10 50, 20 25), (-5 -8, -10 -8, -15 -4))

With geoms[0], you decompose the Multigeometry 

 for line in geoms[0]:
     print line

 LINESTRING (3 4, 10 50, 20 25)
 LINESTRING (-5 -8, -10 -8, -15 -4)

Therefore the script becomes

  if line.geom_type== "MultiLineString":
        ....
   else:
Related Question