Postgis & Pgrouting

Mise en place d’un système de routage et de calcul d’itinéraires intégrable à PostgreSQL/PostGIS : « Pgrouting »

isochrones multi modales Nantes Diouck

        I.               Introduction

Il existe de nombreuses recherches sur l’utilisation des algorithmes de pgRouting pour l’analyse du plus court chemin et le calcul des isochrones. PgRouting est une extension de PostgreSQL et PostGIS. Il fournit divers algorithmes comme All Pairs Short Path, Bi-directional Dijkstra, and A*, Shortest Path y compris Turn Restriction, Driving Distance etc… et d’autres fonctions d’analyses de graphes.               Les fonctions pgRouting sont basées sur les coûts, qui sont dynamiques par opposition coûts de pré-calculés et peuvent donc s’appliquer à différents types de scénarios réels sur un réseau routier.                En dépit de la diversité des possibilités offertes par pgRouting pour l’analyse d’un réseau, et malgré d’intenses recherches sur la toile, il existe peu d’exemples sur la réalisation d’isochrones parfaites ou sur la gestion multimodale des coûts (vitesses/temps) sur un réseau.

La plupart des tutoriels ou des recherches existantes se concentrent essentiellement sur sa facilité d’utilisation, sur la recherche du plus court chemin qui elle-même est basé sur la notion de coût.

Dans ce manuel, nous allons apprendre les différentes étapes de modélisation d’un réseau routier. Pour ce faire j’ai fait le choix de travailler uniquement avec des données gratuites d’OpenStreetMap au lieu des données d’IGN, TéléAtlas et autres qui sont payantes.

Nous allons donc voir aussi comment installer pgRouting, importer les données d’OpenStreetMap , créer un graphe planaire orienté et l’exploiter dans Qgis. Puis nous allons modéliser les coûts (voiture, peton, vélo et distance) en respectant le sens de la circulation. Je vous montrerai également comment créer des zones de chalandises ou isochrones.

        I.               Approche méthodologique

Ceci est un guide méthodologique simplifié pour vous aider à démarrer pas à pas.

1.     Prérequis

Nous allons travailler sur Windows. Donc vous aurez besoins de Postgre/Postgis,  d’un éditeur de texte et d’un serveur Apache/PHP pour le déploiement Web. Pour ceux qui sont sur linux je vous invite à lire la documentation.

2.     Installation de PostgreSQL/PostGIS & pgRouting

Voici les instructions rapides pour l’installation PostgreSQL 9.5, PostGIS 2.0 et pgRouting 2.0. Pour PostgreSQL :

Téléchargez le programme d’installation d’enterprisedb.com ici.

Exécutez le programme d’installation. Vous aurez à choisir un login et mot de passe du super-utilisateur. A la fin du processus d’installation poursuivez avec le module Stack Builder.

Dans Stack Builder, sélectionnez l’installation PostgreSQL 9.5 et installer PostGIS dans la liste des extensions disponibles.    La procédure d’installation de PostGIS vous demandera le mot de passe super-utilisateur, le même que vous avez choisi lors de l’installation.

Puis je vous suggère de laisser le programme d’installation créer la base de données PostGIS « postgis_21_sample»  par défaut que nous aurons besoin plus tard comme modèle.

Contrairement à d’autre version de PostgreSQL, la version 9.5 installe directement l’extension pgRouting et pleins d’autres.

Pour plus de détails sur l’installation je vous invite à consulter la documentation détaillée sur le site de learnOSM.

3.     Création de la base de données « routing »

A la fin de l’installation, ouvrez pgAdmin (interface graphique d’administration de PostgreSQL) et connectez-vous avec votre login et mot de passe.

Vous  allez  désormais    créer    une    nouvelle    base  de données  incluant    les    fonctionnalités    pgRouting    en    utilisant la base de données par défaut « postgis_21_sample» comme modèle. Allez sur le menu et cliquez sur « SQL » et exécutez la requête suivante:

–PostGIS :
CREATE EXTENSION postgis ;
–Postgis_topology:
CREATE EXTENSION postgis_topology;
–pgRouting:CREATE EXTENSION pgrouting;

4.     Importation et intégration de données :
a-Installation d’oms2pgsql
Pour importer des données OpenStreetMap dans notre base de données, nous allons installer quelques outils. L’outil principal est appelé osm2pgsql, qui est un utilitaire qui charge les données XML d’OpenStreetMap en un format que nous pouvons stocker dans la base de données. L’utilitaire est facile à installer. Il suffit de le télécharger et d’ajouter son emplacement à notre variable d’environnement.

  • Pour télécharger la version Windows de osm2pgsql, cliquez sur : http://customdebug.com/osm/osm2pgsql.zip
  • Décompressez le fichier dans « C:\ ». Assurez-vous de créer un dossier permanent, parce que nous allons ajouter cet emplacement aux variables d’environnement.
  • Dans le dossier d’osm2pgsql que vous avez décompressé se trouve le fichier osm2pgsql.exe. C’est un programme que nous allons exécuter pour importer les données. Mais pour que le système Windows puisse le trouver, nous devons ajouter son emplacement aux variables d’environnement. Cliquez sur le menu Démarrer de Windows et tapez « chemin d’accès système ». Vous devriez voir une option « Modifier les variables d’environnement du système. » Cliquez sur l’image.
  • Cliquez sur le bouton nommé « Variables d’environnement ».
  • Trouvez la variable nommée « Path » et cliquez sur « Modifier… » pour ajouter le répertoire où se trouve osm2pgsql.exe.osm2pgsql1.png

