lidar_prosessointi

Tähän on koottu joitakin hyödyllisiä skriptejä LiDAR -aineistojen tehokkaaseen käsittelyyn ilmaisilla avoimen lähdekoodin työkaluilla (esim. PDAL).

Skriptit on tehty erityisesti suuria tiedostomääriä ajatellen - eli yhdellä komennolla käsitellään lähtökohtaisesti kaikki käsiteltävänä oleva aineisto. Lisäksi skriptit mahdollistavat usean tiedoston käsittelemisen samanaikaisesti, mikä parhaimmillaan nopeuttaa prosessointia huomattavasti. Prosessoinnin nopeus voi kuitenkin vaihdella huomattavasti käytetyn tietokoneen ominaisuuksien ja parametrien mukaan.

Tällä hetkellä skriptien valikoima on melko rajattu ja tehty lähinnä arkeologien tarpeita ajatellen. Valikoimaan lisätään mahdollisesti myöhemmin muita skriptejä, mutta voit yhtä hyvin jatkaa aineiston käsittelyä myös muilla ohjelmilla (esim. QGIS)

Käytöstä on pyritty antamaan mahdollisimman yksinkertaiset ja seikkaperäiset ohjeet, joita noudattamalla skriptejä on mahdollista käyttää ilman aiempaa kokemusta ohjelmoinnista tai komentorivin käytöstä. Myös skriptien käyttö on pyritty pitämään mahdollisimman yksinkertaisena ja yhtenäisenä.

Esivalmistelut

Python -ympäristön luominen

Skriptit tarvitsevat toimiakseen Python-ympäristön, johon on asennettu aineistojen käsittelyyn tarvittavat paketit (pdal, gdal, numpy, geopandas, rasterio). Alla on yksityiskohtiaset ohjeet python ympäristön luomiseen Anacondan avulla.

• Asenna Anaconda https://docs.anaconda.com/anaconda/install/
• Käynnistä Anaconda prompt
• Tee uusi Python ympäristö seuraavalla komennolla:

conda create -n lidar -c conda-forge pdal gdal geopandas numpy rasterio affine tqdm

Tällä komennolla Anaconda luo uuden Python ympäristön nimeltä 'lidar' ja asentaa siihen tarvittavat paketit (pdal, gdal, geopandas, numpy & tqdm). Hyväksy ympäristön luominen painamalla 'y'.

Kansiorakenne ja aineiston järjestäminen

Jotta skripti toimisi oikein, käsiteltävät tiedostot tulee järjestää kansioihin tietyllä tavalla. Skripti olettaa, että aineistot on järjestetty työskentelykansion sisällä erillisiin kansioihin aineiston tyypin perusteella seuraavasti:

• lidar (laserkeilausaineisto -laz -tiedostoina)
• dem (korkeusmallit geotiffeinä)

Varsinaiset skriptit (eli .py -päätteiset tiedostot) tallenetaan suoraan työskentelykansioon.

Skriptit käsittelevät oletuksena kaiken kansioihin tallennetun aineiston (jolla on oletettu tiedostopääte). Tästä syystä työskentelykansiota ei kannata käyttää aineistojen säilyttämiseen, eli kopioi työskentelykansioon vain se aineisto, minkä haluat käsitellä. Tuloksia sekä mahdollisia väliaikaistiedostoja varten tarvittavat kansiot luodaan lähtökohtaisesti automaattisesti.

Voit kopioida valmiin kansiorakenteen suoraan Githubista. Skriptit luovat tarvittaessa tuloksia ja väliaikaistiedostoja tarvittavat kansiot itse.

Skriptien suorittaminen

Kun aloitat skriptien käytön, aktivoi aina ensin python-ympäristö ja aseta työskentelykansio.

• Avaa Anaconda prompt (ts. Anacondan oma 'komentorivi') ja aktivoi aiemmin luomasi python-ympäristö komennolla:

conda activate lidar

• Aseta työskentelykansio komennolla cd C:\lisää\oikea\polku\tähän, esim:

cd C:\users\kayttajatunnus\lidar_prosessointi

Tämän jälkeen voit suorittaa skriptejä komennolla: python skriptin_nimi_tähän.py

Tarkemmat ohjeet skriptien käyttöön ja mahdolliset lisävalinnat on esitetty alla.

