introduction

The following exercises are modified from Chapters 3, 4, and 5 of Geocomputation with R by Rovin Lovelace.

prerequisites

library(sf)
library(terra)
library(spData)
library(geodata)
library(spDataLarge)
library(tidyverse)
dem = rast(system.file("raster/dem.tif", package = "spDataLarge"))
landsat = rast(system.file("raster/landsat.tif", package = "spDataLarge"))
spain_dem = geodata::elevation_30s(country = "Spain", path = ".", mask = FALSE)
srtm = rast(system.file("raster/srtm.tif", package = "spDataLarge"))

Exercise 1: Plot the histogram and the boxplot of the dem.tif file from the spDataLarge package

hist(dem)

Exercise 2: Manipulate rasters

reclassify the elevation in three classes: low (<300), medium and high (>500).

plot(dem)

rcl = matrix(c(-Inf, 300, 0, 300, 500, 1, 500, Inf, 2), ncol = 3, byrow = TRUE)

dem_rcl = classify(dem, rcl = rcl)

levels(dem_rcl) = data.frame(id = 0:2, cats = c("low", "medium", "high"))

plot(dem_rcl)

compute the mean elevation for each altitudinal class

elevation_mean = zonal(dem, dem_rcl, fun = "mean")
elevation_mean
##     cats      dem
## 1    low 274.3910
## 2 medium 392.0486
## 3   high 765.2197

Calculate the Normalized Difference Water Index (NDWI; (green - nir)/(green + nir)) of a Landsat image.

ndwi_fun = function(green, nir){
    (green - nir) / (green + nir)
}

ndwi_rast = lapp(landsat[[c(2, 4)]], fun = ndwi_fun)
plot(ndwi_rast)

calculate a correlation between NDVI and NDWI for this area

ndvi_fun = function(nir, red){
  (nir - red) / (nir + red)
}

ndvi_rast = lapp(landsat[[c(4, 3)]], fun = ndvi_fun)
combine = c(ndvi_rast, ndwi_rast)
plot(combine)

layerCor(combine, fun = cor)
##            lyr1       lyr1
## lyr1  1.0000000 -0.9132838
## lyr1 -0.9132838  1.0000000

Use terra::distance() to compute distances from all cells of spain to it’s nearest coastline. According to the documentation, terra::distance() will calculate distance for all cells that are NA to the nearest cell that are not NA

spain_dem = aggregate(spain_dem, fact = 20)
plot(spain_dem)

water_mask = is.na(spain_dem)
water_mask[water_mask == 0] = NA
plot(water_mask)

distance_to_coast = terra::distance(water_mask)
## 
|---------|---------|---------|---------|
=========================================
                                          
distance_to_coast_km = distance_to_coast / 1000

plot(distance_to_coast_km, main = "Distance to the coast (km)")

Try to modify the approach used in the above exercise by weighting the distance raster with the elevation raster; every 100 altitudinal meters should increase the distance to the coast by 10 km. Next, compute and visualize the difference between the raster created using the Euclidean distance (E7) and the raster weighted by elevation.

distance_to_coast_km2 = distance_to_coast_km + ((spain_dem / 100) * 10)

plot(distance_to_coast_km2)

Exercise 3:: Geometry Operations with Rasters

The srtm raster has a resolution of 0.00083 by 0.00083 degrees. Change its resolution to 0.01 by 0.01 degrees using all of the method available in the terra package. Visualize the results.

plot(srtm)

rast_template = rast(ext(srtm), res = 0.01)

srtm_resampl1 = resample(srtm, y = rast_template, method = "bilinear")
srtm_resampl2 = resample(srtm, y = rast_template, method = "near")
srtm_resampl3 = resample(srtm, y = rast_template, method = "cubic")
srtm_resampl4 = resample(srtm, y = rast_template, method = "cubicspline")
srtm_resampl5 = resample(srtm, y = rast_template, method = "lanczos")
srtm_resampl_all = c(srtm_resampl1, srtm_resampl2, srtm_resampl3,
                     srtm_resampl4, srtm_resampl5)
plot(srtm_resampl_all)

Please take this survey to give me feedback on last week and this week!