Virtuaalirasterit

Virtuaalirasterit ovat hyödyllinen GDAL-konsepti suurten rasteriaineistojen hallintaan silloin, kun aineisto on jaettu päällekkäisyyttä sisältämättömiin karttalehtiin. Virtuaalirasterit eivät sen sijaan sovellu esimerkiksi aikasarjojen tai päällekkäisten rasterien, kuten kaukokartoitusruutujen, hallintaan.

Teknisesti virtuaalirasteri on vain pieni xml-tiedosto, joka kertoo GDALille, missä varsinaiset datatiedostot sijaitsevat, mutta käyttäjän näkökulmasta virtuaalirastereita voidaan käsitellä pitkälti kuten mitä tahansa muuta rasterimuotoa. Virtuaalirasterit voivat sisältää rasteridataa missä tahansa GDALin tukemassa tiedostomuodossa. Virtuaalirasterit ovat hyödyllisiä, koska niiden avulla suuria aineistoja voidaan käsitellä ikään kuin ne olisivat yksi tiedosto, jolloin oikeiden tiedostojen etsimistä ei tarvita.

Esimerkiksi Maanmittauslaitoksen 2 m:n ja 10 m:n korkeusmallit ovat saatavilla Roihussa. Nämä aineistot on jaettu useisiin tif-tiedostoihin (karttalehtiin), ja jos haluaisimme esimerkiksi laskea vyöhyketilastoja alueille, jotka sijaitsevat hajallaan eri puolilla Suomea, meidän pitäisi jotenkin selvittää, mikä tiedosto kattaa minkäkin alueen, ja laskea tilastot oikeasta tiedostosta. Lisähaasteita syntyisi, jos alue, jolle haluamme laskea tilastoja, sattuu sijaitsemaan kahden tai useamman karttalehden rajalla. Samankaltaisia reunaefekteihin liittyviä ongelmia syntyy esimerkiksi käytettäessä fokustoimintoja, joissa tarvitaan myös ympäröivien tiedostojen tietoja. Nämä ongelmat voidaan helposti välttää luomalla koko tutkimusalueelle virtuaalirasteri, jolloin GDAL huolehtii automaattisesti edellä mainituista tilanteista.

Virtuaalirastereita on mahdollista käyttää myös niin, että vain pieni xml-tiedosto tallennetaan paikallisesti ja suuret rasteritiedostot sijaitsevat Altaassa, Amazon S3:ssa, julkisella palvelimella tai missä tahansa muussa GDALin virtuaaliohjainten tukemassa paikassa. Data siirretään paikallisesti vain siltä alueelta ja sillä zoomaustasolla, jota pyydetään, kun virtuaalirasteri avataan. Parhaiten toimiva tiedostomuoto rasteridatan tallentamiseen etäpalveluun on Cloud optimized GeoTIFF, mutta myös muut muodot ovat mahdollisia.

Virtuaalirasterien käyttö

Virtuaalirastereita tukevat kaikki GDAL-pohjaiset työkalut, mukaan lukien Pythonin ja R:n paikkatietopaketit, ArcGIS, FME, GrassGIS, MapInfo, QGIS ja SagaGIS. Lisäksi ArcGISillä ja MapInfolla on kummallakin myös oma GDALin virtuaalirasteria muistuttava virtuaalirasterimuotonsa.

Vain osan suuren virtuaalirasterin lukeminen bboxin avulla

Virtuaalirasterien tehokkaan käytön avain on usein tuki sille, että suuresta virtuaalirasterista voidaan lukea vain osa. Tämä onnistuu sekä R:llä että Pythonilla; monissa muissa työkaluissa rajaus täytyy tehdä erillisenä vaiheena ja rajattu data tallentaa tiedostoon.

GDALilla on helppo rajata pieni osa suuresta virtuaalirasterista:

gdal_translate -projwin 614500 6668500 644500 6640500 test.vrt test_clip.tif

R:

library(terra)  
vrt <- rast("test.vrt")  
data = crop(vrt , ext(614500, 644500, 6640500, 6668500))

Python:

