16.4. Lesson: Requêtes Spatiales

Les requêtes spatiales ne sont pas différentes des autres requêtes de base de données. Vous pouvez utiliser la colonne de géométrie comme n’importe quelle colonne de base de données. Avec l’installation de PostGIS dans notre base de données, vous avons ajouté des fonctions additionnelles pour requêter la base.

Objectif de cette leçon: Etudier comment les fonctions spatiales sont implémentées de manière similaire aux autres fonctions non spatiales.

16.4.1. Opérateurs spatiaux

Lorsque vous souhaitez connaître quels sont les points distants de 2 degrés par rapport à un point (X,Y), vous pouvez le faire ainsi:

select *
from people
where st_distance(the_geom,'SRID=4326;POINT(33 -34)') < 2;

Résultat:

 id |     name     | house_no | street_id |   phone_no    |   the_geom
----+--------------+----------+-----------+---------------+---------------
  6 | Fault Towers |       34 |         3 | 072 812 31 28 | 01010008040C0
(1 row)

Note

La valeur the_geom a été tronquée pour des questions d’affichage sur cette page. Si vous voulez voir le point exprimé dans des coordonnées lisibles par l’être humain, utilisez quelquechose de similaire à ce que vous avez vu dans la section précédente « Voir un point en WKT ».

Comment savons-nous que la requête précédente retourne tous les points situés à moins de 2 degrés ? Pourquoi pas 2 mètres ? Ou dans une autre unité ?

Réponse

L’unité utilisée dans la requête en exemple est le degré, car le SCR de la couche est WGS 84. Il s’agit d’un SCR géographique, et son unité est donc le degré. Un SCR projeté, tel que les projections UTM, est en mètres.

Rappelez-vous que lorsque vous rédigez une requête, vous devez connaître l’unité du SCR de la couche. Ceci vous permettra d’écrire une requête correcte renvoyant les résultats tels qu’escomptés.

16.4.2. Index spatiaux

Nous pouvons également définir des index spatiaux. Un index spatial permet de rendre les requêtes spatiales plus rapides. Pour créer un index spatial sur une colonne de géométrie, utilisez:

CREATE INDEX people_geo_idx
  ON people
  USING gist
  (the_geom);

\d people

Résultat:

Table "public.people"
   Column   |         Type          |                Modifiers
 -----------+-----------------------+----------------------------------------
  id        | integer               | not null default
            |                       | nextval('people_id_seq'::regclass)
  name      | character varying(50) |
  house_no  | integer               | not null
  street_id | integer               | not null
  phone_no  | character varying     |
  the_geom  | geometry              |
Indexes:
  "people_pkey" PRIMARY KEY, btree (id)
  "people_geo_idx" gist (the_geom)  <-- new spatial key added
  "people_name_idx" btree (name)
Check constraints:
  "people_geom_point_chk" CHECK (st_geometrytype(the_geom) = 'ST_Point'::text
  OR the_geom IS NULL)
Foreign-key constraints:
  "people_street_id_fkey" FOREIGN KEY (street_id) REFERENCES streets(id)

16.4.3. Try Yourself

Modifiez la table des villes de manière à ce que la colonne de géométrie utilise un index spatial.

Réponse

CREATE INDEX cities_geo_idx
  ON cities
  USING gist (the_geom);

16.4.4. Démonstration des fonctions spatiales de PostGIS

Pour faire une démonstration des fonctions spatiales de PostGIS, nous allons créer une base de données contenant des données (fictives).

Pour commencer, créez une nouvelle base de données (sortez d’abord du shell psql):

createdb postgis_demo

Rappelez-vous d’installer les extensions PostGIS:

psql -d postgis_demo -c "CREATE EXTENSION postgis;"

Ensuite, importez les données fournies dans le répertoire exercise_data/postgis/. Référez-vous à la leçon précédente mais souvenez-vous que vous devez créer une nouvelle connexion PostGIS à la nouvelle base de données. Vous pouvez importer depuis le terminal ou via DB Manager. Importez les fichiers dans les tables suivantes de la base:

• points.shp dans building

• lines.shp dans road

• polygons.shp dans region

Chargez ces trois couches dans QGIS via la boîte de dialogue Ajouter des couches PostGIS, comme d’habitude. Lorsque vous ouvrez les tables d’attributs, vous pouvez remarquer qu’elles ont un champ id et un autre gid, créés lors de l’import PostGIS.

Maintenant que les tables sont importées, nous pouvons utiliser PostGIS pour requêter les données. Retournez dans votre terminal (ligne de commande) et entrez le prompt psql en lançant:

psql postgis_demo

Nous allons présenter quelques-uns des traitements SQL en créant des vues avec eux pour que vous puissiez les ouvrir dans QGIS pour visualiser les résultats.

Sélection par localisation

Récupérer tous les bâtiments de la région KwaZulu:

SELECT a.id, a.name, st_astext(a.the_geom) as point
  FROM building a, region b
    WHERE st_within(a.the_geom, b.the_geom)
    AND b.name = 'KwaZulu';

Résultat:

 id | name |                  point
----+------+------------------------------------------
 30 | York | POINT(1622345.23785063 6940490.65844485)
 33 | York | POINT(1622495.65620524 6940403.87862489)
 35 | York | POINT(1622403.09106394 6940212.96302097)
 36 | York | POINT(1622287.38463732 6940357.59605424)
 40 | York | POINT(1621888.19746548 6940508.01440885)
(5 rows)

