Another rainy season Saturday got me time for some GIS exercises. Inspired by Carson Farmer's
excellent introduction to
pgRouting I decided to set up a routing enabled PostGIS database with OpenStreetMap data and compare it with my
GRASS GIS routing setup.
Contrary to Carson Farmer I didn't compile pgRouting from source but installed it from the
ubuntugis-unstable repository
sudo apt-get install postgres-8.4-pgrouting postgres-8.4-pgrouting-dd postgres-8.4-pgrouting-tsp
and followed the instruction from the mentioned blog.
The osm2pgrouting script creates all necessary tables including table "public.ways", which is probably the most important one. The table schema without indexes and constraints looks like:
Table "public.ways"
Column | Type | Modifiers
--------------+------------------+-----------
gid | integer | not null
class_id | integer | not null
length | double precision |
name | character(200) |
x1 | double precision |
y1 | double precision |
x2 | double precision |
y2 | double precision |
reverse_cost | double precision |
rule | text |
to_cost | double precision |
osm_id | integer |
the_geom | geometry |
source | integer |
target | integer |
The relevant columns for the routing are "length","reverse_cost" and "to_cost". The latter two are the forward and backward costs to use a way and are used in directed shortest path calculations. It's necessary to use an algorithm that can handle direction when one-ways have to be taken in account.
Since OpenStreetMap data are in geographic coordinates the length is in degree, but I prefer to have the length of a way in meter. Thankfully PostGIS provides
ST_Distance_Spheroid to calculate lengths on the (earth) spheorid:
UPDATE ways
SET length = ST_Length_Spheroid(the_geom,'SPHEROID["WGS 84",6378137,298.257223563]')::DOUBLE PRECISION;
Next step is to filter the one-ways paths and set the "reverse_cost" to -1 to force pgRouting to choose these ways only in forward direction. Unfortunately only the name and the osm_id are imported from the OpenStreetMap raw data to the database. Thus I needed a list of osm_ids that represents one-ways. Since I also used a PostgreSQL database to store the GRASS attribute tables while working with OpenStreetMap data, I could easily get a comma-separated list of osm_ids with the following query:
psql -A -t -R ',' -d grass -U grass -c "SELECT osm_id FROM osm_net WHERE oneway = 'yes';" > oneways.sql
Then I updated the "reverse_cost" column:
UPDATE ways SET reverse_cost = -1
WHERE osm_id IN(55909972,55909973,28977387, ... ,102306365,27158132)
Finally I could start with undirected and directed queries and compare the results.
-- Undirected shortest path query with astar
SELECT * FROM astar_sp( 'ways',3365,3412)
-- Directed shortest path query with astar
SELECT * FROM astar_sp_directed( 'ways',3365,3412,true,true)
Verify both calculated shortest paths from node 3365 to node 3412 in QGIS:
The red path is the undirected query and ignores apparently one-ways, contrary to the directed green path that follows one-way directions.
Until now I was considering only length and direction of a way. To get a more life-like model I wanted to take the road classes into account, too. I updated table "public.classes" and estimated for each road class an average speed in m/s in column "cost":
psql -A -d osm_routing -U openstreet -c "SELECT id,type_id,TRIM(trailing FROM name) AS name,cost FROM classes LIMIT 8;"
id|type_id|name|cost
101|1|motorway|33.3333
102|1|motorway_link|27.7778
103|1|motorway_junction|27.7778
104|1|trunk|27.7778
105|1|trunk_link|27.7778
106|1|primary|13.8889
107|1|primary_link|13.8889
108|1|secondary|11.1111
(8 rows)
With these assumed costs for each road class I updated the costs in table "public.ways" and divided the length by the average speed and got the time in seconds that is necessary to pass a certain way.
-- Update to_cost column
UPDATE ways
SET to_cost=(length/(SELECT "cost" FROM classes WHERE id=class_id))::DOUBLE PRECISION
-- Update reverse_cost columns
UPDATE ways
SET reverse_cost=(length/(SELECT "cost" FROM classes WHERE id=class_id))::DOUBLE PRECISION
WHERE reverse_cost <> -1;
Now everything is set up and it needs nothing but another rainy weekend to compare pgRouting results with GRASS GIS results.
Last but not least I'd like to share yet another way how to display results of PostGIS queries without using an intermediate (Shape-)file, namely using the very versatile
Virtual Format from the OGR library. Create a text file with the following lines and open it as a vector layer in QGIS or any other OGR backed application:
<OGRVRTDataSource>
<OGRVRTLayer name="shortest_path">
<SrcDataSource>
PG:dbname='osm_routing' user='openstreet'
</SrcDataSource>
<SrcSQL>
SELECT * FROM astar_sp_directed( 'ways',3365,3414,true,true)
</SrcSQL>
</OGRVRTLayer>
</OGRVRTDataSource>
