153ZODH / 15. cvičení

← předchozí cvičení • stránky předmětu

Objektově orientovaná klasifikace

Osnova

Toto cvičení je zaměřeno na možnosti objektově orientované klasifikace družicových dat. V tomto ohledu navazuje na cvičení věnované klasifikaci, řizené klasifikaci a především segmentaci obrazu.

Seznam příkazů

• i.group
• g.region
• i.segment
• r.to.vect
• v.db.dropcolumn
• v.db.addcolumn
• v.db.update
• v.rast.stats
• r.texture
• v.to.db
• v.db.select
• v.extract
• v.class.mlpy
• v.to.rast

Příprava

Nejprve nastavíme výpočetní region (g.region). Vzhledem k tomu, že je proces objektově orientované klasifikace poměrně výpočetně náročný, zvolíme menší region. V našem případě předpokládáme existenci již definované regionu uloženého jako 'ceske-budejovice-10km' (jde o region s offsetem 10km kolem města České Budějovice, viz návod).

 g.region region=ceske-budejovice-10km -p

projection: 1 (UTM)
zone:       33
datum:      wgs84
ellipsoid:  wgs84
north:      5439839.30457264
south:      5410176.48726338
west:       445117.66335159
east:       480328.80428671
nsres:      30.02309444
ewres:      29.99245395
rows:       988
cols:       1174
cells:      1159912

Před dalším výpočtem je vhodné obrazové kanály družicové scény, které budou vstupovat do procesu klasifikace předzpracovat (aplikovat geometrické a atmosferické korekce). Minimálně nahradit místa bez informace hodnotou NULL - viz návod zde.

Dále je potřeba vytvořit obrazovou skupinu (i.group) slouženou ze vstupních kanálů - vynechán je panchromatický kanál 'B6'.

 i.group group=B input=B1,B2,B3,B4,B5,B7

Ve wxGUI je pro vytvoření skupiny dostupný specializovaný interaktivní dialog (viz obr. níže).

Založení obrazové skupiny ve wxGUI

Segmentace a statistika obrazových kanálů

Segmentaci obrazu provedeme příkazem i.segment.

i.segment group=B out=B_segment thresh=0.02 minsize=15

Výsledek segmentace obrazu jako přípravy pro objektově orientovanou klasifikaci

Vygenerované segmenty převedeme do vektorové reprezentace modulem r.to.vect a přidáme do atributové tabulky nové sloupce 'class' a 'brightness', zároveň odstraníme nepotřebný sloupec 'label'. Nepotřebné sloupečky odstraníme z atributové tabulky příkazem v.db.dropcolumn, nové potom v.db.addcolumn. Hodnoty nastavíme příkazem v.db.update.

# vektorizace segmentů
r.to.vect -v B_segment out=B_segment type=area

# modifikace atributové tabulky
v.db.dropcolumn B_segment col=label
v.db.addcolumn B_segment col="class int"
v.db.addcolumn B_segment col="brightness double precision"
v.db.update B_segment col=brightness value=0

Dále do atributové tabulky uložíme statistiku segmentů spojenou s jednotlivými kanály obrazové scény. Příklad pro první kanál:

v.rast.stats vector=B_segment raster=B1 column=B1

Jas ('brighness') vypočteme jako sumu středních hodnot jednotlivých kanálů:

v.db.update B_segment column=brightness qcolumn="brightness+B1_mean"

Po výpočtu můžeme nepotřebné sloupečky z atributové tabulky odstranit.

v.db.dropcolumn B_segment col=B1_n,B1_min,B1_max,B1_range,B1_cf_var,B1_sum

Statistiku můžeme spočítat dávkově pomocí grafického modeleru anebo pomocí jednoduchého skriptu v jazyce Python:

import grass.script as grass

maps = grass.read_command('i.group', group='B', flags='g', quiet=True).splitlines()

for map in maps:
    band = map.split('@', 1)[0]
    grass.message("Processing <%s>..." % band)
    grass.run_command('v.rast.stats', vector='B_segment',
                      raster=band, column=band)
    grass.run_command('v.db.update', map='B_segment',
                      column='brightness', qcolumn='brightness+%s_mean' % band)
    grass.run_command('v.db.dropcolumn', map='B_segment',
                      column='{band}_n,{band}_min,{band}_max,' \
                          '{band}_range,{band}_cf_var,{band}_sum'.format(band=band))