Ou, si nous créons une vue depuis la requête précédente:

CREATE VIEW vw_select_location AS
  SELECT a.gid, a.name, a.the_geom
    FROM building a, region b
      WHERE st_within(a.the_geom, b.the_geom)
      AND b.name = 'KwaZulu';

Ajoutez la vue comme couche dans QGIS:

Sélection des voisins

Affichez la liste des tous les noms de régions qui touchent la région Hokkaido:

SELECT b.name
  FROM region a, region b
    WHERE st_touches(a.the_geom, b.the_geom)
    AND a.name = 'Hokkaido';

Résultat:

    name
--------------
 Missouri
 Saskatchewan
 Wales
(3 rows)

En tant que vue:

CREATE VIEW vw_regions_adjoining_hokkaido AS
  SELECT b.gid, b.name, b.the_geom
    FROM region a, region b
      WHERE st_touches(a.the_geom, b.the_geom)
      AND a.name = 'Hokkaido';

Dans QGIS:

Vous pouvez noter qu’il manque une région (le Queensland). C’est sans doute dû à une erreur topologique. Des artefacts comme celui-ci peuvent nous alerter sur des problèmes de données. Pour résoudre cette interrogation sans analyser toutes les anomalies des données, nous pouvons utiliser une intersection de tampon à la place:

CREATE VIEW vw_hokkaido_buffer AS
  SELECT gid, ST_BUFFER(the_geom, 100) as the_geom
    FROM region
      WHERE name = 'Hokkaido';

Cette requête créé un tampon de 100 mètres autour de la région d’Hokkaido.

La surface la plus sombre est le tampon:

Sélection à l’aide du tampon:

CREATE VIEW vw_hokkaido_buffer_select AS
  SELECT b.gid, b.name, b.the_geom
    FROM
    (
      SELECT * FROM
        vw_hokkaido_buffer
    ) a,
    region b
    WHERE ST_INTERSECTS(a.the_geom, b.the_geom)
    AND b.name != 'Hokkaido';

Dans cette requête, la vue originelle du tampon est utilisée comme n’importe quelle autre table. Elle dispose d’un alias a et son champ de géométrie, a.the_geom, est utilisé pour sélectionner un polygone dans la table region (alias b) qui l’intersecte. Néanmoins, Hokkaido est exclue de cette sélection car nous n’en voulons pas, nous voulons seulement les régions qui la touche.

Dans QGIS:

Il est également possible de sélectionner les objets à une certaine distance sans avoir besoin de créer un tampon:

CREATE VIEW vw_hokkaido_distance_select AS
  SELECT b.gid, b.name, b.the_geom
    FROM region a, region b
      WHERE ST_DISTANCE (a.the_geom, b.the_geom) < 100
      AND a.name = 'Hokkaido'
      AND b.name != 'Hokkaido';

Cela permet le même résultat sans avoir besoin de l’étape du tampon:

Sélection de valeurs uniques

Affiche la liste des noms de ville uniques pour tous les bâtiments de la région du Queensland:

SELECT DISTINCT a.name
  FROM building a, region b
    WHERE st_within(a.the_geom, b.the_geom)
    AND b.name = 'Queensland';

Résultat:

  name
---------
 Beijing
 Berlin
 Atlanta
(3 rows)

Autres exemples…

CREATE VIEW vw_shortestline AS
  SELECT b.gid AS gid,
        ST_ASTEXT(ST_SHORTESTLINE(a.the_geom, b.the_geom)) as text,
        ST_SHORTESTLINE(a.the_geom, b.the_geom) AS the_geom
    FROM road a, building b
      WHERE a.id=5 AND b.id=22;

CREATE VIEW vw_longestline AS
  SELECT b.gid AS gid,
         ST_ASTEXT(ST_LONGESTLINE(a.the_geom, b.the_geom)) as text,
         ST_LONGESTLINE(a.the_geom, b.the_geom) AS the_geom
    FROM road a, building b
      WHERE a.id=5 AND b.id=22;

CREATE VIEW vw_road_centroid AS
  SELECT a.gid as gid, ST_CENTROID(a.the_geom) as the_geom
    FROM road a
      WHERE a.id = 1;

CREATE VIEW vw_region_centroid AS
  SELECT a.gid as gid, ST_CENTROID(a.the_geom) as the_geom
    FROM region a
      WHERE a.name = 'Saskatchewan';

SELECT ST_PERIMETER(a.the_geom)
  FROM region a
    WHERE a.name='Queensland';

SELECT ST_AREA(a.the_geom)
  FROM region a
    WHERE a.name='Queensland';

CREATE VIEW vw_simplify AS
  SELECT gid, ST_Simplify(the_geom, 20) AS the_geom
    FROM road;

CREATE VIEW vw_simplify_more AS
  SELECT gid, ST_Simplify(the_geom, 50) AS the_geom
    FROM road;

CREATE VIEW vw_convex_hull AS
  SELECT
    ROW_NUMBER() over (order by a.name) as id,
    a.name as town,
    ST_CONVEXHULL(ST_COLLECT(a.the_geom)) AS the_geom
    FROM building a
    GROUP BY a.name;

16.4.5. In Conclusion

Vous avez observer comment faire des requêtes spatiales sur des objets en utilisant les nouvelles fonctions apportées par PostGIS.

16.4.6. What’s Next?

Dans le prochain chapitre, nous allons étudier des structures de géométries plus complexes et voir comment les créer en utilisant PostGIS.