This document outlines how to access a new, online, interface to the
ROSETTA model for predicting hydraulic parameters from soil properties.
The ROSETTA API was developed by Dr. Todd Skaggs (USDA-ARS) and links to
the work of Zhang and Schaap, (2017). The ROSETTA function
in soilDB is one possible approach to using the API
from R.
The ROSETTA model relies on a minimum of 3 soil properties, with increasing (expected) accuracy as additional properties are included:
required, sand, silt, clay: USDA soil texture separates (percentages) that sum to 100%
optional, bulk density at 33 kPa (1/3 bar): mass per volume after accounting for >2mm fragments, units of gram/cm3
optional, volumetric water content at 33 kPa (1/3 bar): roughly “field capacity” for most soils, units of cm3/cm3
optional, volumetric water content at 1500 kPa (15 bar): roughly “permanent wilting point” for most plants, units of cm3/cm3
Soil properties must be described, in order, via vars
argument. The API does not use the names but column ordering must
follow: sand, silt, clay, bulk density, volumetric water content at
33kPa (1/3 bar), and volumetric water content at 1500kPa (15 bar).
Consider using the interactive version, with copy/paste functionality.
Version 1: Schaap, M.G., F.J. Leij, and M.Th. van Genuchten. 2001. ROSETTA: a computer program for estimating soil hydraulic parameters with hierarchical pedotransfer functions. Journal of Hydrology 251(3-4): 163-176. doi: 10.1016/S0022-1694(01)00466-8.
Version 2: Schaap, M.G., A. Nemes, and M.T. van Genuchten. 2004. Comparison of Models for Indirect Estimation of Water Retention and Available Water in Surface Soils. Vadose Zone Journal 3(4): 1455-1463. doi: https://doi.org/10.2136/vzj2004.1455.
Version 3: Zhang, Y., and M.G. Schaap. 2017. Weighted recalibration of the Rosetta pedotransfer model with improved estimates of hydraulic parameter distributions and summary statistics (Rosetta3). Journal of Hydrology 547: 39-53. doi: 10.1016/j.jhydrol.2017.01.004.
Install all of the packages that this tutorial depends on via:
# run these commands in the R console
install.packages('httr', dependencies = TRUE)
install.packages('jsonlite', dependencies = TRUE)
install.packages('latticeExtra', dependencies = TRUE)
install.packages('tactile', dependencies = TRUE)
install.packages('soilDB', dependencies = TRUE)
install.packages('aqp', dependencies = TRUE)The following demonstrates how to push soil property data from SDA through the ROSETTA API.
# required libraries
library(soilDB)
# SQL to be submitted to SDA
q <- "SELECT
-- contextual data
compname, comppct_r,
-- horizon morphology
hzname, hzdept_r, hzdepb_r,
-- parameters used by ROSETTA, may contain NULL
sandtotal_r, silttotal_r, claytotal_r, dbthirdbar_r,
wthirdbar_r/100 AS wthirdbar_decimal, wfifteenbar_r/100 AS wfifteenbar_decimal
-- tables of interest
FROM
mapunit AS mu
-- implicit filtering
INNER JOIN component AS co ON mu.mukey = co.mukey
INNER JOIN chorizon AS ch ON co.cokey = ch.cokey
-- single map unit
WHERE mu.mukey = '373596'
-- order for visual inspection
ORDER BY co.cokey, ch.hzdept_r ASC;"
# submit query to SDA
x <- SDA_query(q)
# attempting to use all possible soil properties
vars <- c('sandtotal_r', 'silttotal_r', 'claytotal_r', 'dbthirdbar_r', 'wthirdbar_decimal', 'wfifteenbar_decimal')
# call ROSETTA API
r <- ROSETTA(x, vars = vars)| compname | comppct_r | hzname | hzdept_r | hzdepb_r | sandtotal_r | silttotal_r | claytotal_r | dbthirdbar_r | wthirdbar_decimal | wfifteenbar_decimal | theta_r | theta_s | alpha | npar | ksat | .rosetta.model | .rosetta.version |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Aluf | 50 | H1 | 0 | 117 | 94.4 | 0.6 | 5.0 | 1.60 | 0.075 | 0.042 | 0.035 | 0.367 | -1.244 | 0.273 | 2.705 | 5 | 1 |
| Aluf | 50 | H2 | 117 | 251 | 61.0 | 19.0 | 20.0 | 1.55 | 0.254 | 0.117 | 0.041 | 0.381 | -1.906 | 0.133 | 1.229 | 5 | 1 |
| Hitilo | 30 | H1 | 0 | 137 | 94.9 | 0.6 | 4.5 | 1.55 | 0.077 | 0.040 | 0.032 | 0.381 | -1.228 | 0.258 | 2.736 | 5 | 1 |
| Hitilo | 30 | H2 | 137 | 178 | 55.4 | 17.6 | 27.0 | 1.58 | 0.252 | 0.181 | 0.080 | 0.388 | -1.387 | 0.111 | 1.324 | 5 | 1 |
| Hitilo | 30 | H3 | 178 | 203 | 60.2 | 18.3 | 21.5 | 1.50 | 0.221 | 0.135 | 0.060 | 0.400 | -1.435 | 0.124 | 1.629 | 5 | 1 |
Develop theoretical water retention curves (van Genuchten model) from soil texture class centroids.
# required libraries
library(aqp)
library(latticeExtra)
library(tactile)
# make a table of soil textures
# sand sand/silt/clay values at geometric centroids
tex <- SoilTextureLevels(simplify = FALSE)
x <- texcl_to_ssc(tex)
x <- cbind(x, tex)| sand | silt | clay | tex |
|---|---|---|---|
| 92 | 5 | 3 | cos |
| 92 | 5 | 3 | s |
| 92 | 5 | 3 | fs |
| 92 | 5 | 3 | vfs |
| 82 | 12 | 6 | lcos |
| 82 | 12 | 6 | ls |
Submit example data to the ROSETTA API and generate theoretical water retention curves.
r <- ROSETTA(x, vars = c('sand', 'silt', 'clay'))
# iterate over results and generate VG model curve
res <- lapply(1:nrow(r), function(i) {
# model bounds are given in kPA of suction
vg <- KSSL_VG_model(VG_params = r[i, ], phi_min = 10^-3, phi_max=10^6)
# extract curve and add texture ID
m <- vg$VG_curve
m$texture <- r$tex[i]
return(m)
})
# flatten to data.frame
res <- do.call('rbind', res)| phi | theta | texture |
|---|---|---|
| 0.001 | 0.379 | cos |
| 0.001 | 0.379 | cos |
| 0.002 | 0.379 | cos |
| 0.002 | 0.379 | cos |
| 0.002 | 0.379 | cos |
| 0.003 | 0.379 | cos |
Plot water retention curves. Note that matric potential (typically negative pressure by convention) are displayed as suction (positive values).
tps <- tactile.theme(
plot.line=list(col='royalblue', lwd=2)
)
xyplot(
phi ~ theta | texture, data = res,
type = c('l', 'g'),
scales = list(alternating=3, x=list(tick.number=6), y=list(log=10, tick.number=6)),
yscale.components = yscale.components.logpower,
ylab = expression(Matric~~Potential~~(-kPa)),
xlab = expression(Volumetric~Water~Content~~(cm^3/cm^3)),
par.settings = tps,
strip = strip.custom(bg=grey(0.85)),
as.table = TRUE,
layout = c(7, 3),
main='Idealized Water Retention'
)