Al traçar una trajectòria a traves d’un graf de trams de carrer de manera que a partir d’un punt inicial es calculi fins allà on arriba una distancia determinada seguint el graf, és molt probable que aquesta distancia no coincideixi amb un nus concret sinó que estigui en un punt intermedi d’un segment.
En aquest article es presenta una forma d’obtenir la cobertura a partir d’un punt inicial fins que s’arriba a una distancia màxima seguint l’arbre del graf de trams de carrer, la contribució principal és el càlcul del darrer segment parcial en totes les trajectòries.
Concretament, en aquest article s’explica de quina manera s’ha resolt el problema que representa obtenir els trams complerts i parcials que conformen el càlcul d’una distància recorreguda pel graf de carrers, mitjançant la funció ‘pgr_withPointsDD’ del pgrouting des de punts exteriors al graf, es a dir, que no formen part dels nodes de la xarxa.
Per poder utilitzar el codi d’aquest article s’ha de tenir instal·lat en el PostgresSQL l’extensió ‘pgrouting’ que ens permet fer anàlisi de xarxes, per mes informació visiteu la pagina web (http://docs.pgrouting.org) i (http://pgrouting.org).
ATENCIO: Els següents passos modificaran la base de dades PostgresSQL, cal anar amb molt de compte ja que podeu inutilitzar bases de dades ja existents.
- Descarregar la base de dades PostgresSQL, feu clic aquí.
- Crear base de dades nova amb nom ‘db_temp’.
- Restaurar la base de dades descarregada.
Un cop restaurada la base de dades PostgresSQL, s’ha de tenir 3 entitats:
- “orig_network” (xarxa de carrers)
- “origin” (punt d’inici)
- “orig_network_vertices_pgr” (nodes de la xarxa de carrers)
Per entendre millor el procediment s’ha plantejat un exemple concret basat en un graf de trams de carrer totalment simètric i ortogonal, però seria aplicable a qualsevol graf.
Es parteix d’una xarxa ortogonal de carrers (“orig_network”) de 2.000 x 2.000 metres, amb trams de carrers de 100m de node a node i el punt inicial (“origin”) representat en vermell, tal i com es pot veure en la figura 1.
Fig 1. Xarxa ortogonal de carrers de 100m.
Procediment pel càlcul de cobertura mitjançant pgrouting.
Per calcular una cobertura de 500m des del punt vermell seguirem les següents passes:
- Fer un duplicat de la xarxa al qual anomenarem ‘network’.
- Crear una taula anomenada ‘poi_tmp’ on tindrà totes les geometries dels punts inicials, en el nostre cas només tindrà un valor.
- Modificar la taula ‘poi_tmp’ afegint els camps (x,y,edge_id,side,fraction, newPoint);
DO $$
#variable_conflict use_variable
DECLARE dist_var INT;
BEGIN
dist_var:=500;DROP TABLE IF EXISTS "network";
CREATE TABLE IF NOT EXISTS "network" as (SELECT * FROM "orig_network");
drop table if exists poi_tmp;
CREATE TABLE if not exists poi_tmp as (SELECT ST_Centroid(tmp."geom") the_geom,tmp."id" as pid from (SELECT * FROM "origin") tmp);
ALTER TABLE poi_tmp ADD COLUMN x FLOAT;
ALTER TABLE poi_tmp ADD COLUMN y FLOAT;
ALTER TABLE poi_tmp ADD COLUMN edge_id BIGINT;
ALTER TABLE poi_tmp ADD COLUMN side CHAR;
ALTER TABLE poi_tmp ADD COLUMN fraction FLOAT;
ALTER TABLE poi_tmp ADD COLUMN newPoint geometry; |
- Actualitzar el camp ‘edge_id’ amb l’identificador del tram més proper.
- Actualitzar el camp ‘fraction’ amb la fracció corresponent del tram abans trobat.
UPDATE "poi_tmp" set "edge_id"=nearest_section."tram_id" from (SELECT distinct on(Poi.pid) Poi.pid As Punt_id,Sg.id as Tram_id, ST_Distance(Sg.the_geom,Poi.the_geom) as dist FROM "network" as Sg,"poi_tmp" AS Poi ORDER BY Poi.pid, ST_Distance(Sg.the_geom,Poi.the_geom), Sg.id) nearest_section where "poi_tmp"."pid"=nearest_section."punt_id";
UPDATE "poi_tmp" SET fraction = ST_LineLocatePoint(e.the_geom, "poi_tmp".the_geom),newPoint = ST_LineInterpolatePoint(e.the_geom, ST_LineLocatePoint(e.the_geom, "poi_tmp".the_geom)) FROM "network" AS e WHERE "poi_tmp"."edge_id" = e.id; |
- Crear la taula ‘tbl_final_nodes_tmp’ que contindrà les ‘id’, els costos agregats i el node de partida de tots els nodes pels quals passa la funció ‘pgr_withPointsDD’.
- Crear la taula ‘geo_final_nodes_tmp’ que contindrà els nodes anteriors amb la geometria corresponent.
DROP FUNCTION IF EXISTS coverage();
DROP TABLE IF EXISTS tbl_final_nodes_tmp;
CREATE TABLE IF NOT EXISTS tbl_final_nodes_tmp AS(SELECT node,agg_cost,start_vid FROM pgr_withPointsDD('SELECT id, source, target, cost, reverse_cost FROM "network" ORDER BY "network".id','SELECT pid, edge_id, fraction, side from "poi_tmp"',array(select "pid"*(-1) from "poi_tmp"),dist_var,driving_side := 'b',details := false));
DROP table if exists geo_final_nodes_tmp;
CREATE TABLE IF NOT EXISTS geo_final_nodes_tmp as (select "orig_network_vertices_pgr".*,"tbl_final_nodes_tmp"."agg_cost", "tbl_final_nodes_tmp"."start_vid", dist_var from "orig_network_vertices_pgr","tbl_final_nodes_tmp" where "orig_network_vertices_pgr"."id" ="tbl_final_nodes_tmp"."node" order by "tbl_final_nodes_tmp"."start_vid" desc,"tbl_final_nodes_tmp"."agg_cost"); |
Després d’executar el codi anterior obtenim la imatge mostrada a la figura 2. El punt blau indica la projecció del punt de partida sobre el tram més proper, que serà d’on començaran tots els camins.
Fig 2. Nodes que queden dins de la cobertura.
El que ens retorna la funció ‘pgr_withPointsDD’ son els nodes que queden dins de la cobertura, però si volem saber els trams que composen la cobertura de 500 metres, veurem que no sabem quins son els trams finals de la cobertura que quedarien incomplerts.
Per tal de resoldre aquest problema, s’ha de fer el següent:
- Crear la taula ‘final_sections_tmp’ on tindrem tots els trams que formen part de la cobertura , inclosos els trams finals dels quals només és vàlida una fracció del mateix.
DROP table IF EXISTS final_sections_tmp;
CREATE TABLE IF NOT EXISTS final_sections_tmp as (select "network"."id","network"."the_geom","geo_final_nodes_tmp"."id" as node,"geo_final_nodes_tmp"."agg_cost" as coste,(dist_var "geo_final_nodes_tmp"."agg_cost") as falta,"geo_final_nodes_tmp"."start_vid" as id_punt, (select case when (dist_var-"geo_final_nodes_tmp"."agg_cost") / ST_Length("network"."the_geom")<=1 then (dist_var-"geo_final_nodes_tmp"."agg_cost") / ST_Length("network"."the_geom") when (dist_var- "geo_final_nodes_tmp"."agg_cost") / ST_Length("network"."the_geom")>1 then (1) end) as fraccio from "network","geo_final_nodes_tmp" where ST_DWithin("geo_final_nodes_tmp"."the_geom", "network"."the_geom",1)=TRUE); |
En la figura 3 es mostra el resultat d’executar el codi anterior.
Fig 3. Trams que formen part de la cobertura de 500 metres.
Per tal d’obtenir les fraccions del trams finals que completen la cobertura de 500 metres, s’ha creat una funció en ‘PL/pgSQL’, la qual realitza les modificacions necessàries.
- Crear la funció ‘Coverage’ que realitzarà el procés.
- Crear la taula ‘fraction_sections_raw’ que contindrà el trams (complerts i fraccions) que composen la cobertura de 500 metres.
CREATE OR REPLACE FUNCTION coverage() RETURNS SETOF final_sections_tmp AS
$BODY$
DECLARE
r final_sections_tmp%rowtype;
m final_sections_tmp%rowtype;
BEGIN
DROP TABLE IF EXISTS fraction_sections_raw;
CREATE TABLE fraction_sections_raw (the_geom geometry, punt_id bigint,id_tram bigint,fraccio FLOAT,node bigint,fraccio_inicial FLOAT,cost_invers FLOAT,cost_directe FLOAT,target bigint,radi_inic FLOAT);
FOR r IN SELECT "final_sections_tmp".* FROM "final_sections_tmp" WHERE "final_sections_tmp"."id" not in (select "edge_id" from "poi_tmp")
LOOP
insert into fraction_sections_raw VALUES(ST_Line_Substring((r."the_geom"),case when (select ST_Line_Locate_Point((r."the_geom"),(select "geo_final_nodes_tmp"."the_geom" from "geo_final_nodes_tmp" where "geo_final_nodes_tmp"."id"=r."node" and "geo_final_nodes_tmp"."start_vid"=r."id_punt")))<0.001 then 0 else 1-r."fraccio"
END,
case when (select ST_Line_Locate_Point((r."the_geom"),(select "geo_final_nodes_tmp"."the_geom" from "geo_final_nodes_tmp" where "geo_final_nodes_tmp"."id"=r."node" and "geo_final_nodes_tmp"."start_vid"=r."id_punt")))<0.001 then r."fraccio" else 1
END),r."id_punt"*(-1),r."id",0,r."node",0,0,0,0);
RETURN NEXT r;
END LOOP;
FOR m IN SELECT "final_sections_tmp".* FROM "final_sections_tmp" WHERE "final_sections_tmp"."id" in (select "edge_id" from "poi_tmp")
LOOP
insert into fraction_sections_raw VALUES(m."the_geom",m."id_punt"*(-1),m."id",0,m."node",0,0,0);
RETURN NEXT m;
END LOOP;
RETURN;
END
$BODY$
LANGUAGE 'plpgsql' ;
PERFORM "the_geom" FROM coverage(); |
Tot seguit es modifica la taula ‘fraction_sections_raw’ per tal d’obtenir els trams finals que corresponen a la cobertura de 500 metres.
- Actualitzar ‘fraccio_inicial’ amb el valor que li pertoca.
- Actualizar ‘’cost_directe’, ‘cost_invers’ amb els valors corresponents.
- Actualitzar ‘fraccio’ amb el valors calculats dels trams finals.
- Modificar la geometria segons el camp ‘fraccio’ calculat.
update "fraction_sections_raw" set "fraccio_inicial"="poi_tmp"."fraction" from "poi_tmp" where "id_tram"="edge_id";
update "fraction_sections_raw" set "cost_directe"="network"."cost", "target"="network"."target", "cost_invers"="network"."reverse_cost" from "network" where "id_tram"="id";
UPDATE "fraction_sections_raw" SET "fraccio"=((case when ("fraction_sections_raw"."fraccio_inicial" * ST_Length("fraction_sections_raw"."the_geom"))>dist_var then (dist_var/ST_Length("fraction_sections_raw"."the_geom")) else "fraction_sections_raw"."fraccio_inicial" end)+(case when ((1-"fraction_sections_raw"."fraccio_inicial") * ST_Length("fraction_sections_raw"."the_geom"))>dist_var then (dist_var/ST_Length("fraction_sections_raw"."the_geom")) else (1-"fraction_sections_raw"."fraccio_inicial") end));
update "fraction_sections_raw" set "the_geom"=final."the_geom"from(select distinct(ST_Line_Substring((m."the_geom"),(case when (select ST_Line_Locate_Point((m."the_geom"), (select "the_geom" from "geo_final_nodes_tmp" where "geo_final_nodes_tmp"."id"=m."node" and "geo_final_nodes_tmp"."start_vid" = m."punt_id"*-1)))<0.01 then 0 else 1-m."fraccio" END),(case when (select ST_Line_Locate_Point((m."the_geom") , (select "the_geom" from "geo_final_nodes_tmp" where "geo_final_nodes_tmp"."id"=m."node" and "geo_final_nodes_tmp"."start_vid" = m."punt_id"*-1)))<0.01 then m."fraccio" else 1 END))) the_geom,m."id_tram"from "fraction_sections_raw" m where m."id_tram" in (select "edge_id" from "poi_tmp")) final where final."id_tram" ="fraction_sections_raw"."id_tram";
insert into "fraction_sections_raw" (select SX."the_geom",PI."pid" as punt_id,SX."id"as id_tram,999 as fraccio,SX."source" as node,PI."fraction" as fraccio_inicial , SX."cost",SX."reverse_cost" from "network" SX inner join (Select "edge_id","pid","fraction" from "poi_tmp") PI on SX."id"=PI."edge_id");
UPDATE "fraction_sections_raw" set "the_geom"=final."the_geom" from (select ST_Line_Substring((SXI."the_geom"),(case when (FT."fraccio_inicial" - (dist_var/ST_Length(SXI."the_geom")))>0 then (FT."fraccio_inicial"-(dist_var/ST_Length(SXI."the_geom"))) else 0 end),(case when (FT."fraccio_inicial" + (dist_var/ST_Length(SXI."the_geom")))<1 then (FT."fraccio_inicial"+(dist_var/ST_Length(SXI."the_geom"))) else 1 end)) as the_geom, FT."punt_id", FT."id_tram", FT."fraccio" from "fraction_sections_raw"FT inner join (select SX."the_geom" as the_geom, SX."id" as tram_xarxa from "orig_network" SX, "poi_tmp" PI where SX."id"=PI."edge_id") SXI on FT."id_tram"=SXI.tram_xarxa where FT."fraccio"=999) final where "fraction_sections_raw"."punt_id"=final."punt_id" and "fraction_sections_raw"."fraccio"=999; |
- Crear la taula ‘fractions_sections_tmp’ que contindrà els trams sense duplicats de la cobertura de 500m.
- Crear la taula ‘driving_distance’ que contindrà la unió de tots els trams.
DROP TABLE IF EXISTS fraction_sections_tmp;
CREATE TABLE fraction_sections_tmp AS (select distinct(the_geom),punt_id,radi_inic from fraction_sections_raw);
DROP TABLE IF EXISTS driving_distance;
CREATE TABLE IF NOT EXISTS driving_distance AS (Select ST_Union(TOT.the_geom) the_geom, TOT."punt_id" from (select the_geom,punt_id,radi_inic from fraction_sections_tmp) TOT group by TOT."punt_id");
END $$; |
Un cop executat tot el codi el resultat és el que es mostra en la figura 4.
Fig 4. Cobertura final de 500 metres.
En aquesta figura es poden veure els segments finals dels trajectes a partir del punt projectat amb la seva autèntica dimensió, implicant per tant només un fragment del tram sencer correspondent, es pot comparar amb la figura 3 per veure la diferència de resultats.
Conclusions
Com a conclusió finals del article, podem dir que s’ha creat una funció per obtenir la cobertura a partir d’un punt amb distància o funció de cost constant a partir d’un punt d’inici projectat sobre del graf, donant un resultat exacte.
Millorant, així, el resultat que es podria obtenir només amb les funcions pròpies de la llibreria del ‘pgrouting’.
La funció dissenyada es aplicable a qualsevol funció de rutatge en la que es parteixi d’uns punts inicials, que formin, o no, part del graf, i obtenir la cobertura desitjada.
El codi sencer es pot descarregar del següent enllaç.