import rasterio
from rasterio.windows import from_bounds
vrt_path = "test.vrt"
with rasterio.open(vrt_path) as src:
    rst = src.read(window=from_bounds(614500, 6640500, 644500, 6668500, src.transform))

Analyysin suorittaminen, joka lukee vain osia virtuaalirasterista

On mahdollista työskennellä erittäin suurten virtuaalirasterien kanssa silloin, kun analyysin ei tarvitse tuottaa samankokoista rasteritulosta. Hyvä esimerkki tästä on vyöhyketilastojen laskeminen polygonialueille, jotka sijaitsevat hajallaan laajalla alueella; katso esimerkki CSC:n training githubista. Tässä tapauksessa rasteridataa luetaan vain niiltä alueilta, jotka menevät päällekkäin polygonien kanssa.

Suurten virtuaalirasterien visuaalinen käyttö

On hyvä huomata, että vaikka koko Suomen kattavan 2 m:n korkeusmallin analysointi perus-.vrt-tiedoston avulla on täysin mahdollista supertietokoneessa, datan tarkastelu esimerkiksi QGISillä ei ole käytännöllistä näin suurelle aineistolle ilman lisäoptimointia. Jos haluat tarkastella suurta virtuaalirasteria sujuvasti, sinun täytyy tehdä muutama asia:

Luo yleiskatselutasot virtuaalirasterillesi gdaladdo-komennolla. Kannattaa huolehtia siitä, etteivät yleiskatselutasot ole niin suuria, että niistä itsestään tulee valtavan suuri tiedosto.
Jos virtuaalirasterisi on todella suuri, on järkevää luoda hierarkkinen virtuaalirasterirakenne, jossa ylimmän tason virtuaalirasteri viittaa pienempiin virtuaalirastereihin, jotka puolestaan viittaavat vielä pienempiin virtuaalirastereihin ja niin edelleen, kunnes alin virtuaalirasteri viittaa varsinaisiin tiedostoihin. Tämän lähestymistavan syy on se, että ilman sitä myös käytettävät yleiskatselutasot kasvavat hyvin suuriksi. Huomaa, että tällaisen hierarkkisen rakenteen käyttö voi aiheuttaa joitakin artefakteja analyysiä ajettaessa, joten sitä kannattaa käyttää vain visualisointitarkoituksiin.
Laske tilastot valmiiksi virtuaalirastereillesi ja lähdetiedostoillesi. Tämä nopeuttaa tiedostojen avaamista esimerkiksi QGISissä. QGISin täytyy ottaa otos datan minimi- ja maksimiarvoista voidakseen asettaa väriskaalan oikein, ja tämä vie aikaa suurilla virtuaalirastereilla. Tämän välttämiseksi voit laskea tilastot etukäteen erilliseen XML-tiedostoon komennolla gdalinfo --stats.
Hyvä niksi QGISissä suurten rasterien kanssa työskenneltäessä on ottaa rasterityökalupalkki käyttöön (View->Toolbars->Raster Toolbar). Sen avulla voit helposti säätää väriskaalaa näytöllä näkyvän alueen mukaan, jolloin saat hyvän kontrastin zoomaustasosta riippumatta.
QGIS näyttää käsittelevän suuria aineistoja varsin hyvin, kun edellä mainitut vaiheet on tehty. Jopa koko Suomen 2 m:n korkeusmallilla zoomaus ja kartan siirtäminen ovat melko sujuvia.

Virtuaalirasterien luominen

Seuraavat työkalut tukevat virtuaalirasterien luomista:

GDAL gdalbuildvrt -komentorivityökalu.
GDAL gdal_translate STAC-haulla
Pythonissa ja R:ssä on kääreet GDALin gdalbuildvrtille; pidempi esimerkki R:lle löytyy StackOverflow-vastauksesta.
QGIS, GrassGIS ja SagaGIS tarjoavat graafisen käyttöliittymän gdalbuildvrtille
lidR tukee lidar-data-analyysin tulosten kirjoittamista suoraan virtuaalirasteriksi

Supertietokoneilla glalbuildvrt sisältyy kaikkiin GDALin sisältäviin moduuleihin, Pythonin BuildVRT geocondaan, QGIS QGIS-moduuliin sekä R:n terra ja lidR r-env -moduuliin.

