Opérateurs, fonctions et requêtes spatiales

Introduction

Le support des données matricielles (raster) est une « nouveauté » de PostGIS 2.

A chaque pixel de la grille est associée une valeur (altitude, température, densité de population, type d’occupation du sol, perméabilité…). Cette valeur est stockée dans une bande.
PostGIS permet d’associer plusieurs bandes à un raster. On peut par exemple disposer d’un raster de températures mensuelles moyennes qui contient une bande par mois de l’année.

PostGIS stocke chaque raster dans une table de la base de données. La vue raster_columns est l’équivalent pour les raster de la vue geometry_columns pour les données vecteur.
Il offre la possibilité de tuiler le raster, c’est à dire de le fractionner en unités plus petites qui permettront d’accélérer les traitement en utilisant des index.
Chacune de ces tuiles un un carré de x pixels de coté. Chaque tuile (tile) orrespond à une ligne de la table. C’est à dire que la table qui contient le raster compte autant de lignes que le raster compte de tuiles.

De nombreuses opérations sont possibles avec les raster. Référez vous à l’aide de PostGIS1 pour avoir une vision complète des possibilités offertes, ainsi qu’à la présentation initiale du format par ses auteurs :

http://trac.osgeo.org/postgis/wiki/WKTRaster/Documentation01

Objectifs et organisation du TP

L’objectif de ce TP est de découvrir le type de données « raster » supporté par PostGIS depuis la version 2.0. Comme pour les données vectorielles utilisées au cours de deux précédents TP, ce TP  va vous initier aux opérateurs, fonctions et requêtes spécifiques à ce type de données. Chacune des questions posées doit être traitée de la façon suivante :

• Si cela n’est pas donné, recherche du ou des opérateurs et fonctions permettant de répondre à la question dans le cours ou dans la documentation de PostGIS.
• Construction de la requête et répondre à la question

Pour une liste des opérateurs et fonctions spatiales ainsi que des exemples d’utilisation, vous pouvez vous référer à cette partie de la documentation PostGIS :

http://postgis.net/docs/reference.html

Une « cheat sheet » relative aux fonctions raster de PostGIS est disponible ici :

http://www.postgis.us/downloads/postgis21_cheatsheet.pdf

Pré-requis

Pour ce TP vous utiliserez la même base de données que pour les TP précédents. Nous utiliserons le client graphique PgAdmin pour l’affichage des tables et des résultats de requêtes.

Références et supports en ligne

Sites traitant de l’utilisation des rasters dans PostGIS

http://postgis.net/docs/manual-2.1/using_raster_dataman.html

Ajout des données Raster à la base de données

Données utilisées