laz2dem.py - Pintamallien teko laserkeilausaineistosta

laz2dem.py tekee lidar/ -kansioon tallennetuista .laz päätteisistä tiedostoista pintamallit ja tallentaa ne kansioon dem/.

Skripti suoritetaan komennolla:

python laz2dem.py

Skriptin suorittamisen lopuksi skripti ilmoittaa käsittelyyn kuluneen ajan sekä yhden tiedoston käsittelyyn keskimäärin kuluneen ajan sekunteina. Voit käyttää ominaisuutta esimerkiksi asetusten vaikutusten testaamiseen.

laz2dem.py on hyödyntää laserkeilausaineiston käsittelyyn PDAL-kirjastoa: https://pdal.io/

Valinnaiset lisäparametrit

Skriptille voi antaa seuraavat valinnaiset parametrit:

• --buffer (default=0)
  • Laserkeilaustiilten käsittelyssä käytettävä bufferi. Käytettäessä bufferia pintamallin muodostamiseen käytetään myös käsiteltävää tiiltä ympäröivät pisteet, jolloin myös pintamalleista tulee alkuperäistä laserkeilaustiiltä laajempia. Bufferin käytöstä on hyötyä esimerkiksi tiettyjen visualisointitekniikoiden kanssa (esim. TPI), jotta vältytään poikkeamilta käsitetävien tiilien reunoilla. Bufferin koko annetaan metreissä. Bufferin käyttö hidastaa jonkin verran käsittelyä.
• --cores (default=4)
  • Rinnakkain käsiteltävien tiedostojen määrä. Tällä asetuksella voit vaikuttaa huomattavasti prosessoinnin nopeuteen. Rinnakkain käsiteltävien tiedostojen määrän tulisi olla korkeintaan yhtä suuri kuin tietokoneen loogisten suorittimien määrä, minkä voit tarkistaa esimerkiksi tehtävienhallinnasta. Huomaa kuitenkin, että tietokoneen ominaisuuksista riippuen maksimia pienempi rinnakkain käsiteltävien tiedostojen määrä voi olla kokonaisuutena nopeampi.
• --crs (default=3067, ts. ETRS-TM35FIN / EPSG:3067)
  • Tällä asetuksella voit asettaa haluamasi koordinaattijärjestelmän. Koordinaattijärjestelmän asettamiseen käytetään EPSG-numerokoodia.
• --resolution (default=1)
  • Laserkeilausaineistosta tuotettavien rasterien resoluutio metreinä. Käytä desimaalierottimena pistettä.

Lisäparametrit annetaan varsinaisen komennon jälkeen, esim:

python pdal_laz2dem.py --buffer=30 --cores=8 --resolution=0.5

dem2rvt.py - Pintamallien visualisointi Relief Visualization Toolboxin (RVT) avulla

dem2rvt.py tarjoaa yksinkertaisen ja tehokkaan käyttöliittymän korkeusmallien visualisointiin Relief Visualization Toolboxin (RVT) avulla. RVT tarjoaa kattavan valikoiman erilaisia visualisointitekniikoita.

• Lisätietoa RVT:stä ja RVT:n avulla tehtävistä visualisoinneista: https://rvt-py.readthedocs.io/en/latest/index.html
• RVT Python library, Žiga Kokalj, Žiga Maroh, Krištof Oštir, Klemen Zakšek and Nejc Čož, 2022. (ZRC SAZU and University of Ljubljana)

Esimerkkejä dem2rvt.py käytöstä:

Yksittäisen visualisoinnin tekeminen

python dem2rvt.py --visualisoinnit=hillshade

Komentoa voi lyhentää käyttäen --visualisoinnit sijaan lyhennettä -v

python dem2rvt.py -v=hillshade

Useampien visualisointien tekeminen yhdellä komennolla on mahdollista käyttäen erottimena puolipistettä

python dem2rvt.py --visualisoinnit=slope;hillshade

Avainsanalla 'kaikki' voi tehdä yhdellä komennonnolla kaikki RVT:n mahdollistamat visualisoinnit

python dem2rvt.py --visualisoinnit=kaikki

Lisäparametrit

Osana komentoa skriptille voi antaa valinnaisia lisäargumentteja. Kaikkiin visualisointeihin vaikuttavia lisäargumentteja ovat:

