37

Use PostGIS topologies to clean-up road networks

 3 years ago
source link: https://blog.mathieu-leplatre.info/use-postgis-topologies-to-clean-up-road-networks.html
Go to the source link to view the article. You can view the picture content, updated content and better typesetting reading experience. If the link is broken, please click the button below to view the snapshot at that time.
Use PostGIS topologies to clean-up road networks

Use PostGIS topologies to clean-up road networks

Wed 03 July 2013

This article gives a few basics to get started with using the PostGIS topology extension.

We will take avandtage of topologies to clean-up a real topological road network, coming from Toulouse OpenData files.

toulouse-voirie.png

Topological

A topology is a general concept, where objects are defined by their relationships instead of their geometries. Instead of lines, we manipulate edges, vertices and faces : might remind you the core concepts of graph theory.

A topological road network is supposed to have their lines (edges) connected at single points (nodes).

In this example dataset, JOSM validator detects not less than 1643 errors :) Broken connections, crossing lines ...

toulouse-voirie-error.pngtoulouse-voirie-error2.pngtoulouse-voirie-error3.png

Let's clean this up !

Installation

On Ubuntu 12.04, you just have to install PostGIS :

sudo apt-add-repository -y ppa:ubuntugis/ppa
sudo apt-get update
sudo apt-get install -y postgresql postgis

The topology extension is installed by default. Just activate it in your database:

CREATE DATABASE "roadsdb";
CREATE EXTENSION postgis;
CREATE EXTENSION postgis_topology;
SET search_path = topology,public;

Data Import

Load your shapefile (using command-line) like usual :

schema="public."
db="roadsdb"
user="postgres"
password="postgres"
host="localhost"
ogr2ogr -f "PostgreSQL" PG:"host=${host} user=${user} dbname=${db} password=${password}" -a_srs "EPSG:2154" -nln ${schema}roads -nlt MULTILINESTRING ROAD_SHAPEFILE.SHP

Create and associate the PostGIS topology:

SELECT topology.CreateTopology('roads_topo', 2154);
SELECT topology.AddTopoGeometryColumn('roads_topo', 'public', 'roads', 'topo_geom', 'LINESTRING');

Convert linestrings to vertices and edges within the topology :

-- Layer 1, with 1.0 meter tolerance
UPDATE roads SET topo_geom = topology.toTopoGeom(wkb_geometry, 'roads_topo', 1, 1.0);

From now on, we have a topology, whose imperfections were corrected. It smoothly merged all dirty junctions, whose defects were at most 1.0 meter wide.

toulouse-voirie-clean.png

You may enconter insertion problems : the tool fails [1] and aborts the whole transaction. Use this snippet to skip errors and go on with the next records:

DO $$DECLARE r record;
BEGIN
  FOR r IN SELECT * FROM roads LOOP
    BEGIN
      UPDATE roads SET topo_geom = topology.toTopoGeom(wkb_geometry, 'roads_topo', 1, 1.0)
      WHERE ogc_fid = r.ogc_fid;
    EXCEPTION
      WHEN OTHERS THEN
        RAISE WARNING 'Loading of record % failed: %', r.ogc_fid, SQLERRM;
    END;
  END LOOP;
END$$;

This is rather frustrating to face topological errors at insertion ! You can try with a lower tolerance, or check that your records have at least valid geometries. Any clarification or help on this would be welcome :)

Visualize and export

In order to visualize your topology vertices in QGIS, browse your database tables, and add the following layers: roads_topo.edge_data and roads_topo.node.

toulouse-voirie-topology.png

You can also export the resulting geometries into a new table :

CREATE TABLE roads_clean AS (
    SELECT ogc_fid, topo_geom::geometry
    FROM roads
);

Or obtain your lovable Shapefile in return :

ogr2ogr -f "ESRI Shapefile" ROAD_CLEAN.SHP PG:"host=${host} user=${user} dbname=${db} password=${password}" -sql "SELECT topo_geom::geometry FROM roads"

If, like Amit you want to split the lines at intersections and assign original attributes, just join roads_topo.edge_data and on the roads table :

SELECT r.lib_off, r.ogc_fid, e.geom
FROM roads_topo.edge e,
     roads_topo.relation rel,
     roads r
WHERE e.edge_id = rel.element_id
  AND rel.topogeo_id = (r.topo_geom).id

Going further...

We could collapse crossing lines and disconnected junctions into a nice and clean network.

Yes ahem, we weren't able to repair every topological error of this dataset using this automatic method. Some inconsistencies, like the following one, are like 6 meters wide ! They are, by the way, perfectly described in OpenStreetMap :

toulouse-voirie-error4.png

We could also play with simplifications using Sandro Santilli's SimplifyEdgeGeom [2] function, it will collapse edges with a higher tolerance ...

SELECT SimplifyEdgeGeom('roads_topo', edge_id, 1.0) FROM roads_topo.edge;

Don't hesitate to share your thoughts and feedback. Concrete use cases and examples are rare about this! And as usual, drop a comment if anything is wrong or not clear :)

[1]SQL/MM Spatial exception, geometry intersects edge, side location conflict, ...[2]Just execute the function SQL code. It's just an elegant wrapper around ST_ChangeEdgeGeom and ST_Simplify.

#postgis, #sql, #topology, #opendata, #toulouse - Posted in the Dev category


© Copyright 2020 by Mathieu Leplatre. mnmlist Theme

Content licensed under the Creative Commons attribution-noncommercial-sharealike License.


About Joyk


Aggregate valuable and interesting links.
Joyk means Joy of geeK