Rastrová data

Podporu pro práci s rastrovými daty poskytuje PostGIS až od verze 2.0 a to v rámci rozšíření PostGIS Raster.

Poznámka

PostGIS Raster je součástí rozšíření postgis, není potřeba ho zvlášt přidávat.

Rastrová data mohou být spravována ve dvou formách jako:

  1. IN-DB, tj. uložena přímo v databázi anebo
  2. OUT-OF-DB, tj. mimo databázi.

Nahraní rastrových dat do databáze

Podobně jako u vektorových dat, viz kapitola Dávkové nahrání dat, lze data nahrát pomocí specializovaného nástroje, který je součástí PostGISu a to raster2pgsql.

Důležité

Knihovna GDAL formát PostGIS Raster podporuje nicméně v současnosti neumožnuje konverzi dat z jiných formátu to PostGIS Raster.

Nejprve si stáhneme testovací data DMT (digitální model terénu), soubor dekomprimujeme unzip dmt.zip.

raster2pgsql

Příkaz raster2pgsql funguje obdobně jako shp2pgsql pro vektorová data ve formátu Esri Shapefile. Pomocí tohoto nástroje lze rastrová data do PostGISu naimportovat v obou režimech IN-DB i OUT-OF-DB.

Režim IN-DB

Načtení rastrových dat pomocí raster2pgsql (IN-DB) z příkazové řádky

raster2pgsql -s 3035 dmt.tif ukol_1.dmt | psql pokusnik 2>err

Poznámka

Tento přikaz nicméně může skončit chybou typu

ERROR:  invalid memory alloc request size 1073741824

která naznačuje, že máme nedostatek paměti pro import tohoto rastru.

Tento problém můžeme např. obejít parametrem -Y, který pošle do PGDump příkaz COPY namísto INSERT.

Přidáme ještě užitečný přepínač -C, který nastaví omezení na importovana data. Jinak by například byl ignorován souřadnicový systém a pod.

raster2pgsql -s 3035 -Y -C dmt.tif ukol_1.dmt | psql pokusnik 2>err

Režim OUT-OF-DB

V tomto režimu nedochází k importu vstupních dat do PostGISu, ale pouze k vytvoření odkazu z databáze na původní data. Rastrová data jsou tedy fyzicky uložena mimo databázi.

Načtení rastrových dat pomocí raster2pgsql (OUT-OF-DB) z příkazové řádky

raster2pgsql -s 3035 -R -C `pwd`/dmt.tif ukol_1.dmt_link | psql pokusnik 2>err

Cesta k soubor musí být uplná, jinak nebude link korektní. Pomohli jsem si unixovým příkazem pwd, který vrátí cestu k aktuálnímu adresáři, ve kterém jsou umístěna importovaná data.

Základní metadata

V sekci raster2pgsql jsme naimportovali rastr DMT ve dvou formách jako IN-DB (tabulka ukol_1.dmt) a OUT-OF-DB (tabulka ukol_1.dmt_link).

SELECT r_table_schema,r_table_name,srid,out_db FROM raster_columns;

 r_table_schema | r_table_name | srid | out_db
----------------+--------------+------+--------
 ukol_1         | dmt          | 3035 | {f}
 ukol_1         | dmt_link     | 3035 | {t}

Tabulka raster_columns ukrývá další užitečné informace.

SELECT scale_x,scale_y,blocksize_x,blocksize_y,same_alignment,
 regular_blocking,num_bands,pixel_types,nodata_values,ST_AsText(extent) as extent
 FROM raster_columns where r_table_name = 'dmt';

scale_x          | 25
scale_y          | -25
blocksize_x      | 19615
blocksize_y      | 11119
same_alignment   | t
regular_blocking | f
num_bands        | 1
pixel_types      | {16BUI}
nodata_values    | {65535}
extent           | POLYGON((4470075 3113850,4960450 3113850,4960450 2835875,4470075 2835875,4470075 3113850))

Poznámka

Záporná hodnota scale_y naznačuje orientaci rastru ze severu na jih.

Kde je:

scale_x prostorové rozlišení ve směru osy x
scale_y prostorové rozlišení ve směru osy y
blocksize_x velikost dlaždice ve směru osy x
blocksize_y velikost dlaždice ve směru osy y
same_alignment všechny dlaždice mají stejné zarovnání
regular_blocking všechny dlaždice mají stejný rozměr a nepřekrývají se
num_bands počet kanálů
pixel_types datový typ buněk kanálů
nodata_values hodnota pro no-data jednotlivých kanálů
extent minimální ohraničující obdélník datové vrstvy

Poznámka

Porovnáme-li velikost dlaždice (blocksize_x a blocksize_y) a velikost vstupního rastru (například pomocí nástroje knihovny GDAL gdalinfo, tak dojdeme, že se rastr naimportoval jako jedna dlaždice.

gdalinfo dmt.tif -noct

Pro rozdělení rastrových dat při importu do více dlaždic slouží parametr -t (<šířka>x<výška>) programu raster2pgsql.

Rozdělení dat do více dlaždic při importu z příkazové řádky

Velikost dlaždice zvolíme 400x400px.

raster2pgsql -s 3035 -Y -C -t 400x400 dmt.tif ukol_1.dmt_tiled | psql pokusnik 2>err

Rastr se v tomto případě naimportuje jako 1400 dlaždic.

SELECT COUNT(*) FROM ukol_1.dmt_tiled;

Příklad

Vejce vesmírných oblud v nadmořské výšce nad 300 metrů jsou oslabena. Využijte toho a zlikvidujte je.

Zadání

Určete nadmořskou výšku bodů s výskytem vajec na základě rastru DMT. Vyberte body s nadmořskou výškou větší než 300 metrů.

Řešení

Geometrie tabulky vesmirne_zrudnice je v systému S-JTSK (EPSG:5514), rastrová data v ETRS-89 (EPSG:3035). V rámci řešení tedy musíme počítat s transformaci dat do společného souřadnicového systému pomocí funkce ST_Transform.

-- nastavevíme cestu
SET search_path TO ukol_1, public;

SELECT v.id,ST_Value(r.rast,v.geom) FROM dmt AS r CROSS JOIN
 (SELECT id,ST_Transform(geom_p, 3035) AS geom FROM vesmirne_zrudice) AS v;

-- optimalizovaná verze dotazu (dmt -> dmt_tiled)
SELECT v.id,ST_Value(r.rast,v.geom) FROM dmt_tiled AS r JOIN
 (SELECT id,ST_Transform(geom_p, 3035) AS geom FROM vesmirne_zrudice) AS v ON
 ST_Intersects(r.rast,v.geom);

Výsledek uložíme do nového sloupečku v tabulce vesmirne_zrudnice a vybereme body s nadmořskou výškou větší než 300 metrů.

ALTER TABLE vesmirne_zrudice ADD COLUMN vyska FLOAT;

UPDATE vesmirne_zrudice SET vyska = value FROM
(
 SELECT v.id AS vid,ST_Value(r.rast,v.geom) AS value FROM dmt_tiled AS r JOIN
  (SELECT id,ST_Transform(geom_p, 3035) AS geom FROM vesmirne_zrudice) AS v ON
  ST_Intersects(r.rast,v.geom)
) AS v WHERE id = vid;

SELECT id FROM vesmirne_zrudice WHERE vyska > 300;