2

I am analyzing Germany's high voltage electricity grid. I extracted the data from OSM, imported it in QGIS and simplified its structure by hand, specifically creating single nodes where power lines are connected in contrast to having a lot of high voltage substations with complicated geometries. Ultimately, I want to export that network and use it in GAMS to run an optimization model on it. But long story short. I have two .shp files, standard WGS 84 CRS.

  1. a file containing the network as lines (about 250)
  2. a file containing all the nodes (about 350) in the network (where the power lines intersect)

Now I want to extract a shortest distance matrix containing the distances between all nodes in file 2, having used file 1 as my "road" network.

I have read through a couple of answers suggesting python, GRASS or PostGIS with pgrouting but none of them worked for me. I'm new to Python, GRASS and PostGIS, which sure does not really help.

(almost forgot. Mac OS 10.7.5, GRASS 7.0, QGIS 2.8.1, PostgreSQL 9.0.13, PostGIS 2.1.7, pgrouting 2.0.0)

  • Could you elaborate on why the answers you have found do not work for you? What is the difficulty? – whuber May 14 '15 at 15:55
  • In GRASS, I can prepare the network with v.net, but neither 'nodes' nor 'connect' seem to be able to implement the network correctly. Selecting two points in the Network Analysis Tool after preparing the network will sometimes give a reasonable answer, sometimes gives wrong answers (in terms of shortest path, takes lots of detours), most of the time says "destination(s) unreachable". In PostGIS/PostgreSQL/pgrouting I couldn't really get anything to work, simple single queries like this: anitagraser.com/2011/02/07/a-beginners-guide-to-pgrouting give a blank answer. – mpetersen May 14 '15 at 17:46
  • I feel more comfortable with GRASS overall, so if it could be solved there somehow, I'd be more than happy! I uploaded my files here: dropbox.com/s/xlda4pid7ae7aal/… I guess the node file might even be redundant, since GRASS evaluates these nodes by itself. – mpetersen May 14 '15 at 17:56
2

Ok, update after another couple of hours of trying. I think I solved the problem, more or less by accident.

Since I digitized/simplified the grid structure by hand, one of the problems supposedly was, that a couple of lines would have intersections between their beginning and end points. It seems like a good idea to let QGIS handle this and use 'Split lines with lines' from the Processing toolbox.

Having done this, I imported the grid and the nodes as shp files into GRASS. The nodes are quite easy to obtain via QGIS's vector analysis tool 'Line Intersections'. Having both layers imported in GRASS, I used v.net with operation:connect and option -c to attach the node layer to the line layer (nodes as layer 2). Alternatively, operation:nodes does almost the same, because the nodes computed are the same as those I imported in my shp file from QGIS. Only difference: the nodes get different IDs, so you can't really preserve a personal ID-scheme to identify certain points in QGIS.
Nevertheless, as the network is now "prepared", the v.net.allpairs function can be used. All my lines are in layer 1, all the nodes are in layer 2. Selecting these, activating geodesic calculation (my shp files were standard WGS 84) and waiting for a few minutes did the job. Having not set any costs, the 'from_cat' column and the 'to_cat' column in the attributes table of the output vector layer are now my old QGIS IDs, and the cost column is the distance between the nodes in km.

I hope this description can be of help for others like me, not too familiar with SQL/Python etc. and looking for a quite clicky solution. I needed the shortest distance table for a GAMS model and thought this would be a piece of cake.

2

NetworkX (http://cheeseshop.python.org/pypi/networkx/) is a Python package with many functions for graph and network analysis. With this package installed you can solve the problem of generating a shortest distance matrix using the Python console in QGIS.

All you need is a network layer (the edges) with a valid topology. The edges need at least 3 attributes: the IDs of start node and end node of each edge, and the length of each edge. In the following code I expect these as the fields with index 1 to 3. You may create this data schema using an appropriate SQL query or a spatial join, followed by (field) calculation of the line length.

Make the network layer the active layer and then in the Python console:

import networkx as nx
aLayer = qgis.utils.iface.activeLayer()
DM = nx.Graph()
DM.add_weighted_edges_from([(feat[1], feat[2], feat[3]) for feat in aLayer.getFeatures()])
print nx.all_pairs_dijkstra_path_length(DM)

The print command in the last line outputs the result as a dictionary of dictionaries. Instead of printing the output to the console this structure can be used for further processing or written to file.

This way many of the algorithms in NetworkX can be used for geoprocessing in QGIS.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.