opgaande_elementen.py 1.36 KB
Newer Older
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
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
"""
Short script to combine opgaande elementen kaarten from 2018 and 2019 into a single map showing all possible
combinations.

First reclass Open (1), close (2) cells to powers of 2:

2018 Open  0
2018 Close 1
2019 Open  2
2019 Close 4

Then make all combinations:

      2019O 2019C
2018O     2     4
2018C     3     5

Opgaande elementenkaarten by Henk Meeuwsen.
Zie W:\PROJECTS\Landschapsmonitor\Documentatie\LM werkwijze openheid 2020.docx.

Code by Hans Roelofsen, WEnR, 5 August 2020
"""

import os
import numpy as np
import rasterio as rio

from sample.utils import lm as lm

# Read source rasters
og18_src = r'c:\apps\temp_geodata\gesloten_2018_20200629time114900.flt'
og19_src = r'c:\apps\temp_geodata\gesloten_2019_20200702time193947.flt'

og18 = rio.open(og18_src)
og19 = rio.open(og19_src)

# Apply reclass dictionaries
og18_rc_dict = {-9999: 0, 1: 0, 2: 1}
og19_rc_dict = {-9999: 0, 1: 2, 2: 4}

og18_reclass = lm.vec_translate(og18.read(1), og18_rc_dict)
og19_reclass = lm.vec_translate(og19.read(1), og19_rc_dict)

# Add two arrays together
og1819 = np.add(og18_reclass, og19_reclass)
print(np.unique(og1819))

# write to file
basename = 'lm_gesloten_18-19_combi'
prof = og18.profile
prof.update(dtype=rio.uint8, driver='gtiff')
with rio.open(os.path.join(r'c:/apps/temp_geodata', '{}.tif'.format(basename)), 'w', **prof) as dst:
    dst.write(og1819.astype(np.uint8), 1)