0

I'm trying to build a distance matrix from all municipalities to all airports in the Canary Islands using pgrouting. I've had problems with disconnected islands and so on using OpenStreetMap data so I've moved towards using official data from the Spanish geographic information centre. This data only contains interurban roads but since some streets are classified as such, I have enough nodes (in theory) to get relevant results for all municipalities. The problem is, on the one hand I'm losing a lot of nodes when performing ST_Drivingdistance even with unrealistically high distance limits and when trying to set the nearest node to the centroid of each municipality, 20% of the municipalities lose the node (I'm not too worried about this because I can always use the median distance from all nodes in each municipality or something of the sort, but it's still annoying).

I'm running the following code.

To create the network

-- Create new index col to avoid out of range error
ALTER TABLE tramos_canarias ADD COLUMN id serial;
select * from tramos_canarias limit 1;

-- force 2d
UPDATE tramos_canarias SET the_geom = ST_Force2D(the_geom)

-- Create network

SELECT pgr_nodeNetwork('tramos_canarias', 0.001, id:= 'id', the_geom:= 'the_geom', table_ending:='noded');
select * from tramos_canarias_noded
limit 5;

SELECT  pgr_createTopology('tramos_canarias_noded', 0.001,
        the_geom:='the_geom', id:='id', source:='source', target:='target');

-- read

select * from tramos_canarias_noded_vertices_pgr
limit 50;

-- analyse

select pgr_analyzegraph('tramos_canarias_noded', 0.000001, the_geom:='the_geom', id:='id', source:='source', target:='target'); 

To determine the closest node to each airport

The same procedure is used with the municipal centroids using SELECT ST_Centroid(the_geom) from municipalities.

-- Create airport table

create table aeropuertos_canarias as (SELECT DISTINCT nombre, cod_iata, wkb_geometry  from rt_aerodromo_p WHERE cod_iata in ('LPA', 'VDE', 'FUE', 'GMZ', 'SPC', 'ACE', 'TFN', 'TFS'));
ALTER TABLE aeropuertos_canarias
   ADD COLUMN nearest_node integer;

select * from aeropuertos_canarias
-- define CRS

SELECT ST_SRID(wkb_geometry) from aeropuertos_canarias

SELECT DISTINCT ST_SRID(the_geom) from ways_vertices_pgr

SELECT UpdateGeometrySRID('aeropuertos_canarias', 'wkb_geometry', 4326)

select ST_Force2D(wkb_geometry) from aeropuertos_canarias

-- Create table with nearest nodes
CREATE TABLE temp AS
   SELECT a.nombre, b.id, min(a.dist)
   FROM
     (SELECT distinct airport.nombre, 
             min(ST_Distance(airport.wkb_geometry::geography, node.the_geom::geography)) AS dist
      FROM aeropuertos_canarias as airport, tramos_canarias_noded_vertices_pgr as node
      GROUP BY airport.nombre) AS a,
     (SELECT airport.nombre, node.id, 
             ST_Distance(airport.wkb_geometry::geography, node.the_geom::geography) AS dist
      FROM aeropuertos_canarias as airport, tramos_canarias_noded_vertices_pgr as node) AS b
   WHERE a.dist = b. dist
         AND a.nombre = b.nombre
   GROUP BY a.nombre, b.id;
select * from temp; 

-- Join nearest nodes
UPDATE aeropuertos_canarias
SET nearest_node = 
   (SELECT id 
    FROM temp
    WHERE temp.nombre = aeropuertos_canarias.nombre);
    
select * from aeropuertos_canarias

This brings the expected results

To create the driving distance table

-- Get all distances
CREATE TABLE matriz_aeropuertos_canarias as (with matriz as (SELECT * FROM pgr_drivingDistance('SELECT id, 
                         source, target, ST_Length(ST_Transform(the_geom,25828)) 
                         AS cost, - ST_Length(ST_Transform(the_geom,25828)) 
                         AS reverse_cost FROM tramos_canarias_noded',
       array(select distinct nearest_node from aeropuertos_canarias), 1000000000, false)
       JOIN aeropuertos_canarias on aeropuertos_canarias.nearest_node = pgr_drivingDistance.from_v)
       SELECT node, nombre, agg_cost from matriz);

select * from matriz_aeropuertos_canarias limit 10;

-- Join geometry of each node

ALTER TABLE matriz_aeropuertos_canarias ADD COLUMN geom geometry;

select * from tramos_canarias_noded_vertices_pgr limit 10;

UPDATE matriz_aeropuertos_canarias SET geom = tramos_canarias_noded_vertices_pgr.the_geom from tramos_canarias_noded_vertices_pgr
                                               WHERE tramos_canarias_noded_vertices_pgr.id = matriz_aeropuertos_canarias.node;


-- Check all airports are taken into account - works

SELECT count(distinct nombre) from matriz_aeropuertos_canarias

-- Check all nodes have geometry - works

SELECT sum(case when matriz_aeropuertos_canarias.geom is NULL then 1 else 0 end) from matriz_aeropuertos_canarias

-- Check all nodes are represented - does not work

select * from matriz_aeropuertos_canarias;

The result is that many nodes are missing when they, in theory, should not be.

Right side, all nodes in the distance matrix. Left side, all nodes in the network

Results of two queries: distance matrix and network

Your Answer

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

Browse other questions tagged or ask your own question.