Ve výsledku atributová tabulka obsahuje pro každý obrazový kanál informaci o střední hodnotě ('mean'), směrodatné odchylce ('stddev') a varianci nebo-li rozptylu ('variance') např.

B1_mean|72.6521739130435
B1_stddev|3.0447202434803
B1_variance|9.27032136105877
B2_mean|32.6086956521739
B2_stddev|1.99432085933559
B2_variance|3.97731568998105
...

Analýza textury panchromatického kanálu

Analýzu textury provedeme příkazem r.texture:

r.texture input=B6 prefix=B6_text method=entr,idm

Jedné metodě odpovídá jedna rastrová mapa, v našem případě vzniknou dvě rastrové mapy: 'B6_text_IDM' a 'B6_text_Entr'.

Analýza textury panchromatického kanálu (vlevo metoda IDM, vpravo metoda ENTR)

Na základě těchto dvou rastrů aktualizujeme atributovou tabulku segmentů, zároveň odstraníme nepotřebné sloupečky. Zůstanou pouze: 'mean', 'stddev' a 'variance'. Statistiku do atributové tabulky uložíme příkazem v.rast.stats, nepotřebné sloupečky odstraníme v.db.dropcolumn.

# pro texturu IDM
v.rast.stats vector=B_segment raster=B6_text_IDM col=B6_text_IDM
v.db.dropcolumn B_segment col=B6_text_IDM_n,B6_text_IDM_min,B6_text_IDM_max,B6_text_IDM_range,B6_text_IDM_cf_var,B6_text_IDM_sum

# pro texturu Entr
v.rast.stats vector=B_segment raster=B6_text_Entr col=B6_text_Entr
v.db.dropcolumn B_segment col=B6_text_Entr_n,B6_text_Entr_min,B6_text_Entr_max,B6_text_Entr_range,B6_text_Entr_cf_var,B6_text_Entr_sum

Geometrické vlastnosti segmentů obrazu

V dalším kroku uložíme do atributové tabulky informace spojené s geometrií segmentů: jejich plochu ('area'), obvod ('perimeter'), kompaktnost ('compact') a fraktálovou dimenzi ('fd'). Sloupečky přidáme příkazem v.db.addcolumn anebo pomocí správce atributových dat, odvozené hodnoty z geometrie nahrajeme do atributové tabulky příkazem v.to.db.

# přidání nových sloupečků do atributové tabulky
v.db.addcolumn B_segment col="area int, perimeter int, compact double precision, fract_dim double precision" --q

# výpočet hodnot odvozených z geometrie segmentů
v.to.db B_segment op=area col=area 
v.to.db B_segment op=perimeter col=perimeter 
v.to.db B_segment op=compact col=compact 
v.to.db B_segment op=fd col=fract_dim

Výsledek si ověříme (v.db.select nebo ve správci atributových dat) jednoduchým atributovým dotazem - vybereme segmenty s výměrou větší než 3 km²:

v.db.select B_segment col=cat,area,perimeter,compact,fract_dim where="area > 3e6"

cat|area|perimeter|compact|fract_dim
645628|3866602|35710|5.122953|1.382287
885913|3072390|37568|6.046235|1.410357

Volba trénovacích ploch

Volbu trénovacích ploch (resp. segmentů) provedeme v digitalizačním nástroji přiřazením příslušného kódu reprezentujího danou třídu (atribut 'class').

Příklad volby trénovacích ploch na základě barevné syntézy 453

Trénovací segmenty vybereme a uložíme do nové vektorové mapy příkazem v.extract:

v.extract input=B_segment where="class > 0" output=tr_areas

Nastavíme vhodnou tabulku barev tak aby byly jednotlivé trénovací plochy dobře rozlišitelné. Interaktivní nástroj pro nastavení tabulky barev je k dispozici z menu wxGUI Vector → Manage colors → Manage color rules interactively.

Nastavení tabulky barev pro vektorovou mapu trénovacích ploch: (1) vodní plochy (2) les (3) zástavba (4) zemědělská půda

Vektorová mapa trénovacích ploch s nastavenou tabulkou barev

Aplikace klasifikátoru

V našem případě použijeme klasifikátor DLDA implementovaný v rámci modulu v.class.mlpy.

Poznámka: Tento modul vyžaduje knihovnu mlpy.

v.class.mlpy input=B_segment training=tr_areas

Na závěř převedeme výsledek klasifikace do rastrové reprezentace pomocí v.to.rast:

v.to.rast input=B_segment output=B_obj_class use=attr attrcolumn=class

Výsledek objektově orientované včetně legendy