Pour vérifier votre installation allez sur le menu démarrer et entrer « cmd » puis dans le terminal entrer « osm2pgsql ». Vous devrez voir ceci :

osm2pgsql2.png

b-      Les données   Ø  Le réseau de voirie

Comme stipulé précédemment, nous allons travailler uniquement avec les données d’OpenStreetMap. Les données sont gratuitement disponibles sur le site de geofabrik.

Etant donné qu’on va travailler sur la commune de Nantes, allez sur « download » et choisissiez Europe, la France puis « Pays de la Loire ». Sinon voici le lien direct de téléchargement.http://download.geofabrik.de/europe/france/pays-de-la-loire.html

Il y a plusieurs formats disponibles mais je vous conseille « pays-de-la-loire-latest.osm.bz2 » car ce dernier est plus complet.  Ø  La commune de Nantes               Importer les communes dans le schéma public. Vous pouvez le télécharger depuis le site de l’open data.https://www.data.gouv.fr/fr/datasets/decoupage-administratif-communal-francais-issu-d-openstreetmap/.

Vous pouvez faire un glissez déposer avec QGIS ou avec OGR

c-      Intégration des données 

Ø  La voirie

Il existe plusieurs façons de le faire évidemment (cf. ETL, shp2pgsql, ogr2ogr ou QGIS) mais nous allons passer par « osm2pgsql », un outil de chargement des données OSM dans PostgreSQL) qui est le plus rapide à ma connaissance pour ces données-là. Pour ce faire voici la commande à passer après son installation :

osm2pgsql3.png

osm2pgsql -v -d routing -U postgres -H localhost -P 5432 -S "C:\osm2pgsql\default.style" C:\osm2pgsql\donnees\pays-de-la-loire-latest.osm\pays-de-la-loire-latest.osm.bz2

Définition

-v mode verbeux-H le host-P le port-U : votre login -W votre mot de passe-d la base de données-S

le style par default téléchargeable ici :https://github.com/openstreetmap/osm2pgsql/blob/master/default.style

Voici le résultat dans le schéma public4 tables sont créées : planet_osm_line => pour les routes détaillées que vous allez utiliser: planet_osm_route : pour les routes simplifiées; planet_osm_point : Adresses,  POI…; planet_osm_polygon : Occupation du sol, bâtis …

Ø  La commune de Nantes

Dans la terminale de Windows (cmd) tapez

"C:\Program Files\PostgreSQL\9.5\bin\shp2pgsql" –V –H localhost –P 5432 –U admin -d -D -I -W CP1252 -s 2154    "C:\Script_Abdou\sda\departement.shp"  public.communes | "C:\Program Files\PostgreSQL\9.5\bin\psql" -d routing

Définiton : -v mode verbeux,-H le host,-P le port,-U : votre login -s la projection, CP1252 l’encodage en Latin et –I pour indexer les données-d la base de données  Pour la description des commandes je vous invite à consulter la documentation ici :http://www.volkerschatz.com/net/osm/osm2pgsql-usage.html

d-      Visualisation des données avec QGIS

Maintenant que vous avez créé votre base de données et y importé les données. Nous allons les visualiser avec QGIS. Pour cela créer une connexion avec QGIS et visualiser les résultats. Je vous invite à consulter la documentation de learnOSM.

Ø  Paramétrez votre de connexion

qgis connexion.png

Ø  Ajouter et explorer les données « planet_osm_line »

Alors comme vous pouvez le constater, nous avons tous le réseau de voirie du Pays de la Loire. Mais ce qui nous intéresse c’est uniquement le réseau routier sans les voies ferrées, le découpage administratif, les parkings, l’hydrographie  etc…  Pour ce faire vous devez filtrer sur « « highway » is not null ». La colonne « highway » classifie les types de routes et est toujours remplie.

qgis connexion 3

Vous pouvez explorer les données et pour plus d’informations détaillées sur la typologie du réseau et des objets graphiques, je vous invite à consulter le WIKI en français ici :http://wiki.openstreetmap.org/wiki/FR:Élements_cartographiques

5.    Création d’un graphe la topologique sur le réseau

a-      Structuration et modélisation topologique

Dans un réseau topologique, les objets sont définis par leurs relations plutôt que par leurs géométries. Il est supposé avoir ses lignes (arrêtes) connectés à des points uniques (nœuds). Je vous invite à revoir la théorie des graphes pour plus d’infos complémentaire.

Nous allons donc construire un réseau topologique avec nos données routières. Cela signifie que pour toute arête dans nos données routières, les terminaisons de chaque arête vont être connectées en un nœud unique et aux autres arêtes qui sont aussi connectées à ce même nœud. Une fois que les arêtes sont connectées aux nœuds nous aurons un graphe qui pourra être utilisé pour le routage avec pgRouting. Pour cela voici la commande à utiliser :