Pour ce TP, nous allons utiliser les données raster du projet WorldClim (http://worldclim.org), disponibles à l’échelle de la planète sur internet :
l’altitude les variables bioclimatiques du projet worldclim :

• BIO5 = Max Temperature of Warmest Month
• BIO6 = Min Temperature of Coldest Month

Ces fichiers raster sont lourds (1.4 Go par variable), nous avons donc réalisé une extraction des données du territoire français (métropole et Corse). Vous trouverez ces données sur l’ENT.

En complément, nous allons ajouter à notre base de données la couche vecteur de l’ensemble des communes de France issues de la base de données GEOFLA de l’IGN :

http://professionnels.ign.fr/geofla

Utilisez Quantum Gis (ajout d’une couche raster) pour visualiser l’ensemble de ces couches. Les raster vont apparaître uniformément gris. Pour avoir une idée de la variabilité de la valeur des différents pixels, double cliquez sur chacune des couches. Dans l’onglet Style, choisissez « Bande grise unique » et modifiez la palette de couleur en choisissant «Pseudo-Couleurs »

Ajout de la composante communale de GEOFLA à la table commune

Rendez vous dans le répertoire qui contient votre couche shp de la base de données geofla de l’iGN avec la commande « cd »

Shell> cd /home/…/geofla

Dans ce répertoire utilisez la commande shp2pgsql, vue au cours des séances précédentes pour ajouter ce shp à la base de données. La table de destination s’appellera « communes », la colonne géométrique sera nommée « geometrie », la projection utilisée sera le Lambert 93 (code EPSG 2154) et elle sera indexée.

Attention, l’encodage de cette couche est le LATIN1. La commande shp2pgsql dispose de l’option « W » qui vous permet de spécifier cet encodage.

Shell> shp2pgsql -s 2154 -g geometrie -I -W LATIN1 COMMUNE.SHP communes_france | psql -h gis -U login nom_base_de_données

Cette couche de données est en Lambert 93. Pour faciliter les calculs futurs, qui seront effectués sur des raster en WGS 84, nous allons ajouter une colonne géométrique en WGS 84 (nommée the_geom) à notre table des communes, la remplir et l’indexer. Vous pouvez utiliser PgAdmin3 ou psql pour écrire et exécuter ces requêtes.

psql> ALTER TABLE communes_france  ADD column the_geom geometry(‘MULTIPOLYGON’,4326);

psql> UPDATE communes_france  SET the_geom=st_transform(geometrie, 4326);
Mise à jour de la colonne « the_geom » de la table commune avec la géométrie contenue dan sla colonne « geometrie » reprojetée en WGS 84.

psql> CREATE INDEX communes_geom_gist ON communes_france USING gist (the_geom);
Création de l’index sur cette nouvelle colonne

Ajout des rasters à la base de données

Comme pour les données vecteur, PostGIS propose un chargeur de données en ligne de commande qui s’appelle raster2pgsql.

Consultez la documentation en ligne de PostGIS ou utilisez la commande suivante pour savoir quelles sont les options disponibles :

Shell> raster2pgsql – ?

Le système de coordonnées des grilles que nous allons intégrer à la base de données est le WGS 84 (code EPSG : 4326), nous créerons des tuiles de 5 pixels de coté et nous indexerons la table.

Nous intégrerons le raster bio_5_fr.tif dans une table nommée temperature_max, le raster bio_6_fr.tif dans une table nommée temperature_min et le raster d’élévation nommé alt_fr.tif dans une table nommée altitude.

Pour cela, rendez vous dans le répertoire contenant ces données avec la commande cd :

Shell> cd /chemin/vers/mes/donnees/raster
Shell> /chemin/vers/raster2pgsql -s 4326 -I -t 5×5 bio_5_fr.tif temperature_max | psql -h gis -U login nom_base_de_données
Shell> /chemin/vers/raster2pgsql -s 4326 -I -t 5×5 bio_6_fr.tif temperature_min | psql -h gis -U login nom_base_de_données
Shell> /chemin/vers/raster2pgsql -s 4326 -I -t 5×5 alt_fr.tif altitude | psql -h gis -U login nom_base_de_données

Requêtes géographiques avec PostGIS

Requêtes d’information sur les objets géographiques

Quelle est la taille des pixels des raster « temperature_min » et « temperature_max’ ? En quelle unité est-elle exprimée ?

Utilisez les fonctions ST_PixelHeight() et ST_PixelWidth()

De combien de tuiles sont composés chacun des raster ?

Quelle est la température minimale du centroïde de la commune de Montpellier ?

Vous utiliserez les fonctions ST_intersection(rast,the_geom) et ST_Centroïd(geometry)
Est-ce cohérent ? -> http://worldclim.org/faq

Quelle est la température minimale moyenne de la commune de Montpellier ?

Vous utiliserez la fonction ST_intersection(rast, the_geom) et l’opérateur d’agrégation AVG()

Quelle est la commune la plus chaude du département de l’Aude ?

Quelle est la commune la moins élevée du département des Pyrénées Orientales ?

Création et utilisation d’un raster multibande

Nous allons, pour simplifier les calculs utiliser la possibilité d’ajouter des bandes aux pixels du raster.
Commençons par ajouter au raster altitude une bande contenant la température minimale du mois le plus froid :

UPDATE altitude SET rast=ST_AddBand(altitude.rast, min.rast)
FROM temperature_min min
WHERE ST_UpperLeftX(min.rast) = ST_UpperLeftX(altitude.rast)
AND ST_UpperLeftY(min.rast) = ST_UpperLeftY(altitude.rast);

Ajoutons ensuite une bande contenant la température la plus élevée du mois le plus chaud :

UPDATE altitude SET rast=ST_AddBand(altitude.rast, max.rast)
FROM temperature_max max
WHERE ST_UpperLeftX(max.rast) = ST_UpperLeftX(altitude.rast)
AND ST_UpperLeftY(max.rast) = ST_UpperLeftY(altitude.rast);

Actualisons les informations relatives à nos tables contenues dans la tables raster_columns :

SELECT AddRasterConstraints(‘altitude’::name, ‘rast’::name);
SELECT AddRasterConstraints(‘temperature_min’::name, ‘rast’::name);
SELECT AddRasterConstraints(‘temperature_max’::name, ‘rast’::name);

Vérifions le nombre de bande de notre raster d’altitude :

SELECT num_bands FROM raster_columns WHERE r_table_name = ‘altitude’;

Requêtes complexes

Affichez l’altitude moyenne, la moyenne de la température la plus basse du mois le plus froid, la moyenne de la température la plus élevée du mois le plus chaud de la commune de Sanilhac-Sagries

Ajouter aux communes l’information relative à l’altitude de leur centroïde

•     ajouter un attribut « altitude » de type numeric à la table « communes_france »
•     mettre à jour la valeur de cette attribut en interrogeant le raster (utilisez la fonction st_value(rast, geometry), le code de la région est 91…
•     visualisez la couche dans Quantum Gis en utilisez une symbologie mettant en évidence la variation d’altitude

Bonus pour aller plus loin

Extraire du MNT les zones basses (altitude < 50 cm) des communes littorales de la région de Montpellier (PALAVAS,  LATTES, MAUGUIO, PEROLS, LA GRANDE-MOTTE, VILLENEUVE LES MAGUELONE, VIC LA GARDIOLE, FRONTIGNAN)

On extrait tout d’abord la géométrie et la valeur de chacun des pixels de la table altitude qui sont contenus dans les communes littorales

SELECT communes_france.nom_comm, (st_intersection(altitude.rast, communes.the_geom)).val AS val,
(st_intersection(altitude.rast, communes.the_geom)).geom AS geom
FROM communes_france, altitude
WHERE communes_france.nom_comm::text ILIKE ANY(array[‘PALAVAS%’,’LATTES%’,’MAUGUIO%’,’PEROLS%’,’%GRANDE-MOTTE%’,’%MAGUELONE’,’%GARDIOLE’,’FRONTIGNAN’]) AND st_intersects(altitude.rast, communes_france.the_geom)

De cet ensemble de géométries et de valeurs, nous allons extraire et regrouper les géométrie des pixels de moins de 50 m d’altitude

WITH valeurs_tous_pixels AS (
SELECT communes_france.nom_comm, (st_intersection(altitude.rast, communes_france.the_geom)).val  AS val, (st_intersection(altitude.rast, communes_france.the_geom)).geom AS geom
            FROM communes_france, altitude
           WHERE communes_france.nom_comm::text ILIKE ANY(array[‘PALAVAS%’,’LATTE%’,
‘MAUGUIO%’,’PEROLS%’,’%GRANDE-MOTTE %’,’%MAGUELONE’,’%GARDIOLE’,’FRONTIGNAN’])
AND st_intersects(altitude.rast, communes_france.the_geom)
        )
SELECT valeurs_tous_pixels.nom_comm, avg(valeurs_tous_pixels.val) AS moyenne_temp_min, st_union(valeurs_tous_pixels.geom) AS the_geom
FROM valeurs_tous_pixels
WHERE valeurs_tous_pixels.val < 50::double precision
GROUP BY valeurs_tous_pixels.nom_comm;

On crée la vue.

CREATE OR REPLACE VIEW zone_inondable AS …

Solution proposée

-- 3.1.3

SELECT nom_comm,
--(ST_intersection(rast, st_centroid(the_geom))).val,
ST_value(rast,1,st_centroid(the_geom))
  FROM communes_france, temperature_min
  WHERE nom_comm = 'MONTPELLIER'
  AND st_intersects(rast,st_centroid(the_geom))     -- pour rendre la requête plus rapide
  ;

-- 3.1.4

SELECT nom_comm, avg(val)/10 AS temp_min_moyenn
FROM (
SELECT nom_comm,
(ST_intersection(rast, the_geom)).val
  FROM communes_france, temperature_min
  WHERE nom_comm = 'MONTPELLIER'
  AND st_intersects(rast,the_geom)
  ) AS valeur_des_differents_raster         -- Les sous-requêtes doivent être "nommées"
GROUP BY nom_comm;

/* Alternative en utilisant les CTE (Comman table expression)
   http://www.postgresql.org/docs/9.0/static/queries-with.html */

WITH valeur_des_differents_raster AS (
SELECT nom_comm,
(ST_intersection(rast, the_geom)).val
  FROM communes_france, temperature_min
  WHERE nom_comm = 'MONTPELLIER'
  AND st_intersects(rast,the_geom)
  )
SELECT nom_comm, avg(val)/10 AS moyenne_temp_min
FROM valeur_des_differents_raster
GROUP BY nom_comm;

-- 3.1.5
SELECT nom_comm, avg(val)/10 AS moyenne_temp_max
FROM (
SELECT nom_comm,
(ST_intersection(rast, the_geom)).val
  FROM communes_france, temperature_max
  WHERE code_dept = '11'
  AND st_intersects(rast,the_geom)
  ) AS valeur_des_differents_raster       
GROUP BY nom_comm
ORDER BY 2 DESC LIMIT 1;

-- 3.1.6
SELECT nom_comm, avg(val) AS moyenne_altitude
FROM (
SELECT nom_comm,
(ST_intersection(rast, the_geom)).val
  FROM communes_france, altitude
  WHERE code_dept = '34'
  AND st_intersects(rast,the_geom)
  ) AS valeur_des_differents_raster       
GROUP BY nom_comm
ORDER BY 2 ASC LIMIT 1;

-- 3.3.1
SELECT nom_comm, avg(altitude) AS moyenne_altitude, avg(temp_min)/10 AS moyenne_temp_min, avg(temp_max)/10 AS moyenne_temp_max
FROM (
SELECT nom_comm,
(ST_intersection(the_geom, rast, 1)).val AS altitude,
(ST_intersection(the_geom, rast, 2)).val AS temp_min,
(ST_intersection(the_geom, rast, 3)).val AS temp_max
  FROM communes_france, altitude
  WHERE lower(nom_comm) = lower('Sanilhac-Sagries')
  AND st_intersects(rast,the_geom)
  ) AS valeur_des_differents_raster       
GROUP BY nom_comm
ORDER BY 2 ASC LIMIT 1;

-- 3.3.2
ALTER TABLE communes ADD COLUMN altitude numeric;

UPDATE communes SET altitude = ST_VALUE(rast, 1, ST_CENTROID(the_geom))
FROM altitude WHERE ST_INTERSECTS(rast, ST_CENTROID(the_geom));

Dans le cadre d’un projet de collecte de données en ligne, nous souhaitions mettre en place depuis longtemps un système d’envoi de mail qui prévient à intervalle régulier les contributeurs, que leurs données ont été examinnée par un « expert » et qu’elles posent question.

L’idée générale était de mettre en place un trigger qui se déclenche aprés chaque modification du statut de validation d’une observation mais cela peut générer beaucoup de messages en cas de « validation » par lot.

Une vue sollicitée à intervalle régulier par une tâche cron sera donc plus appropriée.

pgmail se chargera de l’envoi des mails. C’est une fonction écrite en tcl qui a été créée pour envoyer des mails depuis la base de données.

L’installation de pgmail nécessite l’extension pltclu (u pour untrusted, c’est à dire que les fonctions sytèmes sont accessibles au language.

CREATE EXTENSION pltclu;

On installe ensuite la fonction pgmail dans le schéma qui nous convient le mieux :

CREATE OR REPLACE FUNCTION outils.pgmail(text, text, text, text)
RETURNS integer AS
[…]

Il ne reste plus qu’à ècrire le mail.

On souhaite envoyer un mail et un seul à chaque contributeur pour l’ensemble des données dont le statut de validation à changé.

Nous avons pour cela une table de suivi qui enregistre toutes les versions de chacun des tuples de la table de saisie ainsi que la date de la modification, son auteur et l’opération effectuée (INSERT, UPDATE, DELETE). Cette table et les triggers qui l’alimentent sont directement issus des exemples de la documentation de postgreSQL : http://docs.postgresql.fr/9.2/plpgsql-trigger.html#plpgsql-trigger-audit-example

Nous allons nous intéresser dans notre cas aux colonnes date_operationn validateur, statut_validation et decision_validation

L’idée générale est de detecter les lignes modifiées par un validateur (utilisateur = validateur) à la date du jour et d’en extraire la dernière valeur saisie pour les colonnes validation et decision_validation.

Les fonctions de fenêtre lead et first_value vont nous permettre cela :

http://docs.postgresql.fr/9.2/functions-window.html

La fonction d’aggrégation string_agg, d’aggréger les différentes lignes de résultat :

http://docs.postgresql.fr/9.2/functions-aggregate.html

Il y a probablement plus simple ou efficace comme requête mais celle-ci est une bonne première version fonctionnelle.

A reprendre avec les fonction window !

Voici la présentation réalisée lors de notre Conférence Technique et Thématique « Géomatique » du dernier congrès des CENs qui s’est tenu à Thionville du 1er au 5 octobre 2014.

Regina Obe et les requêtes SQL

Writing better PostGIS queries — REGINA OBE, PARAGON CORPORATION from FOSS4G on Vimeo.

Paul Ramsey

PostGIS Feature Frenzy — Paul Ramsey from FOSS4G on Vimeo.

Vincent Picavet et PostgreSQL 9.4

Gimme some YeSQL! – and a GIS — Vincent Picavet, Oslandia from FOSS4G on Vimeo.

Pour diverses raisons, on se souci de moins en moins de la taille des fichiers que nous produisons (diaporama de plus de 100 Mo, rapports d’un poids similaire).

Cependant, si ce comportement se généralise au sein de l’équipe, la capacité de stockage de notre serveur de fichiers et sa sauvegarde sont rapidement malmenées.

Afin surveiller un peu cela, et pour pouvoir « harceler » un peu les collègues concernés, nous avons mis en place un petit script shell qui utilise les commandes find, ls et awk

Lancé toutes les nuits, il génère un petit fichier de log avec les fichiers répondant au critère :

echo fichiers de + de 5 Mo >>/home/cenlr/fichiers_produits_`date +%Y-%m-%d`.txt find /home/cenlr/ -type f -size +5000k -mtime -1 -exec ls -lh {} \; | awk ‘/^-/{print $5, substr($0, index($0,$9))}’ >> /home/cenlr/fichiers_produits_`date +%Y-%m-%d`.txt

Ressources utilisées :

• http://forum.ubuntu-fr.org/viewtopic.php?pid=11075911#p11075911
• http://www.zem.fr/trouver-les-fichiers-sur-linux-ubuntu/

Le nouveau protocole national et la charte régionale du SINP ont été discutés, partagés et adoptés en concertation au printemps 2013.

Ces deux documents précisent les modalités de la publication et de la circulation des données naturalistes. Ils en précisent aussi le droit d’accès pour chacun des grands types d’acteurs (autorité publique, gestionnaires d’espaces naturels, chercheurs, grand public, bureaux d’études). Ainsi les gestionnaires d’espaces naturels peuvent prétendre à l’accès aux données du SINP les plus précises possibles pour l’accomplissement de ses missions.

Le CEN L-R produit avec le soutien financier de partenaires publics et privés de nombreuses données sur la nature qui nous entoure. Ces informations ont vocation à circuler et à être réutilisées. Par ailleurs, le CEN L-R recueille et gère dans son système d’information des données de bénévoles et de sympathisants qui souhaitent que leurs observations soit utiles à la conservation de la nature et qu’elles ne restent pas en dormance au sein d’une base de données.
Le CEN L- R opérait déjà des exports de données vers ses partenaires mais la cadre proposé par le SINP simplifie et allège la gestion administrative des ces échanges. La charte du SINP, dès lors qu’elle est signée par les partenaires, définit les modalités d’échanges de données entre les acteurs du SINP et remplace les multiples conventions précédemment nécessaires.

A l’automne 2013, le CEN L-R a donc déposé à la DREAL une demande officielle d’adhésion au protocole national et à la charte régionale du SINP.
La démarche qui est aujourd’hui lancée comprend différentes étapes :

• le recensement et la description des grands lots de données
• la mise en forme de ces lots au format national publié fin 2013. Les données seront mise en forme conformément à ce standard, sans floutage géographique
• la mise à disposition de ces lots auprès de la DREAL
• la mise à disposition de ces données auprès des têtes de réseau du SINP signataires de la charte régionale
• la mise en place d’un catalogue de données

Les 3 premières étapes seront accomplies avant le 1er mai 2014. Ce sont prés de 120 000 données collectées sur les espaces naturels de la région qui seront ainsi mises à disposition des acteurs régionaux.
Le calendrier de la mise à disposition des données auprès des têtes de réseau dépendra de l’adhésion de leurs structures animatrices au protocole et à la charte.

A l’avenir, un export régulier de données sera réalisé vers la DREAL et les différentes têtes de réseau. La fréquence et les modalités techniques de ces export seront définies avec chacun des partenaires.

Le CEN L-R informera les naturalistes qui participent à la connaissance de la biodiversité régionale à travers l’utilisation de SICEN, des modalités de diffusion des données ainsi collectées auprès de l’état et des têtes de réseau du SINP.

J’utilise depuis Qgis 1.8 les possibilités de personnalisation du formulaire de renseignement des attributs.

Une nouveauté (une autre!) a fait son apparition avec la sortie de QGis 2, c’est la possibilité d’utiliser des couples clé/valeur stockés dans une table de la base de données.

J’aimerai aller un peu plus loin et appliquer un « filtre d’expression », consistant à ne proposer dans la liste déroulante que les couples clé/valeur des objets présents dans l’emprise courante de ma carte.

Concrètement, j’ai besoin de délimiter du parcellaire et je dispose pour cela du scan du cadastre et des localisants de parcelles fournis par la BD Parcellaire.

Lors du dessin d’une nouvelle parcelle, j’aimerai que les valeurs proposées pour renseigner le numéro de parcelle soient issues des valeurs de numéro de parcelle des localisants contenus dans l’emprise actuelle de ma carte, et non de l’ensemble des localisants de ma bd parcellaire.

L’éditeur de filtre ne propose rien qui soit relatif à mon objet « carte ».

En cherchant un peu partout une solution à mon problème, j’ai bien compris que python serait la solution.

Les 3 ressources ci-dessous m’ont bien aidé :

  – http://nathanw.net/2012/11/10/user-defi … -for-qgis/
  – http://www.3liz.com/blog/rldhont/index. … ns-de-QGIS
  – http://gis.stackexchange.com/questions/ … n-qgis-1-8

J’ai donc créé un script python (userfunctions.py), qui contient une fonction appelée current_canvas_extent, qui retourne une géométrie de type Polygone représentant l’emprise courante de la fenêtre.

from qgis.utils import qgsfunction from qgis.utils import iface from qgis.core import QGis @qgsfunction(0, « Python ») def current_canvas_extent(values, feature, parent): «  » » retourne l etendue courante de la carte «  » » extend = iface.mapCanvas().extent().asWktPolygon() return extend

Cette fonction peut-être modifiée pour ne pas passer par le format WKT mais directement par un objet Geometry de QGis :

from qgis.utils import qgsfunction from qgis.utils import iface from qgis.core import (QGis, QgsGeometry) @qgsfunction(0, « Python ») def current_canvas_extent(values, feature, parent): «  » » retourne l etendue courante de la carte «  » » extend = QgsGeometry.fromRect( iface.mapCanvas().extent()) return extend

J’appelle se script au démarrage de QGis (ajouter à la commande de lancement de QGis l’option « –code chemin\vers\userfunctions.py ») J’ai maintenant accès dans mon éditeur d’expression à cette « variable ».

Et je peux demander à QGis, grâce à l’éditeur d’expression de ne me proposer dans la liste que les parcelles présentes dans l’emprise de ma carte :

Ca n’a jamais été aussi vrai :

Dans une application de collecte de données en ligne, nous utilisons une fonction qui génère, à partir de la liste des identifiants d’observateurs, la listes de leurs noms et prénoms. Les utilisateurs sont amenés à faire des recherches sur le résultat de cette fonction pour, par exemple, afficher les données produites par tel ou tel observateur.

La requête ci-dessous met environ 20 secondes à renvoyer un résultat :

SELECT * FROM saisie.saisie_observation WHERE md.liste_nom_auteur(observateur) ILIKE ‘%BOSS%’

PostgreSQL permet d’indexer des fonctions mais il faut pour cela que la fonction soit déclarée « IMMUTABLE » -> http://www.postgresql.org/docs/9.2/static/sql-createfunction.html

CREATE OR REPLACE FUNCTION md.liste_nom_auteur(text) RETURNS text AS $BODY$ DECLARE var_liste_sql_personne ALIAS for $1; BEGIN RETURN string_agg(nom || ' ' || prenom,' & ') FROM (SELECT regexp_split_to_table(var_liste_sql_personne,'&')::integer as id_personne) t LEFT JOIN md.personne USING(id_personne); END; $BODY$ LANGUAGE plpgsql IMMUTABLE COST 100;

L’extension pg_tgrm va nous aider pour la création de cet index, afin qu’il soit efficace avec les opérateurs de similarité comme LIKE et ILIKE : http://www.postgresql.org/docs/9.2/static/pgtrgm.html

CREATE EXTENSION pg_trgm SCHEMA public VERSION "1.0";

Création de l’index sur md.liste_nom_auteur(observateur) utilisé dans le filtre de la grille

CREATE INDEX saisie_observation_liste_observateurs_idx ON saisie.saisie_observation USING gist(md.liste_nom_auteur(observateur) gist_trgm_ops);

La requête est désormais exécutée en moins de 50 ms !

SELECT * FROM saisie.saisie_observation WHERE md.liste_nom_auteur(observateur) ILIKE '%BOSS%';

L'outil de saisie des données naturaliste est enfin motorisé par les dernières version de PostgreSQL (9.2) et de Postgis (2.0).

Nous allons utiliser les trigger de PostgreSQL et la capacité de PostGIS 2 à gérer les données raster pour renseigner automatiquement l'altitude des données ponctuelles renseignées dans l'interface.

Le MNT de la BD TOPO de l'IGN a été intégré à la base pour l'ensemble de la région dans la table mnt_lr du schéma ign_bd_topo.

La commande utilisée pour créer cette table à partir d'un raster régional est la suivante :

raster2pgsql -s 2154 -C -I -r -M -F -t 100x100 mnt_lr.tif ign_bd_topo.mnt_lr | psql -h localhost -d sicen -U dba

Nous allons donc utiliser un trigger qui se déclenchera pour chaque ligne avant chaque insertion ou mise à jour.

Création de la fonction

Création du déclencheur

CREATE TRIGGER maj_altitude BEFORE INSERT OR UPDATE ON saisie.saisie_observation FOR EACH ROW EXECUTE PROCEDURE saisie.renseigne_altitude();

Connaitre le nombre d’observation par commune

Voici comment arriver à un résultat qu’on recherche souvent avec QGis. Je dispose d’un maillage de mon territoire (maille régulière ou non) et d’une couche d’observation. J’aimerai savoir pour chaque maille, le nombre d’observation qu’elle contient.

L’exemple qui suit est inspiré de cet article : http://datagistips.blogspot.fr/2012/04/le-carroyage-avec-qgis-et-le-plugin.html

Données sources

• les communes de l’Hérault
• Une couche d’observation

Lancement du plugin QMarxan

Paramétrage du plugin

On va renseigner ici

• le nom de la couche que l’ion veut qualifier (notre maillage),
• la couche utilisée pour faire le calcul (ici nos données d’observations),
• le type de statistique que l’on souhaite calculer
• et le nom de la colonne qui receuillera la statistique

Résultat de l’analyse thématique

Les données naturalistes que nous collectons sont consolidées dans le SIG à traveur une interface web de saisie.

Une table de la base de données (saisie.suivi_saisie_observation) recence toutes les opérations d'insertion, de modification ou de suppression de ces donénes.

Nous allons interroger cette table pour calculer la somme mensuelle des données saisie via l'interface ainsi que la somme cumulée depuis sa mise en place en aout 2010.

WITH calendrier AS (
SELECT mois, annee
FROM generate_series(1,12) mois, generate_series(2010,extract(year FROM now())::integer)annee
ORDER BY 2,1),

effort_saisie AS (
SELECT extract(‘month’ FROM date_operation) as mois, extract(‘year’ FROM date_operation) as annee, count(id_obs) as nb_donnees
FROM saisie.suivi_saisie_observation
WHERE operation =’INSERT’
GROUP BY extract(‘month’ FROM date_operation), extract(‘year’ FROM date_operation))

SELECT mois, annee, COALESCE(nb_donnees,0) AS saisie_mensuelle, sum(COALESCE(nb_donnees,0)) OVER (ORDER BY annee,mois) AS saisie_cumulee
FROM calendrier LEFT JOIN effort_saisie USING(mois,annee)
WHERE (’01’||’/’||mois||’/’||annee)::date BETWEEN ‘2010-08-01’::date AND now()::date
ORDER BY 2,1

Comment calculer l’altitude moyenne des communes héraultaises ?

Données utilisées

• l’altitude mondiale fournie par WorldClim : http://biogeo.ucdavis.edu/data/climate/worldclim/1_4/grid/cur/alt_30s_bil.zip
• la couche vecteur des communes de l’IGN (geofla) : http://professionnels.ign.fr/geofla

Outils utilisés

Quantum Gis et GdalTools

→pour le découpage du ratser mondial :

gdal_translate -projwin -5.27 51.50 9.84 41.26 -of GTiff ~/Documents/tp_postgis_raster/alt_30s_bil/alt.bil ~/Documents/tp_postgis_raster/alt_fr

PostGIS 2

→chargement du raster dans postgis

/usr/pgsql-9.2/bin/raster2pgsql -s 4326 -t 5x5 -I ~/Documents/tp_postgis_raster/alt_fr.tif alt_5x5|psql -h localhost -U dba fmin206

→chargement du shape des communes dans postgis

shp2pgsql -d -s 2154 -g the_geom -W LATIN1 -I ~/Documents/tp_postgis_raster/GEOFLA/COMMUNE.SHP communes|psql -h localhost -U dba fmin206

Requête finale

WITH altitudes as (
SELECT insee_com, nom_comm, (ST_intersection(rast, st_transform(the_geom, 4326))).val::integer AS altitude 
FROM alt_5x5 JOIN communes ON ST_Intersects(st_transform(the_geom, 4326), rast) 
WHERE insee_com like '34%') 
SELECT insee_com, nom_comm, avg(altitude) AS altitude_moyenne 
FROM altitudes 
GROUP BY insee_com, nom_comm

Etat initial

PostgreSQL : PostgreSQL 8.4.3 on i686-redhat-linux-gnu, compiled by GCC gcc (GCC) 4.1.2 20080704 (Red Hat 4.1.2-46), 32-bit

PostGIS : POSTGIS= »1.5.0″ GEOS= »3.2.0-CAPI-1.6.0″ PROJ= »Rel. 4.6.1, 21 August 2008″ LIBXML= »2.6.26″ USE_STATS

Etat attendu

PostgreSQL : 9.2.x

PostGIS : 2.x

Documents utilisés :

• http://fedoraproject.org/wiki/EPEL#How_can_I_use_these_extra_packages.3F
• http://trac.osgeo.org/postgis/wiki/UsersWikiPostGIS20CentOS6pgdg
• http://postgis.refractions.net/documentation/manual-2.0/postgis_installation.html#hard_upgrade

Démarche

Sauvegarde des bases de données désirées de l’ancien serveur

on monte le répertoire de sauvagarde

mount -t cifs //192.168.1.240/sauvegardes -o rw,username=admin,password=cen122009,noperm /mnt/svg/

on sauvegarde les bases une à une

pg_dump -h localhost -p 5432 -U dba -Fc -N temp -b -v -f « /mnt/svg/sicen_2013_01_14.backup » sicen

on sauvegarde les rôles

pg_dumpall -g > roles.sql

on monte le répertoire de sauvegarde sur le nouveau serveur

mount -t cifs //192.168.1.240/sauvegardes -o rw,username=admin,password=cen122009,noperm /mnt/svg/

on crée la nouvelle base de données avec PostGIS

on recrée les rôles en utilsant rôles.sql

on lance le script de restauration (hard upgrade)

perl /usr/pgsql-9.2/share/contrib/postgis-2.0/postgis_restore.pl « /home/admin/Documents/sicen_2013_01_14.backup » | psql -h 127.0.0.1 -p 5432 -U dba newdb 2> erreurs.txt

Les erreurs remontées sont dûes à des fonctions qui n’existent plus dans cette version de PostGIS ( cas de area2d(geometry), buffer(geometry, double precision) et centroid(geometry) dont que nous avons recréées à partir des définitions SQL de leurs version modernes (ST_buffer(geometry, double precision)…)
 

Objectif

mettre à disposition de chaque structure productrice de données dans la base de données :

• un export des données qui lui sont attribuées
• à intervalle régulier
• en l'informant par mail de la mise à jour.

#On se rend dans le répertoire où seront mises à disposition les données

cd /chemin/vers/le/repertoire/de/donnees/

#Export des données concernées en shp

pgsql2shp -f export_cenlr.shp -h localhost -u admin odopap "SELECT * FROM export.tous_point_espece_selon_format_esri WHERE id_entite ILIKE('atlasop%') AND structure = 'CENLR'"

#suppression de l'ancienne archive

rm -f export_cenlr.zip

#compression du nouvel export

zip export_cenlr.zip export_cenlr.*

#Suppression des fichiers sources (shp et consorts)

rm -f export_cenlr.*

#Envoi du mail au destinataire des fichiers avec copie à l'expéditeur

echo -e "Bonjour,\n\nL'export des données de l'atlas des observateurs de votre structure a été mis à jour.\n\nL'archive est accessible ici :\nhttp://consultation.libellules-et-papillons-lr.org/chemin/vers/le/repertoire/de/donnees/export_cenlr.zip\nCet export est mis à jour dans ce même dossier tous les lundis.\n\nCordialement,\n\nMathieu Bossaert" | mail -r webmestre@cenlr.org -c webmestre@cenlr.org -s "export des données de l'atlas" cenlr@mon_fai.org

#Configuration de cron avec l'ajout de cette ligne à la crontab

# lancement de la tâche chaque lundi à 1h00

0 1 * * 1 /root/maj_export_structures.sh