Nothing
# create a valid URL filter for SoilWeb API
# arguments are NA by default
.buildFilter <- function(series, bbox, mlra, pedlabsampnum, pedon_id, pedon_key) {
# sanity-check: user must supply some kind of criteria
if(missing(series) & missing(bbox) & missing(mlra) & missing(pedlabsampnum) & missing(pedon_id) & missing(pedon_key))
stop('you must provide some filtering criteria')
# init empty filter
f <- vector()
# process filter components
if(!is.na(series)) {
f <- c(f, paste('&series=', series, sep=''))
}
# note: bbox has already been converted into text representation, suitable for URL
if(!is.na(bbox)) {
f <- c(f, paste('&bbox=', bbox, sep=''))
}
if(!is.na(mlra)) {
f <- c(f, paste('&mlra=', mlra, sep=''))
}
if(!is.na(pedlabsampnum)) {
f <- c(f, paste('&pedlabsampnum=', pedlabsampnum, sep=''))
}
if(!is.na(pedon_id)) {
f <- c(f, paste('&pedon_id=', pedon_id, sep=''))
}
if(!is.na(pedon_key)) {
f <- c(f, paste('&pedon_key=', pedon_key, sep=''))
}
# combine filters
f <- paste(f, collapse='')
return(f)
}
# get extended data as JSON -> list, from SoilWeb API
.getExtended_SoilWeb <- function(f) {
# KSSL geochem, XRD, glass
x <- URLencode(paste0('https://casoilresource.lawr.ucdavis.edu/soil_web/kssl/query.php?gzip=1&format=json&what=extended', f))
# list of dataframe objects; note: missing data are returned as FALSE
ext <- .soilDB_curl_get_JSON(x, gzip = TRUE, quiet = TRUE)
# done
return(ext)
}
# get simplified NASIS morphologic data as JSON -> list, from SoilWeb API
.getMorphologic_SoilWeb <- function(f) {
# NASIS morphology
x <- URLencode(paste0('https://casoilresource.lawr.ucdavis.edu/soil_web/kssl/query.php?gzip=1&format=json&what=nasis_morphologic', f))
# list of dataframe objects; note: missing data are returned as FALSE
m <- .soilDB_curl_get_JSON(x, gzip = TRUE, quiet = TRUE)
# done
return(m)
}
# get simplified KSSL site+hz data as JSON -> list, from SoilWeb API
.getKSSL_SoilWeb <- function(f) {
# KSSL site + horizon
x <- URLencode(paste0('https://casoilresource.lawr.ucdavis.edu/soil_web/kssl/query.php?gzip=1&format=json&what=site_hz', f))
# list of dataframe objects; note: missing data are returned as FALSE
site_hz <- .soilDB_curl_get_JSON(x, gzip = TRUE, quiet = TRUE)
# report missing data
if (all(c(isFALSE(site_hz[['site']]),
isFALSE(site_hz[['horizon']])))) {
# no data
return(NULL)
}
# done
return(site_hz)
}
# single record getter function, called by new vectorized interface
.fetchSingle_KSSL <- function(f, returnMorphologicData=FALSE, returnGeochemicalData=FALSE, simplifyColors=FALSE) {
# basic site + hz data
sh <- .getKSSL_SoilWeb(f)
## TODO: this may no longer be a valid approach, once vectorization is complete
# account for missing data
if(is.null(sh)) {
return(NULL)
}
# local copies
s <- sh[['site']]
h <- sh[['horizon']]
# sometimes columns containing all NA are treated as logical, they are not
logical.idx <- which(sapply(h, is.logical))
for(i in logical.idx) { h[[i]] <- as.numeric(h[[i]]) }
# upgrade to SoilProfileCollection
suppressMessages(depths(h) <- pedon_key ~ hzn_top + hzn_bot)
site(h) <- s
# NASIS morphoogy
if(returnMorphologicData) {
m <- .getMorphologic_SoilWeb(f)
}
# get geochemical, optical and XRD data from extended query
if(returnGeochemicalData) {
ext <- .getExtended_SoilWeb(f)
# local copies, these are FALSE if data are missing
geochem <- ext$geochem
optical <- ext$optical
xrd_thermal <- ext$xrd_thermal
# cleanly return missing data as 0-length data.frame
# having a 'labsampnum' (character) column for joins
ul <- 'dummy.labsampnum'
if(isFALSE(geochem)) {
geochem <- data.frame(labsampnum = ul)[0, , drop = FALSE]
}
if(isFALSE(optical)) {
optical <- data.frame(labsampnum = ul)[0, , drop = FALSE]
}
if(isFALSE(xrd_thermal)) {
xrd_thermal <- data.frame(labsampnum = ul)[0, , drop = FALSE]
}
}
## TODO: clean this up
if(returnMorphologicData & returnGeochemicalData)
return(list(SPC=h, morph=m, geochem=geochem, optical=optical, xrd_thermal=xrd_thermal))
else if(returnGeochemicalData)
return(list(SPC=h, geochem=geochem, optical=optical, xrd_thermal=xrd_thermal))
else if(returnMorphologicData)
return(list(SPC=h, morph=m))
else
return(h)
}
# fully vectorized in all arguments except BBOX
#' Get Kellogg Soil Survey Laboratory Data from SoilWeb snapshot
#'
#' Download soil characterization and morphologic data via BBOX, MLRA, or soil
#' series name query, from the KSSL database.
#'
#' This interface has largely been superseded by the Soil Data Access
#'snapshot of the Laboratory Data Mart, available via `fetchLDM()`.
#'
#' Series-queries are case insensitive. Series name is based on the "correlated
#' as" field (from KSSL snapshot) when present. The "sampled as"
#' classification was promoted to "correlated as" if the "correlated as"
#' classification was missing.
#'
#' When \code{returnMorphologicData} is TRUE, the resulting object is a list.
#' The standard output from \code{fetchKSSL} (\code{SoilProfileCollection}
#' object) is stored in the named element "SPC". The additional elements are
#' basic morphologic data: soil color, rock fragment volume, pores, structure,
#' and redoximorphic features. There is a 1:many relationship between the
#' horizon data in "SPC" and the additional dataframes in \code{morph}. See
#' examples for ideas on how to "flatten" these tables.
#'
#' When \code{returnGeochemicalData} is TRUE, the resulting object is a list.
#' The standard output from \code{fetchKSSL} (\code{SoilProfileCollection}
#' object) is stored in the named element "SPC". The additional elements are
#' geochemical and mineralogy analysis tables, specifically:
#' geochemical/elemental analyses "geochem", optical mineralogy "optical", and
#' X-ray diffraction / thermal "xrd_thermal". \code{returnGeochemicalData} will
#' include additional dataframes \code{geochem}, \code{optical}, and
#' \code{xrd_thermal} in list result.
#'
#' Setting \code{simplifyColors=TRUE} will automatically flatten the soil color
#' data and join to horizon level attributes.
#'
#' Function arguments (\code{series}, \code{mlra}, etc.) are fully vectorized
#' except for \code{bbox}.
#'
#' @param series vector of soil series names, case insensitive
#' @param bbox a single bounding box in WGS84 geographic coordinates e.g.
#' \code{c(-120, 37, -122, 38)}
#' @param mlra vector of MLRA IDs, e.g. "18" or "22A"
#' @param pedlabsampnum vector of KSSL pedon lab sample number
#' @param pedon_id vector of user pedon ID
#' @param pedon_key vector of KSSL internal pedon ID
#' @param returnMorphologicData logical, optionally request basic morphologic
#' data, see details section
#' @param returnGeochemicalData logical, optionally request geochemical,
#' optical and XRD/thermal data, see details section
#' @param simplifyColors logical, simplify colors (from morphologic data) and
#' join with horizon data
#' @param progress logical, optionally give progress when iterating over
#' multiple requests
#' @return a \code{SoilProfileCollection} object when
#' \code{returnMorphologicData} is FALSE, otherwise a list.
#' @note SoilWeb maintains a snapshot of these KSSL and NASIS data. The SoilWeb
#' snapshot was developed using methods described here:
#' \url{https://github.com/dylanbeaudette/process-kssl-snapshot}. Please use
#' the link below for the live data.
#' @author D.E. Beaudette and A.G. Brown
#' @seealso \code{\link{fetchOSD}}
#' @references \url{http://ncsslabdatamart.sc.egov.usda.gov/}
#' @keywords utilities
#' @examplesIf curl::has_internet()
#' @examples
#' \donttest{
#' library(aqp)
#'
#' # search by series name
#' s <- fetchKSSL(series='auburn')
#'
#' # search by bounding-box
#' # s <- fetchKSSL(bbox=c(-120, 37, -122, 38))
#'
#' # how many pedons
#' length(s)
#'
#' # plot
#' plotSPC(s, name='hzn_desgn', max.depth=150)
#'
#' ##
#' ## morphologic data
#' ##
#'
#' # get lab and morphologic data
#' s <- fetchKSSL(series='auburn', returnMorphologicData = TRUE)
#'
#' # extract SPC
#' pedons <- s$SPC
#'
#' # if (requireNamespace("farver")) {
#' # ## automatically simplify color data (requires farver)
#' # s <- fetchKSSL(series='auburn', returnMorphologicData = TRUE, simplifyColors=TRUE)
#' # # check
#' # par(mar=c(0,0,0,0))
#' # plot(pedons, color='moist_soil_color', print.id=FALSE)
#' # }
#' }
#'
#' @export fetchKSSL
fetchKSSL <- function(series=NA, bbox=NA, mlra=NA, pedlabsampnum=NA, pedon_id=NA, pedon_key=NA, returnMorphologicData=FALSE, returnGeochemicalData=FALSE, simplifyColors=FALSE, progress=TRUE) {
if(!requireNamespace('jsonlite', quietly=TRUE))
stop('please install the `jsonlite` packages', call.=FALSE)
# convert BBOX into text representation
if(!missing(bbox)) {
# invalid BBOX
if( length(bbox) != 4) {
stop('invalid BBOX')
}
# convert BBOX into text representation
# not vectorized, would require a different kind of input
bbox <- paste(bbox, collapse=',')
}
# create argument matrix
arg <- expand.grid(
series = series,
bbox = bbox,
mlra = mlra,
pedlabsampnum = pedlabsampnum,
pedon_id = pedon_id,
pedon_key = pedon_key
)
if (all(is.na(arg))) {
stop("must specify series, bbox, mlra, pedlabsampnum, pedon_id, or pedon_key argument",
call. = FALSE)
}
# number of unique argument sets (rows in matrix)
n.arg.sets <- nrow(arg)
# list to store results
res <- vector(mode = 'list', length = n.arg.sets)
# disable progress bar when number of argument sets (calls to .fetchSingle_KSSL) is 1
if(n.arg.sets == 1) {
progress <- FALSE
}
# allow toggling of progress bar
if(progress) {
pb <- txtProgressBar(min = 0, max = n.arg.sets, style = 3)
}
# iterate over argument set
for(i in 1:n.arg.sets) {
# build single URL filter
f <- with(
arg[i, ],
.buildFilter(series, bbox, mlra, pedlabsampnum, pedon_id, pedon_key)
)
# process a single request
req.i <- .fetchSingle_KSSL(f, returnMorphologicData, returnGeochemicalData)
# setting a list element to NULL effectively removes it
if(!is.null(req.i)) {
res[[i]] <- req.i
}
if(progress) {
setTxtProgressBar(pb, i)
}
}
if(progress) {
close(pb)
rm(pb)
}
## TODO: enforce unique-ness in results: unique.SPC and unique.data.frame on extended data
## make composite SPC and optionally additional parts
# simple request, result is a list of SPCs
if(!returnMorphologicData & !returnGeochemicalData) {
h <- aqp::pbindlist(res)
# NO site/hz data, stop here
if(is.null(h)) {
message('query returned no data')
return(NULL)
}
} else {
# complex request, result is a list of lists
# SPC
h <- aqp::pbindlist(lapply(res, '[[', 'SPC'))
# NO site/hz data, stop here
if(is.null(h)) {
message('query returned no data')
return(NULL)
}
## TODO enforce unique-ness on SPC here
## NOTE: simpler with purrr::transpose()
# morph
if(returnMorphologicData) {
# add new tables here
v <- c("phcolor", "phfrags", "phpores", "phstructure", "pediagfeatures", "rmf")
m <- vector(mode = 'list')
# iterate over tables and unwind - rbind - store
for(i in v) {
m[[i]] <- do.call('rbind', lapply(res, function(j) j$morph[[i]]))
## TODO enforce unique-ness on data.frames here
}
}
# geochem: may be missing for a lot of labsampnum
if(returnGeochemicalData) {
# safely combine + fill missing columns
geochem <- do.call('rbindlist', args = list(l = lapply(res, '[[', 'geochem'), fill = TRUE))
optical <- do.call('rbindlist', args = list(l = lapply(res, '[[', 'optical'), fill = TRUE))
xrd_thermal <- do.call('rbindlist', args = list(l = lapply(res, '[[', 'xrd_thermal'), fill = TRUE))
# convert back to data.frame
geochem <- as.data.frame(geochem)
optical <- as.data.frame(optical)
xrd_thermal <- as.data.frame(xrd_thermal)
## TODO enforce unique-ness on data.frames here
}
}
# do NOT set KSSL-specific horizon identifier
# labsampnum is NOT unique 0.1% of the time AGB 2019/11/14
# hzidname(h) <- "labsampnum"
# set KSSL-specific hzdesgn/hztexcl fields
hzdesgnname(h) <- "hzn_desgn"
hztexclname(h) <- "lab_texture_class"
## set metadata
# TODO: check before clobbering / consider standard var name
metadata(h)$origin <- 'KSSL via Soilweb / fetchKSSL'
metadata(h)$created <- Sys.time()
# cleaning up the results
if(returnMorphologicData & simplifyColors) {
if(inherits(m$phcolor, 'data.frame')) {
# simplify color data: 1 row / horizon, from morphologic data tables
x.colors <- simplifyColorData(m$phcolor, id.var = 'labsampnum', wt='colorpct')
# safely LEFT JOIN with @horizons
suppressMessages(horizons(h) <- x.colors)
}
}
# report object size
res.size <- round(object.size(res) / 1024 / 1024, 2)
# some feedback via message:
message(paste(length(h), ' pedons loaded (', res.size, ' Mb transferred)', sep=''))
## TODO: clean this up
if(returnMorphologicData & returnGeochemicalData)
return(list(SPC=h, morph=m, geochem=geochem, optical=optical, xrd_thermal=xrd_thermal))
else if(returnGeochemicalData)
return(list(SPC=h, geochem=geochem, optical=optical, xrd_thermal=xrd_thermal))
else if(returnMorphologicData)
return(list(SPC=h, morph=m))
else
return(h)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.