Pasq.fr

Parce-qu'il y a forcément du sens à tout ce bordel !

Geocodage directement dans POSTGIS

Rédigé par Alain 6 commentaires
c'est ici sur la carte

Comment utiliser les API de la BAN (Base d'adresse Nationale) directement dans une base POSTGIS et géocoder vos adresses sans y penser ?

C'est quand je suis tombé sur cet article (https://www.crunchydata.com/blog/geocoding-with-web-apis-in-postgres) que je me suis dit qu'il fallait le faire avec notre chère BAN. Alors court, intense et fortement utile pour géocoder directement dans la base postgis grâce à une fonction.

Ceci est possible grâce au API fournies par le site adresse.gouv.fr et sa page d'aide que je vous invite à consulter pour connaître les spécifications de la recherche.

Prérequis

  • une base postgresql (nan, sans déconner...) > version 11 (ah, vous voyez),
  • l'extension plpython3u à installer avec, chez moi, j'ai du l'installer dans le pg_catalog, et je ne sais pas pourquoi... essayer sans SCHEMA d'abord et si ça marche pas, ajouter le :

CREATE EXTENSION plpython3u SCHEMA "pg_catalog"

Créer la "function"

Dans votre administration de base données préférée (fichier .sql à la fin) :

CREATE OR REPLACE FUNCTION geocode(address text)

RETURNS geometry

AS $$

import requests

try:

payload = {'q' : address}

base_geocode = 'https://api-adresse.data.gouv.fr/search/'

r = requests.get(base_geocode, params = payload)

coords = r.json()['features'][0]['geometry']['coordinates']

lon = coords[0]

lat = coords[1]

geom = f'SRID:4326;POINT({lon} {lat})'

except Exception as e:

geom = None

return geom

$$

LANGUAGE 'plpython3u';

Pour l'explication rapide, le code récupère les coordonnées à l'adresse de l'api, la réponse étant en .json, on extrait les données "coordinates" dans un tableau pour faire les coordonnées. Celles-ci sont servies en WGS84 (espg 4326), et on retourne le champs "geom" (si c'est le nom de votre champs geometry, ou des fois the_geom).

Utilisation

On peut tester cette fonction : 

SELECT ST_AsText(geocode('4 rue du bois d''amour brest'));

ce qui me retourne :  POINT(-4.486943 48.388149)

 J'ai choisi cette rue parce que 1 le nom est sympa et 2 nécessite un double ''

Il faut faire attention à cela et nettoyer la base autant que possible avant d'utiliser la fonction. Un regex peut sûrement améliorer la chose, faut que j'y pense !

On peut alimenter une table avec les points directement, avec par exemple :

INSERT INTO saisie (adresse, geom) VALUES('4 rue du bois d''amour brest', st_asewkt(geocode('4 rue du bois d''amour brest')));

Ce qui me crée bien :

adresse                    |geom                       |
---------------------------+---------------------------+
4 rue du bois d'amour brest|POINT (-4.486943 48.388149)|

point visible :

point placé

en faisant les tests, on peut taper des fois n'importe quoi : INSERT INTO saisie (adresse, geom) VALUES ('test débile',st_asewkt(geocode('pendach les oies 29100 ailleurs'))); Et bien, il me trouve une adresse quand même !!!

test débile

L’intérêt étant de l'utiliser une base constituée où vous avez des adresses à placer. Je vous laisse l'utilité pour vous.

Et le géocodage inverse, alors

Effectivement, si on peut l'un, on peut faire l'autre, en modifiant un peu le code de la fonction :

CREATE OR REPLACE FUNCTION geocode_reverse(lon text,lat text)

RETURNS TEXT

AS $$

import requests

try:

payload = {'lon':lon, 'lat':lat}

base_geocode = 'https://api-adresse.data.gouv.fr/reverse/'

r = requests.get(base_geocode, params = payload)

ad = r.json()['features'][0]['properties']['label']

except Exception as e:

ad = 'erreur'

return ad

$$

LANGUAGE 'plpython3u';

Le .json récupéré est le même format que pour la recherche par adresse, ce qui est très pratique. Ici, j'ai choisi de faire ressortir le libellé entier de l'adresse appelé "label", mais j'aurais pu sortir que la ville, ou le score de matching. alors testons :

SELECT geocode_reverse('2.37','48.357');

on obtient bien : 23 Grande Rue 91720 Prunay-sur-Essonne

En reprenant ma base, on pourrait très bien avoir ( j'ajoute une colonne pour le test) et je lance, sans oublier de caster mes coordonnées en TEXT :

UPDATE saisie SET "sortie reverse" = geocode_reverse(ST_X(geom)::text,ST_Y(geom)::text);
  adresse geom sortie reverse
  4 rue du bois d'amour brest POINT (-4.486943 48.388149) 4 rue du Bois d'Amour 29200 Brest
  test débile POINT (-4.259477 48.064154)  Bois du Juch 29100 Le Juch

Ici, c'était facile, on était déjà en coordonnées WGS84, mais si vous êtes en lambert, un petit ST_TRANSFORM, cf dessous. Et on peut effectivement l'utiliser en triggers (déclencheurs).

J'avais promis le fichier SQL : creation_geocodage_BAN_postgis.sql

Grand test

Essayons cela à grande échelle, j'ai généré une table des centroïdes des communes de Métropole, soit 34 836 lignes, je vais tester pour trouver leur adresse.

UPDATE centrecommune SET "adresse" = geocode_reverse(ST_X(ST_Transform (geom, 4326))::text,ST_Y(ST_Transform (geom, 4326))::text);

Et voilà, test fini, j'ai 8206 erreurs, et cela a pris 1h37 sur mon vieux pc (avec la fibre quand même), il n'y a rien de significatif, normal qu'il n'y ait pas d'adresse au milieu des forêts et des champs, c'était juste un test.

sortie adresse

Alors non, je n'ai pas testé grandeur nature pour le géocodage, je n'ai pas de fichier d'importance avec adresse sous la main, un jour au taff peut-être.

Si vous avez des idées d'optimisation ou de meilleure façon de faire, merci de laisser un commentaire. 


6 commentaires


Écrire un commentaire

Quelle est le deuxième caractère du mot 0mnxj ?