• --vrt
  • Määrittää tehdäänkö visualisoinnista lopuksi virtuaalirasteri
  • Virtuaalirasteri mahdollistaa suuren tiedostomäärän käsittelyn yhtenä tasona paikkatieto-ohjelmissa.
  • Vaihtoehdot: True -> tehdään virtuaalirasteri / False -> ei tehdä virtuaalirasteria
  • Oletuksena True
• --vrt_pyramidit (default=True)
  • Jos visualisoinnille tehdään virtuaalirasteri, --vrt_pyramidit määrittää tehdäänkö virtuaalirasterista pyramidit, eli matalaresoluutioiset esikatselukuvat.
  • Pyramidit sujuvoittavat etenkin suurempien aineistojen käyttöä.
• --gdal_datatype
  • määrittää tallennettavien visualisointien datatyypin käyttäen gdal numerokoodia (oletuksena 6, eli float32)

Lisäparametreja käytetään seuraavasti

#Tekee vinovalovarjosteen, mutta ei laske lopuksi virtuaalirasteria
python dem2rvt.py --visualisoinnit=kaikki --vrt=False
#Tekee vinovalovarjosteen ja sille virtuaalirasterin, mutta ei pyramideja
python dem2rvt.py --visualisoinnit=kaikki --vrt_pyramidit=False 

Vinovalovarjoste / hillshade

Laskee korkeusmallista vinovalovarjosteen ja tallentaa tuloksen kansioon hillshade/

python dem2rvt.py --visualisoinnit=hillshade

Lisäparametrit

Vinovalovarjosteen lisäparametrit ja oletusasetukset

• --ve_factor
  • korkeuserojen korostamiseen käytettävä kerroin. Oletuksena 3
• --sun_azimuth
  • Valonlähteen suunta asteina. Oletuksena 315
• --sun_elevation
  • Valonlähteen korkeus asteina. Oletuksena 45

Esimerkikki lisäparamterien käytöstä:

python dem2rvt.py --visualisoinnit=hillshade --ve_factor=2 --sun_azimuth=115 --sun_elevation=20

Monisuuntainen vinovalovarjoste / multiple direction hillshade (mdhs)

Laskee korkeusmallista useasta suunnasta valaistun vinovalovarjosteen ja tallentaa sen rgb-rasterina kansioon multi_hillshade/

python dem2rvt.py --visualisoinnit=multi_hillshade
python dem2rvt.py --visualisoinnit=mdhs

Lisäparametrit

Monisuuntaisen vinovalovarjosteen lisäparametrit ja oletusasetukset

• --ve_factor
  • korkeuserojen korostamiseen käytettävä kerroin. Oletuksena 1
• --sun_elevation
  • Valonlähteen korkeus asteina. Oletuksena 45
• --nr_directions
  • valaisusuuntien määrä, oletuksena 3

Rinteenkaltevuus / slope

Laskee korkeusmallista rinteenkaltevuuden (slope) ja tallentaa sen kansioon slope/

python dem2rvt.py --visualisoinnit=slope

Simple local relief model (slrm)

Laskee korkeusmallista simple relief modelin (slrm) ja tallentaa sen kansioon simple_relief_model. Slrm lasketaan vertaamalla jokaisen solun korkeutta sitä ympäröivien solujen korkeuksien keskiarvoon halutulla säteellä.

python dem2rvt.py --visualisoinnit=simple_local_relief_model
python dem2rvt.py --visualisoinnit=slrm

Multi-scale relief model (msrm)

Laskee korkeusmallista multi-scale relief modelin (msrm) ja tallentaa sen kansioon multi-scale_relief_model.

python dem2rvt.py --visualisoinnit=multi-scale_relief_model
python dem2rvt.py --visualisoinnit=msrm

Sky-view factor (svf)

Laskee korkeusmallista sky-view factorin tai anisotropic sky-view factorin.

Sky-view factor tarkempi kuvaus ja viittaukset: Zakšek, K., Oštir, K., Kokalj, Ž. 2011. Sky-View Factor as a Relief Visualization Technique. Remote Sensing 3: 398-415 https://doi.org/10.3390/rs3020398

#sky-view factor
python dem2rvt.py --visualisoinnit=sky-view_factor
python dem2rvt.py --visualisoinnit=svf
#anisotropic sky-view factor
python dem2rvt.py --visualisoinnit=anisotropic_sky-view_factor
python dem2rvt.py --visualisoinnit=asvf

