Přístup k rastrovým datům

Přístup k rastrovým datům umožňuje PyGRASS ve třech režimech:

  • RasterRow (náhodné čtení po řádcích, sekvenční zápis)
  • RasterRowIO (čtení po řádcích z vyrovnávací paměti, sekvenční zápis)
  • RasterSegment (náhodné čtení a zápis po dlaždicích)

Další informace v dokumentaci PyGRASS.

Varování

GRASS při čtení rastrových dat vždy data převzorkuje podle aktuálního výpočetního regionu. Manipulaci s regionem má v PyGRASS na starost třída Region anebo lze přímo použít modul g.region.

Statistika rastrových dat

V následující ukázce vypíšeme statistiku rastru:

  1. Před načtením dat je nastaven výpočetní region (řádek 11-13).
  2. Rastrová data jsou načtena pomocí třídy RasterRow (řádek 15-16).
  3. Jednotlivé řádky a sloupce rastru jsou procházeny cyklem for (řádky 20-21).
  4. Na konci skriptu nezapomeneme rastrovou mapu korektně uzavřít 34.
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
#!/usr/bin/env python3

import numpy as np
from grass.pygrass.raster import RasterRow
from grass.pygrass.modules import Module

from grass.pygrass.gis.region import Region

name = 'dmt@PERMANENT'

reg = Region()
reg.from_rast(name)
reg.set_current()

rast = RasterRow(name)
rast.open('r')

min = max = None
count = ncount = 0
for row in rast:
    for value in row:
        if np.isnan(value):
            ncount += 1
        else:
            if min is None:
                min = max = value
            else:
                if min > value:
                    min = value
                elif max < value:
                    max = value
        count += 1

rast.close()

print("min={:.2f} max={:.2f} count={} (no-data: {})".format(
    min, max, count, ncount)
)

Skript ke stažení zde.

Výpis může vypadat následovně:

min=53.80 max=1530.51 count=138116 (no-data: 59244)

Poznámka

Tento skript berte jako ilustrační, rozhodně jej nelze považovat za optimální cestu pro zjištění extremních hodnot v rastru. Porovnejte s modulem r.univar a verzí skriptu založené na knihovně NumPy (17-18).

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
#!/usr/bin/env python3

import numpy as np
from grass.pygrass.raster import RasterRow
from grass.pygrass.gis.region import Region

name = 'dmt@PERMANENT'

reg = Region()
reg.from_rast(name)
reg.set_current()

with RasterRow(name) as rast:
    array = np.array(rast)

print("min={:.2f} max={:.2f} count={} (no-data: {})".format(
      array.min(), array.max(), array.size,
      np.count_nonzero(np.isnan(array)))
)

Skript ke stažení zde.

Dotazování na rastrová data

Skript vypisuje pro definiční body obcí v ČR jejich nadmořské výšky odvozené z digitálního modelu terénu (rastrová mapa dmt).

  1. Před načtením rastrových dat na řádcích 10-12 je na základě rastrové mapy dmt nastaven výpočetní region.
  2. Rastrová mapa dmt je načtena třídou RasterRow (řádka 14-15).
  3. Jelikož se jedná u vstupní vektorové mapy o data bodová, tak stačí mapu otevřít bez topologie (řádky 17-18).
  4. Souřadnice definičních bodů obcí jsou převedeny na souřadnice rastru funkcí coor2pixel (řádek 21)
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
#!/usr/bin/env python3

from grass.pygrass.raster import RasterRow
from grass.pygrass.vector import Vector
from grass.pygrass.gis.region import Region
from grass.pygrass.utils import coor2pixel

name = 'dmt@PERMANENT'

reg = Region()
reg.from_rast(name)
reg.set_current()

dmt = RasterRow(name)
dmt.open('r')

obce = Vector('obce_bod')
obce.open('r')

for o in obce:
    x, y = coor2pixel(o.coords(), region)
    value = dmt[int(x)][int(y)]
    print ('{:40}: {:.0f}'.format(o.attrs['nazev'], value))

obce.close()
dmt.close()

Skript ke stažení zde.

Výpis může vypadat následovně:

...
Kopidlno                                : 225
Neratov                                 : 223
Podhorní Újezd a Vojice                 : 336
...

Poznámka

Rychlost implementace můžete porovnat s modulem v.what.rast.

v.what.rast -p map=obce_bod@ruian raster=dmt@PERMANENT