Ø  Création du schéma topologique en Lambert  93

SELECT topology.CreateTopology(‘routing’, 2154);

Voir pgr_createTopology pour plus d’informations. Cette commande va vous créer le schéma topologique « routing » avec différentes tables notamment les points de départ et d’arrivées pour chaque tronçon « start_node et end_node ».

Maintenant il vous faut importer vos données dans le schéma avec la projection en Lambert 93 :

Ø  Ajouter une colonne géometrique  en Lambert 93 dans la table « planet_osm_line »

ALTER TABLE  public.planet_osm_line ADD COLUMN geom geometry (LINESTRING,2154);

Ø  Mettre à jour la colonne géométrique

UPDATE  public.planet_osm_line SET geom=st_transform(way,2154)

Ø  Création de la table route qui intersecte la commune de Nantes

DROP TABLE routing.route ;
CREATE TABLE routing.route as SELECT osm_id, access, "addr:housename", "addr:housenumber", "addr:interpolation",        admin_level, aerialway, aeroway, amenity, area, barrier, bicycle,        brand, bridge, boundary, building, construction, covered, culvert,        cutting, denomination, disused, embankment, foot, "generator:source",        harbour, highway, historic, horse, intermittent, junction, landuse,        layer, leisure, lock, man_made, military, motorcar, name, "natural",        office, oneway, operator, place, population, power, power_source,        public_transport, railway, ref, religion, route, service, shop,        sport, surface, toll, tourism, "tower:type", tracktype, tunnel,        water, waterway, wetland, width, wood, z_order, st_intersection( a.geom,b.geom) as geom  FROM public.planet_osm_line a, public.communes bWHERE ST_intersects(a.geom,b.geom)AND codgeo=’44109’AND highway IS NOT NULL; 

Ø  Creer un index spatial
 CREATE INDEX idx_ routing_route_geom  ON routing.route  USING gist  (geom);

Comme dans QGIS « highway is not null »  permet de filtrer uniquement sur la voirie. Prévoir une pause de 1mn pour Nantes si vous avez un pc qui tourne bien. Sinon faite le test sur une petite commune.

b–      Créer et associer la topologie PostGIS

Cette commande vous permettra de couper les tronçons à chaque intersection. C’est-à-dire la création d’un graphe planaire avec des nœuds et des arcs.

Ø  Ajouter une colonne topologique ‘topo’