Topographic openness

Laskee pintamallista topographic positiivisen- ja negatiivisen topographic opennessin ja tallentaa ne erillisiin kansioihin.

python dem2rvt.py --visualisoinnit=openness

Multi-scale topographic position (mstp)

Laskee pintamallista multi-scale topographic positionin (mstp).

python dem2rvt.py --visualisoinnit=multi-scale_topographic_position
python dem2rvt.py --visualisoinnit=mstp

Local dominance

Laskee pintamallista local dominancen.

python dem2rvt.py --visualisoinnit=dominance

dem2tpi.py - Topographic position index

dem2tpi.py laskee dem/ -kansioon tallennetuista .tif päätteisistä korkeusmalleista topographic position indexin (TPI) ja tallentaa sen uuteen kansioon. Uusi kansio nimetään TPI:n laskennassa käytettyjen parametrien perusteella.

TPI on hyödyllinen tekniikka lähiympäristöön korkeampien tai matalampien maastonmuotojen visualisointiin. Topographic position index lasketaan vertaamalla korkeusmallin jokaisen solun korkeutta sen ympäristön solujen korkeuksien keskiarvoon. Lähiympäristön tarkasteluun käytettävää sädettä muuttamalla TPI:n avulla voidaan korostaa hyvin eri kokoisia maastonmuotoja. Esimerkiksi 5 metrin sädettä käyttämällä TPI korostaa paikallisia mittakaavaltaan pieniä maastonmuotoja, mutta ei juurikaan reagoi suurempiin maastonmuotoihin. Vastaavasti suurempaa sädettä käytettäessä TPI korostaa mittakaavaltaan suurempia maastonmuotoja.

Skripti suoritetaan komennolla:

python dem2tpi.py --radius=7

Pakolliset parametrit

• --radius
  • määrittää TPI:n laskentaan käytettävän säteen metreinä. Esimerkiksi arkeologisesti kiinnostavien pinnanmuotojen visualisointiin radius on tyypillisesti hyvä asettaa välille 5-15.
  • Huomioithan, että suurempi radius kasvattaa prosessointiin kuluvaa aikaa.
  • HUOM! Tämä versio skriptistä soveltuu lähinnä suhteellisen pienten ja paikallisten korkeuserojen visualisointiin, koska tulokset skaalataan vakioidulla kaavalla tietylle vaihteluvälille. Muihin tarpeisiin on parempi käyttää RVT:n simple local relief model -visualisointia, joka on käytännössä sama asia.

Skriptille voi antaa seuraavat valinnaiset parametrit:

• --cores (default=4)
  • Rinnakkain käsiteltävien tiedostojen määrä. Tällä asetuksella voit vaikuttaa huomattavasti prosessoinnin nopeuteen. Rinnakkain käsiteltävien tiedostojen määrän tulisi olla korkeintaan yhtä suuri kuin tietokoneen loogisten suorittimien määrä, minkä voit tarkistaa esimerkiksi tehtävienhallinnasta. Huomaa kuitenkin, että tietokoneen ominaisuuksista riippuen maksimia pienempi rinnakkain käsiteltävien tiedostojen määrä voi olla kokonaisuutena nopeampi.
• --crs (default=3067, ts. ETRS-TM35FIN / EPSG:3067)
  • Tällä asetuksella voit asettaa haluamasi koordinaattijärjestelmän. Koordinaattijärjestelmän asettamiseen käytetään EPSG-numerokoodia.

Skriptin suorittamisen lopuksi skripti ilmoittaa käsittelyyn kuluneen ajan sekä yhden tiedoston käsittelyyn keskimäärin kuluneen ajan sekunteina. Voit käyttää ominaisuutta esimerkiksi asetusten vaikutusten testaamiseen.

dem2tpi.py perustuu Zoran Čučkovićin numpy-kirjastoa hyödyntävään tapaan laskea TPI (https://landscapearchaeology.org/2021/python-tpi/). Oleellisimpana erona alkuperäiseen tässä esitetty versio skriptistä mahdollistaa suoraan usean tiedoston rinnakkaisen käsittelyn ja skaalaa tulokset vakioidulla kaavalla kokonaisluvuiksi väliltä +/- 100.