Virtuaalirasterin luominen GDALin gdalbuildvrtillä

Jos datasi on jaettu alikansioihin, helpoin vaihtoehto on luoda virtuaalirasteri syötetiedostoluettelon avulla:

gdalbuildvrt -input_file_list file_list.txt virtual_raster.vrt

Katso lisävalinnat gdalbuildvrtin dokumentaatiosta, esimerkiksi no-data-arvon asettamiseen.

Tiedostoluettelon luominen gdalbuildvrtille

Tiedostoluettelon tulisi sisältää mieluiten täydet polut, mutta paikallisille tiedostoille voidaan käyttää myös suhteellisia polkuja.

Linux

find /dataset/project_2019680/mml/dem2m/ -name "*.**tif**" > file_list.txt

Windows

dir \data\dem2m\*.**tif** /S /B > file_list.txt

Rasteritiedostot Altaassa / jossakin muussa S3-palvelussa

Jos teet tämän supertietokoneella, lataa allas-moduuli.

Listaa tiedostonimet siinä muodossa kuin ne ovat ämpärissä käyttäen rclonea tai jotakin muuta työkalua:

rclone lsf --include '*.**tif**' allas:<your_bucket_name > file_list.txt

Lisää seuraavaksi tiedostoluetteloon täydet polut siinä muodossa kuin GDAL niitä vaatii käyttäen vsicurl-, vsis3- tai vsiswift-ajureita. Katso pidemmät selitykset GDAL-ajureista ja Altaasta CSC:n GDAL-sivulta.

sed -i -e 's-^-/**vsicurl**/https://a3s.fi/<your_bucket_name>/-' file_list.txt

sed -i -e 's-^-/**vsis3**/<your_bucket_name>/-' file_list.txt

sed -i -e 's-^-/**vsiswift**/<your_bucket_name>/-' file_list.txt

Määritä tunnistetietosi GDALille ennen kuin suoritat komennon gdalbuildvrt.

Virtuaalirasterin luominen STAC-haun tuloksista

STAC (Spatio-Temporal Asset Catalog) on tapa kuvata (rasteri)aineistoja siten, että dataa voidaan hakea ajan ja sijainnin perusteella. Esimerkiksi Paituli STAC sisältää useita suomalaisia aineistoja ja selittää STAC-käsitteitä tarkemmin.

GDAL tukee STAC-hakuja ja virtuaalirasterien luomista hakutulosten perusteella.

Esimerkiksi seuraava komento luo virtuaalirasterin Paituli STACin hakutulosten perusteella. Kysely pyytää tiedostoja corine_land_cover_at_geocubes-kokoelmasta määritellyllä aikavälillä ja alueella (bbox) käyttäen Assetteja (=tiedostoja), joiden nimi on COG.

Linuxissa ja Macissa:

gdal_translate "STACIT:\"https://paituli.csc.fi/geoserver/ogc/stac/v1/search?collections=corine_land_cover_at_geocubes&datetime=2017-01-05/2019-02-14&bbox=19.5,61.5,28.7,63.0\":asset=COG" -of VRT corine.vrt

Windowsissa (lisää ^ ennen &-merkkejä):

gdal_translate "STACIT:\"https://paituli.csc.fi/geoserver/ogc/stac/v1/search?collections=corine_land_cover_at_geocubes^&datetime=2017-01-05/2019-02-14^&bbox=19.5,61.5,28.7,63.0\":asset=COG" -of VRT corine.vrt

Huomaa, että GDAL poistaa virtuaalirasterista Assetit, jotka uudemmat Assetit peittävät kokonaan. Esimerkiksi Paituli STACissa CORINE-dataa on saatavilla vuosilta 2012 ja 2018. Jos yllä olevaa hakua muutetaan kattamaan myös vuosi 2012, STAC-haku löytäisi myös CORINE 2012 -Itemit, mutta niitä ei silti sisällytettäisi luotuun VRT:hen. Jos vuoden 2012 data tarvitaan mukaan virtuaalirasteriin, aikaväliä täytyy säätää niin, ettei vuosi 2018 sisälly siihen.