SELECT topology.AddTopoGeometryColumn('routing’, 'routing', 'route', 'topo', 'LINESTRING');

Ø  Convertir les lignes brisées en noeuds et arrêtes au sein de la topologie 

UPDATE routing.route SET topo = topology.toTopoGeom(geom, 'routing’, 1, 0.0001);

Ou « 0.0001 » est le seuil de tolérance du réseau à adapter.

À partir de maintenant, nous avons un réseau topologique, où les imperfections ont été corrigées. On a fusionné les jonctions imparfaites.              

Vous pouvez également vérifier votre réseau avec les fonctions d’Analyse de graphe.

Ø  Maintenant nous allons enrichir notre table en rajoutant les noms et  les type de voiries.

ALTER TABLE routing.edge_data  add COLUMN  type   character varying;
ALTER TABLE routing.edge_data  add COLUMN  nom  character varying;

Ø  Créer une table intermédiaire « route_edge_data » pour stocker les noms et types.

drop table routing. route_edge_data;
create table routing.route_edge_data as
SELECT e.edge_id,r.highway as type , r.oneway,r.name as nom  FROM routing.edge e,     routing.relation rel,     routing.route rWHERE e.edge_id = rel.element_id  AND rel.topogeo_id = (r.topo).id;

Ø  Mettre à jour attributs depuis la table « route_edge_data »

UPDATE routing.edge_data  a SET (type, nom) = (r.type, r.nom)FROM routing.route_edge_data  rWHERE a.edge_id=r.edge_id;

6. Les algorithmes de pgRouting

Plusieurs algorithmes sont utilisés pour parvenir aux résultats attendus. j’en décrirai que deux (2) :

  • L’algorithme de Dijkstra :

Dans l’étude des graphes, l’algorithme de Dijkstra est utilisé pour retrouver le chemin le moins coûteux. Cet algorithme permet de calculer le chemin le plus court pour se rendre dans une ville, d’un point à un  autre ou d’une ville à une autre quand d’une région. Le  réseau routier devra être modélisé en un graphe connexe – un graphe non orienté sur lequel il existe un chemin de tout point vers un autre avec des poids ou valeurs pouvant être des distances (pour le calcul de plus chemin), le temps de parcourt estimé (pour le chemin le plus rapide), la consommation de carburant et les frais de routes (pour le chemin le plus économique)   affectés à chaque arc – dont le poids lié à un arc est positif ou nul.

Du nom de son inventeur, le mathématicien et informaticien néerlandais Edsger Dijkstra, l’algorithme a été publié en 1959.

  • L’algorithme A* :

L’algorithme de recherche A* (qui se prononce A étoile, ou A star à l’anglaise) sert également dans la recherche de chemin optimal dans un graphe, entre un nœud initial et un nœud final. Il est très simple d’usage et est souvent présenté en exemple typique d’algorithme de planification dans le domaine de l’intelligence artificielle. Cet algorithme a été élaboré par Peter E. Hart, Nils Nilsson et Bertram Raphel en 1968, pour que la solution trouvée en premier lieu dans les calculs, soit l’une des meilleurs possibles. C’est ce qui explique qu’il soit célèbre dans des applications privilégiant la durée du calcul par rapport à la vraisemblance des résultats.

Un système de routage a plusieurs avantages :

  • La réduction des coûts du transport car le chemin parcouru est celui le plus rentable ;
  • l’augmentation de la fluidité du trafic, et donc gain de temps et aussi et surtout ;
  • la réduction considérable des émissions de CO2 (car à chemin minimal, énergie ou carburant minimal) destructeur de l’environnement si cher à l’espèce humaine à l’heure où tous les efforts sont canalisés vers la sauvegarde et la protection de cette précieuse ressource.

Une fois la modélisation du graphe terminé, vous allez pouvoir tester votre réseau. Pour cela voici une rapide présentation de quelques fonctionnalités de pgRouting que nous pouvons utiliser.

  • Autres fonctions

Pour plus d’information voici la liste essentielle des fonctions de pgRouting :

Fonctions Nom Traduction
pgr_apspJohnson All Pairs Shortest Path, Johnson’s Algorithm Retourne tous les coûts pour chaque paire de nœuds dans le graphe
pgr_apspWarshall All Pairs Shortest Path, Floyd-Warshall Algorithm Retourne tous les coûts pour chaque paire de nœuds dans le graphe
pgr_astar Shortest Path A* Retourne le plus court chemin en utilisant l’algorithme A*.
pgr_bdAstar Bi-directional A* Shortest Path Retourne le plus court chemin en utilisant l’algorithme bidirectionnel A*.
pgr_bdDijkstra Bi-directional Dijkstra Shortest Path Retourne le plus court chemin en utilisant l’algorithme bidirectionnel Dijkstra
pgr_dijkstra Shortest Path Dijkstra Retourne le plus court chemin en utilisant l’algorithme Dijkstra
pgr_kDijkstra destination Shortest Path Mutliple pgr_kdijkstraCost Retourne les coûts pour les K plus courts chemins utilisant l’algoithme Dijkstra.
pgr_kdijkstraPath – Retourne les chemins pour les K plus courts chemins en utilisant l’algorithme Dijkstra
pgr_ksp K-Shortest Path Retourne les “K” plus courts chemins.
pgr_tsp Traveling Sales Person pgr_tsp – Retourne la meilleure route à partir d’un point de départ via une liste de nœuds.
pgr_tsp – Retourne le meilleur ordre de route quand passée une matrice de distance.
pgr_makeDistanceMatrix – Retourne une matrice de distance Euclidienne à partir de points fournis par le résultat sql.
pgr_trsp Turn Restriction Shortest Path (TRSP) Retourne le plus court chemin avec support pour les restrictions de virage.
pgr_drivingDistance Driving Distance Retourne la distance de conduite à partir d’un point de départ.

Ø  Calcul du plus court chemin avec Dijkstra.

L’algorithme de Dijkstra a été la première implémentation disponible dans pgRouting. Il ne nécessite pas d’autres attributs que les champs source et target et les attributs de cost (coût). Il peut être utilisé sur des graphes orientés ou non. Vous pouvez spécifier que votre réseau à un coût de parcours inverse (reverse cost) ou non. Pour être en mesure d’utiliser un coût de parcours inverse, vous devez ajouter une colonne de reverse_cost. Dans un premier temps, nous allons utiliser que la distance des tronçons de voirie comme coûts pour calculer le chemin le plus court. On reviendra sur les notions de coûts inverses pour respecter le sens de la circulation.

Le coût pour les distances c’est très simple. C’est la longueur de la route en mètre. Voici la requête :

Ø  Ajout de la colonne  colonne tps_distance pour les distances

ALTER TABLE routing.edge_data  add COLUMN tps_distance    double precision;

Ø  Mise çà jour des distances

UPDATE routing.edge_data a  SET tps_distance=st_length(st_transform(geom,2154))/1000;

Ø  Visualisation du résulat avec Qgis

Si vous avez bien suivi le tutorial vous devriez avoir ça dans QGIS. La couche du réseau s’appelle « edge_data ».  Et quand vous cliquez sur ce tronçon, vous devriez avoir les mêmes attributs. Il est de type « cycleway ». Le nom de la route c’est « le chemin de halage de la Loire ». La longueur de son tronçon est de 1.898 km.

nantes.png

        III.       Calcul du plus court chemin avec Dijkstra

Ø  Calcul du plus court chemin  en isodistance

Important ! Les attributs source (start_node) et target (end_node) sont vos points de départs et d’arrivées On n’aura pas les mêmes. Et tps_distance votre coût estimé ici en distance (en mètre).

Par exemple qu’elle est la distance la plus courte entre le  centre commerciale de Beaulieu et le château de Nantes. Voici la requête à passer:

Ø  le plus court chemin

SELECT seq, id1 AS node, id2 AS edge, cost,type, geom FROM pgr_dijkstra('SELECT edge_id as id, start_node as source, end_node as target,  tps_distance as  cost from routing.edge_data',            1187, 8168, false, false)             as diJOIN routing.edge_data ptON di.id2 = pt.edge_id;

Evidemment vous allez me dire mais comment je vais faire pour connaitre les nœuds de départ et d’arrivées ? Pas de soucis je les ai identifiés visuellement sous QGIS. Il faut également retenir que vos attributs peuvent être différents des miens. Dans Qgis importez les tables « node, edge_data », sur un fond OSM  ou Google cherchez les nœuds les plus près de ces lieux et exécuter la requête.

requete.png

Pour résumer on a le nœud de départ « 1187» du centre commercial  qui va jusqu’au nœud « 1439» sur une distance de 0.150 km soit 150 m   et prendre le nœud 939  sur 406m et ainsi de suite jusqu’à atteindre le nœud le plus près du château de Nantes.

Pour vous faciliter la tâche je vais vous épargner ça on va allez directement sur QGIS pour s’amuser avec notre réseau.

Ø  Calcul du plus court chemin avec pgRouting Layer de Qgis

Parmi les nombreuses extensions que compte QGIS il y en a un spécialement conçu pour pgRouting et c’est « pgRouting Layer ». Allez dans QGIS, dans le menu « Extension », cliquez sur installez les extensions. Et dans la barre de recherche tapez sur pgr. Vous allez voir pgRouting Layer. Cliquer sur installer. 

pgrouting layer.png

Importez la table du réseau routier « edge_data » et activer pgRouting Layer depuis le menu « Base de données » si ce n’est pas fait. Je vous invite à consulter la documentation QGIS pour l’importation sur depuis lien http://planet.qgis.org/planet/tag/pgrouting/.                A présent renseignez les champs comme suit et calculer votre itinéraire.

pgrouting layer distance 1.png

Comme vous pouvez le voir c’est rapide et interactif. Et n’hésitez pas à exporter les résultats. Pour vérifier votre requête faites un clic droit et cliquez sur propriété. Et dans « Général » vous pouvez récupérer la requête et le passez dans pgAdmin. Elle sera la même que celle  qu’on avait formulé avant à quelques détails prés.

Ø  Le plus court chemin à pied

D’après ma recherche bibliographique, il est généralement établit  qu’un piéton marche en moyenne à une vitesse de 4 à 5km/h. Pour estimer te temps de parcours pour un tronçon donné, la règle est la suivante : Temps = distance /vitesse  sachant qu’on connait déjà la distance des tronçons  (tps_distance).

Ø  Donc on ajoute l’attribut tps_pieton dans notre table « edge_data »  

ALTER TABLE routing.edge_data  ADD COLUMN tps_pieton   double precision; Ø  MAJ tps_pieton   
UPDATE routing.edge_data  SET tps_pieton =tps_distance/5;

Cependant on sait aussi qu’un piéton ne doit pas emprunter les autoroutes et voix expresses donc nous allons les restreindre. Mais il faut également avoir en tête  qu’il peut circuler dans les deux sens.

Ø  MAJ restreinte tps_pieton  

UPDATE routing.edge_data  SET tps_pieton =-1 WHERE type IN ('trunk','trunk_link','motorway','motorway_link')

Ø  Le plus court chemin à vélo

En moyenne, la vitesse des cyclistes est estimée entre 13 et 15 km selon les études. Pour aller plus loin nous allons calculer la vitesse moyenne des cyclistes en fonction des types de routes empruntées.

Pour ce faire après d’intense recherche bibliographique, je suis tombé sur un article intéressant qui m’a conduit à un projet Gituhb :https://github.com/Project-OSRM/osrm-backend/tree/develop/profiles

Vous y trouverez donc les vitesses en fonction des types de routes que j’ai légèrement adaptées.

Type Vitesse
cycleway 15
primary 15
primary_link 15
secondary 15
secondary_link 15
tertiary 15
tertiary_link 15
residential 15
unclassified 15
living_street 15
road 15
service 15
track 15
path 15
footway 12
pedestrian 12
Steps 12
Autres (trunk,trunk_link,motorway,motorway_link…) -1

Ø  Donc on ajoute l’attribut tps_velo dans notre table « edge_data »

ALTER TABLE routing.edge_data  ADD COLUMN tps_velo double precision;

Ø  MAJ tps_velo

UPDATE routing.edge_data SET tps_velo =tps_distance /15 ;
UPDATE routing.edge_data SET tps_velo=tps_distance /12 WHERE type IN ('footway','pedestrian' ) ;
UPDATE routing.edge_data SET tps_velo =tps_distance /2 WHERE type ='steps' ;

Ø  Restrictions

UPDATE routing.edge_data SET tps_velo  =-1 WHERE tps_velo  IS NULL ;

Voici le résultat comparé avec Google. Vous pouvez constater aussi que le sens de la circulation n’a pas été pris en compte. On y reviendra plus tard.

Dijkstra velo.png

Dijkstra velo google.png

 Ø  Le plus court chemin en voiture

  Concernant OSM voici une documentation complètes sur les limitations de vitesses :http://wiki.openstreetmap.org/wiki/OSM_tags_for_routing/Maxspeed.

Pour affecter le temps passé par tronçon nous allons reprendre la documentation de Github. Les vitesses je les ai optimisées pour simuler les conditions de circulation en temps réelle.

Type Vitesse
motorway 90
motorway_link 45
trunk 85
trunk_link 40
primary 65
primary_link 30
secondary 55
secondary_link 25
tertiary 40
tertiary_link 20
residential 25
unclassified 25
living_street 10
road 25
Autres (path, footway, pedestrian, service, Steps…) -1

Pour ce faire on ajoute l’attribut tps_voiture dans notre table « edge_data »

  Ø  Couts tps_voiture

ALTER TABLE routing.edge_data  ADD COLUMN tps_voiture double precision;

Ø  MAJ tps_voiture

update routing.edge_data set tps_voiture  =tps_distance /90 where type =  'trunk' ;
;update routing.edge_data set tps_voiture  =tps_distance /45 where type = 'trunk_link';  
update routing.edge_data set tps_voiture  =tps_distance /85 where type ='motorway';
 update routing.edge_data set tps_voiture  =tps_distance /40 where type = 'motorway_link'  ;
update routing.edge_data set tps_voiture  =tps_distance /65 where type ='primary'  ;
update routing.edge_data set tps_voiture  =tps_distance /30 where type = 'primary_link' 
;update routing.edge_data set tps_voiture  =tps_distance /55 where type = 'secondary'  ;
update routing.edge_data set tps_voiture  =tps_distance /25 where type = 'secondary_link' ;
update routing.edge_data set tps_voiture  =tps_distance /40 where type = 'tertiary' ;
update routing.edge_data set tps_voiture  =tps_distance /20 where type = 'tertiary_link' ;
update routing.edge_data set tps_voiture  =tps_distance /25 where type = 'residential' ;
update routing.edge_data set tps_voiture  =tps_distance /25 where type = 'road'  ;
update routing.edge_data set tps_voiture  =tps_distance /25 where type = 'unclassified';
update routing.edge_data set tps_voiture  =tps_distance /10 where type = 'living_street';

 Ø  Coûts inverses

update routing.edge_data set tps_voiture  =-1 WHERE tps_voiture IS NULL ;

Voici les résultats entre mon lieu de travail et mon domicile. Vous noterez également que le sens de la circulation n’est pas encore pris en compte dans notre modèle.

Dijkstra voiture direction.png

Dijkstra voiture google.png

   

     IV. Calcul du plus court chemin avec les coûts inverses
En complément du cost (coût de la traversée du tronçon dans le sens de numérisation), le reverse_cost (coût inverse) sert essentiellement pour les graphes orientés, c’est-à-dire possédant un champ (comme sens de la circulation) indiquant une différence de parcours selon ce champ. Si ce n’est pas le cas, alors les 2 champs  cost et reverse_cost prennent la même valeur. Dans la plupart des cas, la valeur choisie pour représenter ces coûts est la longueur des tronçons. C’est une valeur que toute base de données possède, soit directement en tant qu’attribut, soit en champ calculable avec la géométrie des tronçons. Néanmoins, il est possible de choisir d’autres valeurs (en fonction des données), d’effectuer des règles de calculs ou de pondérer, etc. Voici la règle à respecter :

   1. Dans le cas d’un tronçon de route à double sens (aucune restriction de circulation), les  cost et reverse_cost sont les mêmes.

   2. Dans le cas d’un tronçon de route à sens unique (information remontée par le champ « sens par exemple », et relatif au sens de numérisation), le parcours en sens contraire au sens de circulation doit être empêché.

Par convention, il a été choisi de remplir avec la valeur suffisamment importante (1 000 000) ou une valeur négative (-1).

Calcul de coût entre 2 points

calcul entre 2 points.png

  1. Cas pratique :

               Calculons le cout du trajet aller et retour entre les 2 points A et D.

 graphe orienté.png

Exemple de graphe orienté

Dans le cas présent, tous les tronçons sont à double-sens, sauf le tronçon « BD » (circulation impossible de B vers D). Si l’on applique la règle de calcul précédente cela donne :

Coûts et parcours pour un trajet aller retour entre 2 points

Coûts et parcours pour un trajet aller/retour entre 2 points

2.      Applications sur notre réseau de voirie

Pour orienter notre graphe, nous allons rajouter le sens de circulation la route. Avec les données OSM, le sens de la circulation se trouve dans la colonne « oneway » que nous allons récupérer dans un champ qu’on va appeler  « sens ».

Pour enrichir notre table en rajoutant le sens de la circulation voici la requête.

drop table routing.route_edge_data;
create table routing.route_edge_data asSELECT e.edge_id,r.oneway as sens  FROM routing.edge e,     routing.relation rel,     routing.route rWHERE e.edge_id = rel.element_id  AND rel.topogeo_id = (r.topo).id;

Ø  Ajouter les colonnes à remplir

ALTER TABLE routing.edge_data  add COLUMN  sens character varying;

Ø  Maj attributs depuis la table route

UPDATE routing.edge_data  a SET (sens) = (r. sens)FROM routing.route_edge_data  rWHERE a.edge_id=r.edge_id;

Avec ça vous allez pouvoir affecter des couts inverses. Pour ce faire nous allons aussi modifier les coûts directs. Après une documentation approfondie voici les règles à suivre pour respecter le sens de la circulation dans OSM.

Régle 1 Quand  » oneway  » = ‘yes’ or « oneway » = ‘1’ alors   cost = distance/vitesseQuand  » oneway  » = ‘yes’ or « oneway » = ‘1’ alors   reverse_cost= -1

 Régle 2 Quand  » oneway  » = ‘no’   alors   cost = distance/vitesseQuand  » oneway  » = ‘no’   alors   reverse_cost= distance/vitesse

Régle 3 Quand  » oneway  » = ‘is null’ or « oneway » = ‘no;yes’ alors   cost = distance/vitesseQuand  » oneway  » = ‘is null’ or « oneway » = ‘no;yes’ alors   reverse_cost= -1

Régle 4 Interdiction de circuler autres que les véhiculesQuand le type n’est pas réservé pour les voiture, cost et reverse cost = ‘-1’

Voici les requêtes à passer pour attribuer les coûts pour chaque type de tronçons. D’abord on ajoute de la colonne  rc_voiture  pour les coûts inverses. Pour info les coûts inverses vont nous servir ici pour prendre en compte le sens de la circulation.

ALTER TABLE routing.edge_data  add COLUMN rc_voiture  double precision;
Puis executez la requete en respectant les regles de sens de circulations

Exemple de l’itinéraire aller-retour  entre  mon domicile et mon lieu de  travail en voiture.

Dijkstra voiture direction 1.png

Dijkstra voiture direction 2.png

Et avec Google

Dijkstra voiture direction g1.png

Dijkstra voiture direction g2.png

V.Calcul du plus court chemin avec L’algorithme A* :

L’algorithme de recherche A* (qui se prononce A étoile, ou A star à l’anglaise) sert également dans la recherche de chemin optimal dans un graphe, entre un nœud initial et un nœud final.

               On va   préparer les colonnes pour le plus court chemin en utilisant l’algorithme A*.

—Ajout des colonnes

ALTER TABLE routing.edge_data ADD COLUMN x1 double precision;
ALTER TABLE routing.edge_data  ADD COLUMN y1 double precision;
ALTER TABLE routing.edge_data  ADD COLUMN x2 double precision;
ALTER TABLE routing.edge_data  ADD COLUMN y2 double precision;

–Mise à jour des colonnes A*

UPDATE routing.edge_data  SET x1 = ST_x(ST_PointN(geom, 1));
UPDATE routing.edge_data  SET y1 = ST_y(ST_PointN(geom, 1));
UPDATE routing.edge_data  SET x2 = ST_x(ST_PointN(geom, ST_NumPoints(geom)));
UPDATE routing.edge_data  SET y2 = ST_y(ST_PointN(geom, ST_NumPoints(geom)));

 Voici un exemple d’utilisation. Vous pouvez également le tester directement dans Qgis.

SELECT seq, id1 AS node, id2 AS edge, cost,geom  FROM pgr_astar('              
  SELECT edge_id as id,start_node as source , end_node as target,tps_distance as cost, x1, y1, x2, y2  FROM osm_topo.edge_data ', 30, 60, false, false) a  , osm_topo.edge_data b  where a. id1=b.start_node ;  30, 60, false, false);

        VI.               L’ACCESSIBILITE MULTIMODALE EN ISOCHRONES ET ISODISTANCES

isochrones.jpg

En cartographie, une courbe isochrone est une courbe géométrique délimitant les points accessibles par un véhicule – terrestre ou aérien – en un temps donné (par exemple, la zone pouvant être desservie en moins de 30 minutes ou 5km (isodistance) par un livreur de pizza ou un dépanneur de matériel informatique).

1.      ObjectifsØ  Créer des isochrones et isodistances à partir du réseau  pour les piétons, cyclistes et automobilistes. Ø  Calculer les profils sociodémographiques (emplois, populations, logements, ménages..) par isochrones

2.      Analyse comparative des méthodes de création d’isochrones                              D’abord, nous allons calculer et visualiser les nœuds et arcs qui se situent à 500m de Nantes métropole. Pour ce faire nous allons utiliser l’algorithme driving distance de pgRouting.

Ø  Driving distance pour les nœuds

--Drop table routing.pgr_drivingdistance_pt;
 create table routing_nantes.pgr_drivingdistance_pt as       
 SELECT   id1 AS node_id ,  cost,geomFROM pgr_drivingdistance(
'SELECT edge_id as id, start_node as source, end_node as target, tps_distance as    cost from routing.edge_data',12616, 0.5, false, false) as diJOIN routing.node ptON di.id1 = pt.node_id;

Ø  Driving distance pour les tronçons

--Drop table routing.pgr_drivingdistance_lgn;
create table routing_nantes.pgr_drivingdistance_lgn as         
SELECT   id1 AS node_id ,pt.tps_distance as cost,geomFROM pgr_drivingdistance(
'SELECT edge_id as id, start_node as source, end_node as target,  tps_distance as  cost from routing.edge_data',12616, 0.5, false, false) as diJOIN routing.edge_data ptON di.id1 = pt.start_node;

Ø  Résultats

arcs et noeuds avec pgr_drivingdistance.jpg

 

Comme vous pouvez le constater les deux résultats sont légèrement différents. ·         Les nœuds retournent exactement la distance entre 0 et 499m ·

Tandis que les arcs dépassent de plusieurs mètres sur les extrémités.  Cette erreur est liée aux intersections entre les nœuds extrêmes et les arcs au-delà des 500m. Maintenant nous allons créer des isochrones à partir de ces deux objets et voir lequel reflète mieux la réalité.

a-      Isochrone avec pgr_alphAShape               Cette fonction retourne un tableau avec des lignes (x, y) qui décrivent les sommets d’une forme alpha. Nous allons donc le coupler avec l’algorithme de driving distance pour créer un polygone à partir des sommets des nœuds et des arcs. Pour plus d’information je vous renvoie vers cette documentation: https://anitagraser.com/2011/09/25/a-closer-look-at-alpha-shapes-in-pgrouting/

Suite….

Publicités

13 réflexions sur “Postgis & Pgrouting

  1. mld dit :

    Bonjour,
    Je bloque depuis quelques semaines maintenant sans trouvé de réponse à mon problème…

    Lors de la création de l’extension pgrouting pgadmin me retourne l’erreur suivante :
    ERROR: could not load library « C:/PROGRA~1/POSTGR~1/pg96/../pg96/lib/postgresql/libpgrouting-2.4.dll »: unknown error 126

    Ok soucis pour postgis et postgis_topology.

    D’où pourrait venir le problème ?

    PostgreSQL : 9.6.0 on x86_64-pc-mingw64, compiled by gcc.exe (Rev5, Built by MSYS2 project) 4.9.2, 64-bit
    Postgis : 2.3.0rc1
    Windows 10 64bits

    pgrouting-pg96-binaries-2.4.0-rc1w64gcc48.zip (téléchargé ici : http://winnie.postgis.net/download/windows/pg96/buildbot/)

    J'aime

  2. mld dit :

    Merci pour le tuyau, si je ne trouve pas de solution je crois bien qu’on va finir par faire ça…
    Pour un info on un a Postgresql sous distribution BigSQL sur un serveur Windows.

    J'aime

  3. Bonjour,
    Je souhaite réaliser des isochrones d’après la méthode dite ‘Isochrone avec St_MakePolygon et St_ExteriorRing’
    La requete ci dessous permet d’obtenir l’isochrone pour toutes les routes à moins de 500m de la poste. Elle fonctionne dès lors que l’on connait l’identifiant (edge_id) du point qui nous intéresse.

    Code:

    with buffer_itineraire as (
    SELECT et.edge_id,
    st_buffer(et.geom,10,’endcap=square join=round’) as geom
    , 1 as factor
    FROM pgr_drivingdistance(
    ‘SELECT edge_id as id, start_node as source, end_node as target, tps_distance as cost from routing.edge_data’,
    8215, 0.5, false, false
    ) firstPath
    CROSS JOIN
    (SELECT id1,cost from pgr_drivingDistance(
    ‘SELECT edge_id as id,start_node as source , end_node as target,tps_distance as cost FROM routing.edge_data’,
    8215,0.2,false,false)
    ) secondPath
    INNER JOIN routing.edge_data et
    ON firstPath.id1 = et.start_node
    AND secondPath.id1 = et.end_node)
    SELECT ST_MakePolygon(ST_ExteriorRing(ST_GeometryN(st_union(geom),1))) as geom FROM buffer_itineraire;

    Je cherche à exécuter cette requete pour tous les points qui ont l’information ‘Agroalimentaire’, puis ‘Services’, puis ‘Autres’ contenue dans un champ ‘type_comm’, sans avoir à chercher manuellement l’edge_id pour chacun des points question (il y en a 37).

    Pourriez vous m’aider à résoudre cette requête?

    Dans l’attente de vos réponses, questions, je vous souhaite une bonne journée.

    Marine.

    J'aime

Laisser un commentaire

Entrez vos coordonnées ci-dessous ou cliquez sur une icône pour vous connecter:

Logo WordPress.com

Vous commentez à l'aide de votre compte WordPress.com. Déconnexion / Changer )

Image Twitter

Vous commentez à l'aide de votre compte Twitter. Déconnexion / Changer )

Photo Facebook

Vous commentez à l'aide de votre compte Facebook. Déconnexion / Changer )

Photo Google+

Vous commentez à l'aide de votre compte Google+. Déconnexion / Changer )

Connexion à %s