# NOT RUN {
library(terra)
res <- mukey.wcs(list(aoi = c(-116.7400, 35.2904, -116.7072, 35.3026), crs = "EPSG:4326"),
db = 'gNATSGO', res = 30)
MUKEY <- unique(values(res))
prp <- setNames(
get_SDA_property(
c("ph1to1h2o_r", "claytotal_r"),
"weighted average",
mukeys = MUKEY,
top_depth = 0,
bottom_depth = 25,
include_minors = TRUE,
miscellaneous_areas = FALSE
)[, c("mukey", "ph1to1h2o_r", "claytotal_r")],
c("ID", "pH1to1_0to25", "clay_0to25")
)
levels(res) <- prp
res2 <- catalyze(res)
res2
plot(res2[['pH1to1_0to25']])
# }
Run the code above in your browser using DataLab