Danube Hack a pôdna mapa z WMS do vektorů

Tento víkend jsem se účastnil jako mentor akce DanubeHack a měl jsem tam možnost řešit celou řadu zajímavých úloh. Někdy se snad dostanu ke svojí malé aplikaci, která vyroběla api postavené nad Djangem a RedHatím OpenShiftem. Dneska bych se s vámi rád podělil o to, jak jsme převáděli rastrovou mapu staženou z WMS na vektory pomocí GRASSu.

Půdní mapu si vyžádal jeden z týmů a našli dokonce i funkční WMS z Pôdneho portálu Slovenska.

Půdní mapa Slovenska. Zdroj: http://www.podnemapy.sk/portal/prave_menu/podna_mapa/podna_mapa.aspx

Nejdřív jsem si založil GRASS LOCATION a MAPSET a nastavil na hranice Slovenska. Potom jsem stáhnul data z WMS služby – použil jsem formát GeoTIFF, protože JPEG díky ztrátové kompresy potrhá barvy v pixelech a místo rozumných cca 20 kategorií by jich tam bylo asi 3000. Samozřejmě jsem nastavil použitelné rozlišení regionu na 50m (původní mapa je 1:400~000):

GRASS> g.region -p
projection: 99 (Krovak)
zone:       0
datum:      hermannskogel
ellipsoid:  bessel
north:      -1132700
south:      -1334770
west:       -591450
east:       -165430
nsres:      50.00494927
ewres:      50.00234742
rows:       4041
cols:       8520
cells:      34429320

a stáhnu data pomocí r.in.wms

GRASS> r.in.wms url="http://sscri.vupop.sk/arcgis/services/vupop_wms/MapServer/WMSServer" wms_version="1.3.0"  srs=102067 format=tiff layers=0 out=podnamapa --o --v

Pudní mapa z WMS

Jak můžete vidět, mapa je ošklivá – plná kartografických šraf a nápisů – nic, co byste mohli převézt na rozumně vypadající vektory. Přemýšlel jsem nad vhodným filtrem a nakonec z toho vypadl modul r.neighbors, který prožene přes raster okno veliké MxM pixelů, jako operátor jsem použil “nejčastěji se vyskytující hodnotu” a velikost okna zkoušel, až jsem došel k 9×9:

GRASS> r.neighbors in=podnamapa out=podnamapa_clean method=mode size=9 --o

r.neighboors na šrafy

Tím jsem se zbavil většiny šraf a nápisů (hranice se při tomhle měřítku taky trochu šouply, ale při tomhle měřítku a času, kdo by to řešil?) Pár jich ale zbylo. Při vypsání počtu pixelů za jednotlivé kategorie jsem prostě odhadnul, že ta malá čísla jsou ty špatné hodnoty – pixely, které tam nemají co dělat a které neprošly filtrem:

GRASS> r.report map=podnamapa_clean units=k

+-----------------------------------------------------------------------------+
|                         RASTER MAP CATEGORY REPORT                          |
|LOCATION: sk                                         Mon Oct 19 22:30:41 2015|
|-----------------------------------------------------------------------------|
|          north:    -1132700    east:     -165430                            |
|REGION    south:    -1334770    west:     -591450                            |
|          res:   50.00494927    res:  50.00234742                            |
|-----------------------------------------------------------------------------|
|MASK: hranice in danubehack, categories 1                                    |
|-----------------------------------------------------------------------------|
|MAP: 9x9 neighborhood: mode of podnamapa2 (podnamapa2_clean in danubehack)   |
|-----------------------------------------------------------------------------|
|                      Category Information                        |  square  |
|    #|description                                                 |kilometers|
|-----------------------------------------------------------------------------|
|    0| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|     5.753|
| 6342| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|     0.010|
| 8126| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|   739.818|
| 9664| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|     0.615|
|10847| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|  3828.899|
|13113| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|    49.410|
|13741| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|     0.690|
|14888| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|     0.050|
|16405| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|     0.105|
|17056| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|     0.020|
|18105| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|10,297.153|
|20250| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|11,647.305|
|20312| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|  2813.458|
|24319| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|   901.181|
|24570| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|  4829.897|
|25427| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|     2.243|
|26005| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|    68.512|
|26047| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|  1369.785|
|27482| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|  3191.988|
|32512| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|   168.625|
|32762| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|  5694.238|
|32767| . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .|  3409.603|
|    *|no data. . . . . . . . . . . . . . . . . . . . . . . . . . .|37,066.504|
|-----------------------------------------------------------------------------|
|TOTAL                                                             |86,085.861|
+-----------------------------------------------------------------------------+

Pomocí r.reclass jsem se jich zbavil:

GRASS> r.reclass in=podnamapa_clean out=podnamapa_clean2 rules=- << EOF
0 = NULL
6342 = NULL
8126 = 8126
9664 = NULL
10847 = 10847 
13113 = 13113 
13741 = NULL
14888 = NULL
16405 = NULL
17056 = NULL
18105 = 18105 
20250 = 20250 
20312 = 20312 
24319 = 24319 
24570 = 24570 
25427 = 25427 
26005 = 26005 
26047 = 26047 
27482 = 27482 
32512 = 32512 
32762 = 32762 
32767 = 3767
* = NULL
EOF

Převod raster na vektor je už hračka

GRASS> r.to.vect in=podnamapa_clean2 out=podnamapa type=area

Půdní mapa vektory

a aby nebyly na polygonech oškolivé pixely, tak jsem to ještě prohnal generalizací

GRASS> v.generalize in=podnamap out=podnamapa_generalize --o method=douglas threshold=100

Generalizovaná mapa

Výsledná mapa nemá kategorie, to si ale už účastníci Hackathonu doplnili z původního datového zdroje ručně.

podnamapa_vector

A že jste to vy, vrazil jsem to ve formátu TopoJSON na GitHub, abyste ti to mohli prohlídnout a stáhnout (doplnění kategorií za domácí úkol).

Výsledek mi zabral asi tak hodinu zkoušení. Má to chyby: např. díry po té reklasifikaci, že, ale šlo o čas a “dostatečně dobré” řešení.

Post navigation