BASICSimulationModel: proposed changes to cryic/frigid logic

Newhall model gains support for additional spatially explicit inputs:
- O horizon (presence/absence)
- Saturation (presence/absence)
- soilAirOffset (degrees C or degrees F)
Cryic Soil Temperature Regime criteria now properly incorporate O horizon and saturation into thresholds

~~DRAFT~~ PR for discussion and testing of changes to cryic/frigid determination in BASICSimulationModel

it appears one of the cryic criteria was misinterpreted at some point, will try to translate it literally and see if results change
as background (more to come on this) climate simulations run in the Sierra Nevada anecdotally skip from mesic directly to cryic. Need to evaluate which of the criteria are triggering this.
May fuel some future proposals for clarifying how to apply cryic criteria 1a and 1b on p.31 in KST (12th edition)

library(jNSMR)
#> jNSMR (0.3.0) -- R interface to the classic Newhall Simulation Model
#> Added JAR file (newhall-1.6.5.jar) to Java class path.
library(soilDB)
library(terra)
library(rgeedim)

x <- vect(fetchSDA_spatial(c("CA067", "CA628", "CA731", "CA729", "CA732"),
                      by.col = "areasymbol",
                      geom.src = "sapolygon"))
#> Using 1 chunks...
#> Chunk #1 completed (n = 5; 4.2 secs)
#> Done in 5 secs; mean/chunk: 4.2 secs; mean/symbol: 1 secs.

# get PRISM climate and SSURGO-derived available water storage data
predictors <- c(newhall_prism_subset(x), newhall_issr800_subset(x))

# take focal window to smooth/remove NoData in AWC
predictors$awc <- focal(predictors$awc, 5, fun="mean", na.rm=TRUE)

# use Lang et al. (2023) canopy height <https://www.nature.com/articles/s41559-023-02206-6>
gd_initialize()

gd_image_from_id('users/nlang/ETH_GlobalCanopyHeight_2020_10m_v1') |>
  gd_download(region = as.polygons(x, ext = TRUE),
              bands = list("b1"),
              filename="canopy_height.tif",
              scale = 100)
#> [1] "canopy_height.tif"
h <- classify(rast("canopy_height.tif"), cbind(NA, 0))

# inspect
plot(h)


# assume if canopy is 2x minimum for tree height (4m) there is an O horizon
predictors$hasOHorizon <- project(h, predictors) > 8

plot(predictors$hasOHorizon)


# add elevation
gd_image_from_id('USGS/NED') |>
  gd_download(region = as.polygons(x, ext = TRUE),
              bands = list("elevation"),
              filename="elev.tif",
              scale = 100)
#> [1] "elev.tif"
e <- rast("elev.tif")
predictors$elev <- project(e, predictors)

# run batches of models on gridded data
res <- newhall_batch(crop(predictors, ext(predictors)), cores=8, core_thresh = 1e5)
#> Batch 1 of 2 (n=99660 on 8 cores) done in 1.3 mins
#> Batch 2 of 2 (n=61152 on 8 cores) done in 40 secs

y <- fetchSDA_spatial(SDA_spatialQuery(predictors, "areasymbol")$areasymbol,
                      by.col = "areasymbol",
                      geom.src = "sapolygon")
#> Using 6 chunks...
#> Chunk #1 completed (n = 10; 5.3 secs)
#> Chunk #2 completed (n = 10; 9.4 secs)
#> Chunk #3 completed (n = 10; 3 secs)
#> Chunk #4 completed (n = 10; 5.1 secs)
#> Chunk #5 completed (n = 10; 8.1 secs)
#> Chunk #6 completed (n = 4; 0.3 secs)
#> Done in 32.2 secs; mean/chunk: 5.2 secs; mean/symbol: 0.6 secs.
# inspect one modeled soil taxonomic property
plot(
  res$dryDaysAfterSummerSolstice,
  main = "dryDaysAfterSummerSolstice",
  col = grDevices::hcl.colors(50, "cividis")
)
plot(as.lines(vect(y)), add = TRUE, col = "darkorange")


# apply xeric threshold (dry >1.5 months after summer solstice)
plot(
  res$dryDaysAfterSummerSolstice > 45,
  main = "dryDaysAfterSummerSolstice",
  col = grDevices::hcl.colors(2, "cividis")
)
plot(as.lines(vect(y)), add = TRUE, col = "darkorange")


plot(res$temperatureRegime,
     col = grDevices::hcl.colors(7, "cividis"))
plot(as.lines(vect(y)), add = TRUE, col = "darkorange")

ncss-tech / jNSMR

BASICSimulationModel: proposed changes to cryic/frigid logic #6