/
_06-ex.Rmd
134 lines (103 loc) · 4.16 KB
/
_06-ex.Rmd
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
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
Some of the exercises use a vector (`zion_points`) and raster dataset (`srtm`) from the **spDataLarge** package.
They also use a polygonal 'convex hull' derived from the vector dataset (`ch`) to represent the area of interest:
```{r 06-ex-e0, message=FALSE, include=TRUE}
library(sf)
library(terra)
library(spData)
zion_points_path = system.file("vector/zion_points.gpkg", package = "spDataLarge")
zion_points = read_sf(zion_points_path)
srtm = rast(system.file("raster/srtm.tif", package = "spDataLarge"))
ch = st_combine(zion_points) %>%
st_convex_hull() %>%
st_as_sf()
```
E1. Crop the `srtm` raster using (1) the `zion_points` dataset and (2) the `ch` dataset.
Are there any differences in the output maps?
Next, mask `srtm` using these two datasets.
Can you see any difference now?
How can you explain that?
```{r 06-ex-e1}
plot(srtm)
plot(st_geometry(zion_points), add = TRUE)
plot(ch, add = TRUE)
srtm_crop1 = crop(srtm, vect(zion_points))
srtm_crop2 = crop(srtm, vect(ch))
plot(srtm_crop1)
plot(srtm_crop2)
srtm_mask1 = mask(srtm, vect(zion_points))
srtm_mask2 = mask(srtm, vect(ch))
plot(srtm_mask1)
plot(srtm_mask2)
```
E2. Firstly, extract values from `srtm` at the points represented in `zion_points`.
Next, extract average values of `srtm` using a 90 buffer around each point from `zion_points` and compare these two sets of values.
When would extracting values by buffers be more suitable than by points alone?
```{r 06-ex-e2}
zion_points_buf = st_buffer(zion_points, dist = 90)
plot(srtm)
plot(st_geometry(zion_points_buf), add = TRUE)
plot(ch, add = TRUE)
zion_points_points = extract(srtm, vect(zion_points))
zion_points_buf = extract(srtm, vect(zion_points_buf))
plot(zion_points_points$srtm, zion_points_buf$srtm2)
```
E3. Subset points higher than 3100 meters in New Zealand (the `nz_height` object) and create a template raster with a resolution of 3 km for the extent of the new point dataset.
Using these two new objects:
- Count numbers of the highest points in each grid cell.
- Find the maximum elevation in each grid cell.
```{r 06-ex-e3}
nz_height3100 = dplyr::filter(nz_height, elevation > 3100)
new_graticule = st_graticule(nz_height3100, datum = "EPSG:2193")
plot(st_geometry(nz_height3100), graticule = new_graticule, axes = TRUE)
nz_template = rast(ext(nz_height3100), resolution = 3000, crs = crs(nz_height3100))
nz_raster = rasterize(vect(nz_height3100), nz_template,
field = "elevation", fun = "length")
plot(nz_raster)
plot(st_geometry(nz_height3100), add = TRUE)
nz_raster2 = rasterize(vect(nz_height3100), nz_template,
field = "elevation", fun = max)
plot(nz_raster2)
plot(st_geometry(nz_height3100), add = TRUE)
```
E4. Aggregate the raster counting high points in New Zealand (created in the previous exercise), reduce its geographic resolution by half (so cells are 6 by 6 km) and plot the result.
- Resample the lower resolution raster back to the original resolution of 3 km. How have the results changed?
- Name two advantages and disadvantages of reducing raster resolution.
```{r 06-ex-e4}
nz_raster_low = raster::aggregate(nz_raster, fact = 2, fun = sum, na.rm = TRUE)
res(nz_raster_low)
nz_resample = resample(nz_raster_low, nz_raster)
plot(nz_raster_low)
plot(nz_resample) # the results are spread over a greater area and there are border issues
plot(nz_raster)
```
```{asis}
Advantages:
- lower memory use
- faster processing
- good for viz in some cases
Disadvantages:
- removes geographic detail
- adds another processing step
```
E5. Polygonize the `grain` dataset and filter all squares representing clay.
```{r 06-ex-e5}
grain = rast(system.file("raster/grain.tif", package = "spData"))
```
- Name two advantages and disadvantages of vector data over raster data.
- When would it be useful to convert rasters to vectors in your work?
```{r 06-ex-e5-2}
grain_poly = as.polygons(grain) %>%
st_as_sf()
levels(grain)
clay = dplyr::filter(grain_poly, grain == "clay")
plot(clay)
```
```{asis}
Advantages:
- can be used to subset other vector objects
- can do affine transformations and use sf/dplyr verbs
Disadvantages:
- better consistency
- fast processing on some operations
- functions developed for some domains
```