16.4. Interrogazioni spaziali

Le interrogazioni spaziali non sono diverse dalle altre interrogazioni del database. Puoi usare la colonna della geometria come qualsiasi altra colonna del database. Con l’installazione di PostGIS nel nostro database, abbiamo funzioni aggiuntive per interrogare il nostro database.

L’obiettivo di questa lezione: Vedere come le funzioni spaziali vengono applicate in modo simile alle funzioni «normali» non spaziali.

16.4.1. Operatori spaziali

Se vuoi sapere quali punti sono entro una distanza di 2 gradi da un punto (X, Y) puoi farlo con:

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

Risultato:

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

Nota

Il valore di the_geom sopra è stato troncato per ragioni di spazio in questa pagina. Se vuoi vedere il punto in coordinate leggibili dall’uomo, prova qualcosa di simile a quello che hai fatto nella sezione precedente «Visualizza un punto come WKT».

Come sappiamo che l’interrogazione fatta restituisce tutti i punti entro 2 gradi? Perché non 2 metri? O qualsiasi altra unità?

Soluzione

Le unità utilizzate dalla query di esempio sono i gradi, perché il SR utilizzato dal layer è il WGS 84. Si tratta di un SR geografico, il che significa che le sue unità sono in gradi. Un SR proiettato, come le proiezioni UTM, è in metri.

Ricorda che quando scrivi una query, devi conoscere le unità in cui si trova il SR del layer. In questo modo potrai scrivere una query che restituisca i risultati attesi.

16.4.2. Indici spaziali

Puoi anche definire indici spaziali. Un indice spaziale rende le interrogazioni spaziali molto più veloci. Per creare un indice spaziale sulla colonna della geometria, usa:

CREATE INDEX people_geo_idx
  ON people
  USING gist
  (the_geom);

\d people

Risultato:

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

Modifica la tabella cities in modo che la sua colonna della geometria sia spazialmente indicizzata.

Soluzione

CREATE INDEX cities_geo_idx
  ON cities
  USING gist (the_geom);

16.4.4. Dimostrazione delle funzioni spaziali postGIS

Per dimostrare le funzioni spaziali di PostGIS, creerai un nuovo database contenente alcuni dati (fittizi).

Per iniziare, crea un nuovo database (esci prima dalla console di psql):

createdb postgis_demo

Ricordati di installare le estensioni PostGIS:

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

Quindi, importa i dati forniti nella cartella exercise_data/postgis/. Fai riferimento alla lezione precedente per le istruzioni, ma ricorda che dovrai creare una nuova connessione PostGIS al nuovo database. È possibile importare dal terminale o tramite DB Manager. Importa i file nelle seguenti tabelle del database:

• points.shp in building

• lines.shp in road

• polygons.shp in region

LivelliLivelliCarica questi tre vettori del database in QGIS tramite la finestra di dialogo Aggiungi layer PostGIS , come al solito. Quando apri le tabelle degli attributi, noterai che hanno sia un campo id che un campo gid creati dall’importazione PostGIS.

Ora che le tabelle vengono importate, puoi usare PostGIS per interrogare i dati. Torna al tuo terminale (riga di comando) e inserisci il terminale di psql eseguendo:

psql postgis_demo

Vedrai alcune di queste dichiarazioni creando delle viste, in modo da poterle aprire in QGIS e vedere i risultati.

Seleziona per posizione

Ottieni tutti gli edifici nella regione di 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';

Risultato:

 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)

Oppure, se crei una vista:

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';

Aggiungi la vista come un vettore e visualizzalo in QGIS:

Seleziona vicino

Mostra un elenco di tutti i nomi delle regioni adiacenti alla regione di Hokkaido:

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

Risultato:

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

Come vista:

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';

In QGIS:

Nota la regione mancante (Queensland). Ciò potrebbe essere dovuto a un errore di topologia. Fatti come questo possono avvisarti di potenziali problemi nei dati. Per risolverlo senza rimanere intrappolati nelle anomalie che i dati potrebbero avere, potremmo invece utilizzare un buffer intersecato:

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

Questo crea un buffer di 100 metri attorno alla regione di Hokkaido.

L’area più scura è il buffer:

Seleziona usando il buffer:

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';

In questa interrogazione, la vista del buffer originale viene utilizzata come qualsiasi altra tabella. Viene assegnato l’alias a, e il suo campo di geometria,:kbd:a.the_geom, viene utilizzato per selezionare qualsiasi poligono nella tabella region (alias b) che lo interseca. Tuttavia, Hokkaido stesso è escluso da questa affermazione di selezione, perché non lo vogliamo; vogliamo solo le regioni adiacenti.

In QGIS:

Puoi anche selezionare tutti gli oggetti all’interno di una data distanza, senza il passo aggiuntivo di creare un buffer:

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';

Raggiungendo lo stesso risultato, senza necessità del passaggio del buffer temporaneo:

Seleziona valori unici

Mostra un elenco di nomi di città unici per tutti gli edifici nella regione del Queensland:

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

Risultato:

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

Ulteriori esempi …

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

Hai visto come interrogare oggetti spaziali usando le nuove funzioni del database di PostGIS.

16.4.6. What’s Next?

Ora esaminerai le strutture di geometrie più complesse e come crearle utilizzando PostGIS.