Package 'aqp'

Description

The aqp (Algorithms for Quantitative Pedology) package for R was developed to address some of the difficulties associated with processing soils information, specifically related to visualization, aggregation, and classification of soil profile data. This package is based on a mix of S3/S4 functions and classes, and most functions use basic dataframes as input, where rows represent soil horizons and columns define properties of those horizons. Common to most functions are the requirements that horizon boundaries are defined as depth from 0, and that profiles are uniquely defined by an id column. The aqp package defines an S4 class, "SoilProfileCollection", for storage of profile-level metadata, as well as summary, print, and plotting methods that have been customized for common tasks related to soils data.

Details

Demos: demo(aqp)

Project homepage

Author(s)

Dylan E. Beaudette [email protected], Pierre Roudier, Andrew G. Brown

See Also

⁠depths<-()⁠, SoilProfileCollection(), sp1, sp2, sp3, sp4, sp5, sp6

Description

You can access the contents of a SoilProfileCollection by profile and horizon "index", i and j, respectively: spc[i, j, ...]. Subset operations are propagated to other slots (such as diagnostics or spatial) when they result in removal of sites from a collection.

• i refers to the profile position within the collection. By default the order is based on the C SORT order of the variable that you specified as your unique profile ID at time of object construction. Note that if your ID variable was numeric, then it has been sorted as a character.

• j refers to the horizon or "slice" index. This index is most useful when either a) working with slice'd SoilProfileCollection or b) working with single-profile collections. j returns the layer in the specified index positions for all profiles in a collection.

• ... is an area to specify an expression that is evaluated in the subset. Currently supported

  • .LAST (last horizon in each profile): return the last horizon from each profile. This uses i but ignores the regular j index.

  • .FIRST (first horizon in each profile): return the last horizon from each profile. This uses i but ignores the regular j index.

  • .HZID (horizon index): return the horizon indices corresponding to i+j+... ("k") constraints

  • .NHZ (number of horizons): return the number of horizons in the profiles resulting from i+j+... ("k") constraints

Usage

## S4 method for signature 'SoilProfileCollection'
x[i, j, ..., drop = TRUE]

Arguments

Description

Get the data from a column accessed by name. Column names other than profile ID are not shared between site and horizons. Bonus: [[ gives access to all site and horizon level variables in tab complete for RStudio using the magrittr pipe operator!

Usage

## S4 method for signature 'SoilProfileCollection'
x[[i, j]]

Arguments

Examples

data(sp2)
depths(sp2) <- id ~ top + bottom
site(sp2) <- ~ surface

# get with [[
sp2[['surface']]

# get using "unknown" expression:
#  "2nd + 3rd horizon column names"
for(i in horizonNames(sp2)[2:3])
 print(sp2[[i]])

data(sp5)

# some column names to work with
rgb.columns <- c("R25","G25","B25")

res <- lapply(rgb.columns, function(x) {

  # [[ allows you to access column names in a loop
  round(sp5[[x]] * 255)

})

# rename scaled results
names(res) <- paste0(rgb.columns,"_scl")

# add horizon ID to results
result <- data.frame(hzID = hzID(sp5), do.call('cbind', res))
head(result)

# join result back into horizons
horizons(sp5) <- result

Description

Add or change the data from a column accessed by name. Column names other than profile ID are not shared between site and horizons. The benefit of using double bracket setter over $ is that name can be calculated, whereas with $, it must be known a priori and hard coded.

When using the double bracket setter the length of input and output matching either the number of sites or number of horizons is used to determine which slot new columns are assigned to.

Usage

## S4 replacement method for signature 'SoilProfileCollection'
x[[i]] <- value

Arguments

Description

Get the data from a column accessed by name x$name. Column names other than profile ID are not shared between site and horizons.

Usage

## S4 method for signature 'SoilProfileCollection'
x$name

Arguments

Examples

data(sp1)

depths(sp1) <- id ~ top + bottom

# get data from a column by name (prop)
sp1$prop

Description

Set the data in a column accessed by name spc$name. Column names other than profile ID are not shared between site and horizons.

When using $<-, the length of input and output matching either the number of sites or number of horizons is used to determine which slot new columns are assigned to. Use site(x)$name <- value or horizons(x)$name <- value to be explicit about which slot is being accessed.

Usage

## S4 replacement method for signature 'SoilProfileCollection'
x$name <- value

Arguments

Description

Fix old-style organic horizon depths or depths with a non-standard datum by the "depth accumulation" method.

Usage

accumulateDepths(
  x,
  id = NULL,
  hzdepths = NULL,
  hzname = NULL,
  hzdatum = 0,
  seqnum = NULL,
  pattern = "O",
  fix = TRUE
)

Arguments

Details

The "depth accumulation" method calculates thicknesses of individual horizons and then cumulative sums them after putting them in id + top depth order. The routine tries to determine context based on hzname and pattern. The main transformation is if a top depth is deeper than the bottom depth, the depths are reflected on the Z-axis (made negative). The data are then id + top depth sorted again, the thickness calculated and accumulated to replace the old depths.

This function uses several heuristics to adjust data before transformation and thickness calculation:

Regex matching of horizon designation patterns and similar

• matches of pattern where both top and bottom depth NA -> ⁠[0,1]⁠ ⁠[top,bottom]⁠ depth

• REMOVE horizons that do not match pattern where both top and bottom depths NA

Over-ride hzname handling with the sequence column argument seqnum

• if seqnum column specified "first record with NA hzname" is considered a pattern match if seqnum == 1

Trigger "fixing" with the fix argument:

• Add 1 cm to bottom-most horizons with NA bottom depth

• Add 1 cm thickness to horizons with top and bottom depth equal

• Add 1 cm thickness to horizons with NA top depth and bottom depth 0

Value

A horizon-level data.frame, suitable for promoting to SPC with ⁠depths<-⁠, or a SoilProfileCollection, depending on the class of x.

Examples

# example using hzdatum argument
data(sp4)
depths(sp4) <- id ~ top + bottom
hz <- accumulateDepths(sp4,
                       id = "id",
                       hzdepths = c("top", "bottom"),
                       hzname = "name",
                       hzdatum = 5 * 1:length(sp4))
plot(hz)
  
# example using old-style O horizons
hz <- read.table(text = "peiidref hzdept hzdepb hzname seqnum phiid
                    1        11      0      5      A      2   295
                    2        11      1      0     Oe      1   294
                    3        11      5     13     C1      3   296
                    4        11     13     58     C2      4   297
                    5        11     58    152     C3      5   298
                    6        13      0      5      A      2   303
                    7        13      1      0     Oe      1   302
                    8        13      5     25     Bw      3   304
                    9        13     25     61      C      4   305
                    10       13     61     NA      R      5   306
                    11      136      0     13     A1      3   695
                    12      136      1      0     Oe      2   694
                    13      136      2      1     Oi      1   693
                    14      136     13     61     C1      4   696
                    15      136     61     76     C2      5   697")
                    
depths(hz) <- peiidref ~ hzdept + hzdepb

hz_fixed <- accumulateDepths(hz,
                             id = "peiidref",
                             hzdepths = c("hzdept", "hzdepb"),
                             hzname = "hzname")
is_valid <- checkHzDepthLogic(hz_fixed)$valid

test0 <- subset(hz_fixed, !is_valid)
test1 <- subset(hz_fixed, is_valid)

origO <- subset(hz, grepl("O", hzname))
fixedO <- subset(hz_fixed, grepl("O", hzname))

par(mfrow=c(2,1), mar=c(0,0,3,2))
plotSPC(origO, max.depth = 25)
plotSPC(fixedO, max.depth = 25)

Description

Add depth brackets to soil profile sketches.

Usage

addBracket(
  x,
  label.cex = 0.75,
  tick.length = 0.05,
  arrow.length = 0.05,
  offset = -0.3,
  missing.bottom.depth = NULL,
  ...
)

Arguments

Details

x may contain multiple records per profile. Additional examples can be found in this tutorial.

Note

This is a low-level plotting function: you must first plot a SoilProfileCollection object before using this function.

Author(s)

D.E. Beaudette

See Also

addDiagnosticBracket, plotSPC

Examples

# sample data
data(sp1)

# add color vector
sp1$soil_color <- with(sp1, munsell2rgb(hue, value, chroma))

# promote to SoilProfileCollection
depths(sp1) <- id ~ top + bottom

# plot profiles
par(mar = c(0, 0, 0, 1))
plotSPC(sp1, width = 0.3)

# extract min--max depths associated with all A horizons
# result is a single-row data.frame / profile
combinedBracket <- function(i) {
  h <- horizons(i)
  idn <- idname(i)
  this.id <- h[[idn]][1]
  
  idx <- grep('^A', h$name)
  
  res <- data.frame(
    id = this.id,
    top = min(h$top[idx]), 
    bottom = max(h$bottom[idx], na.rm=TRUE)
  )
  names(res)[1] <- idn
  
  return(res)
}

# return matching horizon top / bottom depths for A or C horizons
# result is a 0 or more row data.frame / profile
individualBrackets <- function(i) {
  h <- horizons(i)
  idn <- idname(i)
  this.id <- h[[idn]][1]
  
  idx <- grep('^A|^C', h$name)
  
  res <- data.frame(
    id = this.id,
    top = h$top[idx], 
    bottom = h$bottom[idx]
  )
  names(res)[1] <- idn
  
  return(res)
}

# combined brackets
b1 <- profileApply(sp1, combinedBracket, frameify = TRUE)

# individual brackets
b2 <- profileApply(sp1, individualBrackets, frameify = TRUE)

# plot in reverse order
plotSPC(sp1, plot.order = rev(1:length(sp1)), width = 0.25)

# note that plotting order is derived from the call to `plotSPC(sp1)`
addBracket(b1, col='red', offset = -0.35)

# plot in reverse order
plotSPC(sp1, plot.order = rev(1:length(sp1)), width = 0.25)

# note that plotting order is derived from the call to `plotSPC(sp1)`
addBracket(b2, col='red', offset = -0.35)

Description

Annotate diagnostic features within a sketch of soil profiles.

Usage

addDiagnosticBracket(
  s,
  kind,
  feature = "featkind",
  top = "featdept",
  bottom = "featdepb",
  ...
)

Arguments

Details

Additional examples can be found in this tutorial.

Note

This is a low-level plotting function: you must first plot a SoilProfileCollection object before using this function.

Author(s)

D.E. Beaudette

See Also

addBracket, plotSPC

Description

Symbolize volume fraction on an existing soil profile collection plot.

Usage

addVolumeFraction(
  x,
  colname,
  res = 10,
  cex.min = 0.1,
  cex.max = 0.5,
  pch = 1,
  col = "black"
)

Arguments

Details

This function can only be called after plotting a SoilProfileCollection object. Details associated with a call to plotSPC() are automatically accounted for within this function: e.g. plot.order, width, etc..

Note

It may be necessary to adjust both res, cex.min, and cex.max for optimal legibility.

Author(s)

D.E. Beaudette

See Also

plotSPC()

Description

Summarize soil color data, weighted by occurrence and horizon thickness.

Usage

aggregateColor(
  x,
  groups = "genhz",
  col = "soil_color",
  colorSpace = "CIE2000",
  k = NULL,
  profile_wt = NULL,
  mixingMethod = c("estimate", "exact")
)

Arguments

Details

Weights are computed by: w_i = sqrt(sum(thickness_i)) * n_i where w_i is the weight associated with color i, thickness_i is the total thickness of all horizons associated with the color i, and n_i is the number of horizons associated with color i. Weights are computed within groups specified by groups.

Value

A list with the following components:

Author(s)

D.E. Beaudette

See Also

generalize.hz()

Examples

# keep examples from using more than 2 cores
data.table::setDTthreads(Sys.getenv("OMP_THREAD_LIMIT", unset = 2))

# load some example data
data(sp1, package='aqp')

# upgrade to SoilProfileCollection and convert Munsell colors
sp1$soil_color <- with(sp1, munsell2rgb(hue, value, chroma))
depths(sp1) <- id ~ top + bottom
site(sp1) <- ~ group

# generalize horizon names
n <- c('O', 'A', 'B', 'C')
p <- c('O', 'A', 'B', 'C')
sp1$genhz <- generalize.hz(sp1$name, n, p)

# aggregate colors over horizon-level attribute: 'genhz'
a <- aggregateColor(sp1, groups = 'genhz', col = 'soil_color')

# check results
str(a)

## Not run: 
# aggregate colors over site-level attribute: 'group'
a <- aggregateColor(sp1, groups = 'group', col = 'soil_color')

# aggregate colors over site-level attribute: 'group'
# discretize colors to 4 per group
a <- aggregateColor(sp1, groups = 'group', col = 'soil_color', k = 4)

# aggregate colors over depth-slices
s <- dice(sp1, c(5, 10, 15, 25, 50, 100, 150) ~ soil_color)
s$slice <- paste0(s$top, ' cm')
s$slice <- factor(s$slice, levels=guessGenHzLevels(s, 'slice')$levels)
a <- aggregateColor(s, groups = 'slice', col = 'soil_color')

  # optionally plot with helper function
  if(require(sharpshootR))
    aggregateColorPlot(a)

# a more interesting example
  data(loafercreek, package = 'soilDB')
  
  # generalize horizon names using REGEX rules
  n <- c('Oi', 'A', 'BA','Bt1','Bt2','Bt3','Cr','R')
  p <- c('O', '^A$|Ad|Ap|AB','BA$|Bw', 
         'Bt1$|^B$','^Bt$|^Bt2$','^Bt3|^Bt4|CBt$|BCt$|2Bt|2CB$|^C$','Cr','R')
  loafercreek$genhz <- generalize.hz(loafercreek$hzname, n, p)
  
  # remove non-matching generalized horizon names
  loafercreek$genhz[loafercreek$genhz == 'not-used'] <- NA
  loafercreek$genhz <- factor(loafercreek$genhz)
  
  a <- aggregateColor(loafercreek, 'genhz')
  
  # plot results with helper function
  par(mar=c(1,4,4,1))
  aggregateColorPlot(a, print.n.hz = TRUE)
  
  # inspect aggregate data
  a$aggregate.data

## End(Not run)

Description

Estimate the most-likely depth to contact within a collection of soil profiles. Consider getSoilDepthClass followed by group-wise percentile estimation as a faster alternative.

Usage

aggregateSoilDepth(
  x,
  groups,
  crit.prob = 0.9,
  name = hzdesgnname(x),
  p = "Cr|R|Cd",
  ...
)

Arguments

Details

This function computes a probability-based estimate of soil depth by group. If no grouping variable exists, a dummy value can be used to compute a single estimate. The crit.prob argument sets the critical probability (e.g. 0.9) at which soil depth within a group of profiles is determined. For example, a crit.prob of 0.95 might result in an estimated soil depth (e.g. 120cm) where 95% of the profiles (by group) had depths that were less than or equal to 120cm.

Value

A data.frame is returned, with as many rows as there are unique group labels, as specified in groups.

Author(s)

D.E. Beaudette

See Also

estimateSoilDepth() slab()

Examples

data(sp1)
depths(sp1) <- id ~ top + bottom
site(sp1) <- ~ group

# set horizon designation in SPC
hzdesgnname(sp1) <- 'name'

aggregateSoilDepth(sp1, 'group', crit.prob = 0.9)

Description

This function is used to support relative positioning of soil profiles by plotSPC, based on transect or gradient values typically associated with a site level attribute (e.g. elevation). Gradient values specified in x are translated to the range used by plotSPC (usually ⁠1, length(SPC)⁠) specified in x.min and x.max.

Usage

alignTransect(x, x.min, x.max, fix = TRUE, ...)

Arguments

Details

See the Pair-Wise Distances by Generalized Horizon Labels tutorial for additional examples.

Value

list containing:

• grad: values of x in ascending order

• order: ordering vector of x

• relative.pos: elements of x translated to the new relative scale defined by x.min and x.max

Examples

data("sierraTransect")

# split transects
g <- subset(sierraTransect, transect == 'Granite')
a <- subset(sierraTransect, transect == 'Andesite')

g.p <- alignTransect(g$elev, x.min = 1, x.max = length(g), fix = FALSE)
a.p <- alignTransect(a$elev, x.min = 1, x.max = length(a), fix = FALSE)

op <- par(mar=c(2,0,0,2), mfrow=c(2,1))

plotSPC(g, width=0.25, name.style='center-center', 
        cex.names=0.75, 
        relative.pos = g.p$relative.pos, plot.order = g.p$order)

axis(1, at = g.p$relative.pos, labels = g.p$grad, line = -1.5)

plotSPC(a, width=0.25, name.style='center-center', 
        cex.names=0.75, 
        relative.pos = a.p$relative.pos, plot.order = a.p$order)

axis(1, at = a.p$relative.pos, labels = a.p$grad, line = -1.5)


par(op)

Description

Generic function to allocate soil properties to different classification schemes.

Usage

allocate(
  ...,
  to = c("FAO Salt Severity", "FAO Black Soil", "ST Diagnostic Features"),
  droplevels = FALSE
)

Arguments

Details

This function is intended to allocate a set of soil properties to an established soil classification scheme, such as Salt Severity or Black Soil. Allocation is semantically different from classification. While classification is the 'act' of developing a grouping scheme, allocation is the assignment or identification of measurements to a established class (Powell, 2008).

Usage Details

Each classification scheme (to argument) uses a different set of arguments.

• ⁠FAO Salt Severity⁠

  • EC: electrical conductivity column name, dS/m

  • pH: pH column name, saturated paste extract

  • ESP: exchangeable sodium percentage column name, percent

• ⁠FAO Black Soils⁠

  • object: a data.frame or SoilProfileCollection

  • pedonid: pedon ID column name, required when object is a data.frame

  • hztop: horizon top depth column name, required when object is a data.frame

  • hzbot: horizon bottom depth column name, required when object is a data.frame

  • OC: organic carbon column name, percent

  • m_chroma: moist Munsell chroma column name

  • m_value: moist Munsell value column name

  • d_value: dry Munsell value column name

  • CEC: cation exchange capacity column name (NH4OAc at pH 7), units of cmol(+)/kg soil

  • BS: base saturation column name (NH4OAc at pH 7), percent

  • tropical: logical, data are associated with "tropical soils"

• ⁠ST Diagnostic Features⁠

  • object: a data.frame or SoilProfileCollection

  • pedonid: pedon ID column name, required when object is a data.frame

  • hzname: horizon name column, required when object is a data.frame

  • hztop: horizon top depth column name, required when object is a data.frame

  • hzbot: horizon bottom depth column name, required when object is a data.frame

  • texcl: soil texture class (USDA) column name

  • rupresblkcem: rupture resistance column name

  • m_value: moist Munsell value column name

  • m_chroma: moist Munsell chroma column name

  • d_value: dry Munsell value column name

  • BS: base saturation column name (method ??), percent

  • OC: organic carbon column name, percent

  • n_value: ??

  • featkind: ??

Value

A vector or data.frame object.

Note

The results returned by allocate(to = "ST Diagnostic Features") currently return a limited set of diagnostic features that are easily defined. Also, the logic implemented for some features does not include all the criteria defined in the Keys to Soil Taxonomy.

Author(s)

Stephen Roecker

References

Abrol, I., Yadav, J. & Massoud, F. 1988. Salt-affected soils and their management. No. Bulletin 39. Rome, FAO Soils.

FAO. 2006. Guidelines for soil description. Rome, Food and Agriculture Organization of the United Nations.

FAO. 2020. DEFINITION | What is a black soil? (online). (Cited 28 December 2020). http://www.fao.org/global-soil-partnership/intergovernmental-technical-panel-soils/gsoc17-implementation/internationalnetworkblacksoils/more-on-black-soils/definition-what-is-a-black-soil/es/

Powell, B., 2008. Classifying soil and land, in: McKenzie, N.J., Grundy, M.J., Webster, R., Ringrose-Voase, A.J. (Eds.), Guidelines for Survey Soil and Land Resources, Australian Soil and Land Survey Handbook Series. CSIRO, Melbourne, p. 572.

Richards, L.A. 1954. Diagnosis and Improvement of Saline and Alkali Soils. U. S. Government Printing Office. 166 pp.

Soil Survey Staff, 2014. Keys to Soil Taxonomy, 12th ed. USDA-Natural Resources Conservation Service, Washington, D.C.

Examples

# Salt Severity
test <- expand.grid(
  EC  = sort(sapply(c(0, 0.75, 2, 4, 8, 15, 30), function(x) x + c(0, -0.05, 0.05))),
  pH  = c(8.1, 8.2, 8.3, 8.4, 8.5, 8.6),
  ESP = sort(sapply(c(0, 15, 30, 50, 70, 100), function(x) x + c(0, 0.1, -0.1)))
)
test$ss      <- with(test, allocate(EC = EC, pH = pH, ESP = ESP, to = "FAO Salt Severity"))
table(test$ss)

# Black Soil Category 1 (BS1)
test <- expand.grid(
  dept = seq(0, 50, 10),
  OC   = sort(sapply(c(0, 0.6, 1.2, 20, 40), function(x) x + c(0, -0.05, 0.05))),
  chroma_moist  = 2:4,
  value_moist   = 2:4,
  value_dry     = 4:6,
  thickness     = 24:26,
  CEC           = 24:26,
  BS            = 49:51,
  tropical      = c(TRUE, FALSE)
)
test$pedon_id <- rep(1:21870, each = 6)
test$depb     <- test$dept + 10

bs1 <- allocate(test, pedonid = "pedon_id", hztop = "dept", hzbot = "depb", 
                OC = "OC", m_chroma = "chroma_moist", m_value = "value_moist", 
                d_value = "value_dry", CEC = "CEC", BS = "BS", 
                to = "FAO Black Soil"
)

table(BS1 = bs1$BS1, BS2 = bs1$BS2)


# SoilProfileCollection interface

data(sp3)
depths(sp3) <- id ~ top + bottom
hzdesgnname(sp3) <- 'name'

# fake base saturation
horizons(sp3)$bs <- 75

plotSPC(sp3)

allocate(
  sp3, 
  to = 'FAO Black Soil', 
  OC = 'tc', 
  m_chroma = 'chroma', 
  m_value = 'value', 
  d_value = 'value',
  CEC = 'cec',
  BS = 'bs'
)

# make a copy and edit horizon values
x <- sp3
x$value <- 2
x$chroma <- 2
x$cec <- 26
x$tc <- 2

x$soil_color <- munsell2rgb(x$hue, x$value, x$chroma)

plotSPC(x)

allocate(
  x, 
  to = 'FAO Black Soil', 
  OC = 'tc', 
  m_chroma = 'chroma', 
  m_value = 'value', 
  d_value = 'value',
  CEC = 'cec',
  BS = 'bs'
)


# Soil Taxonomy Diagnostic Features
data(sp1)
sp1$texcl = gsub("gr|grv|cbv", "", sp1$texture)
df <- allocate(object = sp1, pedonid = "id", hzname = "name", 
               hzdept = "top", hzdepb = "bottom", texcl = "texcl", 
               to = "ST Diagnostic Features"
)
aggregate(featdept ~ id, data = df, summary)

Description

This is an accessor and replacement method for the aqp_df_class entry in the metadata slot. This entry is used internally by methods that interact with data.frame objects and slots to ensure that the same class used to promote to the SoilProfileCollection initially is used throughout the process.

Usage

## S4 method for signature 'SoilProfileCollection'
aqp_df_class(object)

## S4 replacement method for signature 'SoilProfileCollection'
aqp_df_class(object) <- value

Arguments

Description

Returns the top depth of the argillic horizon as a numeric vector.

Usage

argillic.clay.increase.depth(p, clay.attr = "clay")

Arguments

Details

Uses crit.clay.argillic to determine threshold clay increase, and get.increase.matrix to determine where increase is met within a vertical distance of 30 cm.

Value

A numeric vector containing top depth of argillic horizon, if present, or NA.

Author(s)

Andrew Gene Brown

See Also

getArgillicBounds, get.increase.matrix, crit.clay.argillic

Examples

data(sp1, package = 'aqp')
depths(sp1) <- id ~ top + bottom
site(sp1) <- ~ group

p <- sp1[1]
attr <- 'prop' # clay contents 
foo <- argillic.clay.increase.depth(p, clay.attr = attr)
foo

Description

SoilProfileCollections can be coerced to other R object types using as(spc, 'type').

Possible endpoints include: list, data.frame, SpatialPointsDataFrame and SpatialPoints.

Usage

## S4 method for signature 'SoilProfileCollection'
as.data.frame(x)

Arguments

Value

list

data.frame

tbl_df

data.table

SpatialPointsDataFrame

sf

SpatialPoints

Examples

# load example data stored as SoilProfileCollection
data(sp5)

# sp5
str(sp5)

# list output
str(as(sp5, 'list'))

# data.frame output
str(as(sp5, 'data.frame'))

# Spatial Objects
# make some random coordinate data for each profile
sp5$x <- sp5$y <- rnorm(length(sp5))
initSpatial(sp5, crs = "OGC:CRS84") <- ~ x + y

if (requireNamespace("sf")) {

  # sf output
  str(as(sp5, 'sf'))

  # SpatialPointsDataFrame output
  str(as(sp5, 'SpatialPointsDataFrame'))

  # SpatialPoints output
  str(as(sp5, 'SpatialPoints'))
  
}

Description

Calculate Redness Index after Barron & Torrent (1986) "Use of the Kubelka—Munk Theory to Study the Influence of Iron Oxides on Soil Colour" using Munsell colors converted to LAB. DOI: 10.1111/j.1365-2389.1986.tb00382.x. Accepts vectorized inputs for hue, value and chroma, produces vector output.

Usage

barron.torrent.redness.LAB(hue, value, chroma)

Arguments

Value

A numeric vector of horizon redness index (higher values = redder).

Author(s)

Andrew G. Brown

References

Barron, V. and Torrent, J. (1986), Use of the Kubelka-Munk theory to study the influence of iron oxides on soil colour. Journal of Soil Science, 37: 499-510. doi:10.1111/j.1365-2389.1986.tb00382.x

Description

Simulate realistic sand/silt/clay values (a composition) using multivariate Normal distribution or Dirichlet distribution. Simulations from the multivariate Normal distribution are based on the compositional mean and variance-covariance matrix. Simulations from the Dirichlet distribution are based on maximum likelihood estimation of alpha parameters.

Usage

bootstrapSoilTexture(ssc, method = c("dirichlet", "normal"), n = 100)

Arguments

Details

Simulations from the multivariate normal distribution will more closely track the marginal distributions of sand, silt, and clay–possibly a better fit for "squished" compositions (TODO elaborate). However, these simulations can result in extreme (unlikely) estimates.

Simulations from the Dirichlet distribution will usually be a better fit (fewer extreme estimates) but require a fairly large number of records in ssc (n >= 30?) for a reliable fit.

Additional examples will be added to this tutorial.

Value

a list containing:

• samples - data.frame of simulated sand, silt, clay values

• mean - compositional mean

• var - compositional variance-covariance matrix

• D.alpha - (fitted) alpha parameters of the Dirichlet distribution, NULL when method = 'normal'

Author(s)

D.E. Beaudette

References

Aitchison, J. (1986) The Statistical Analysis of Compositional Data Monographs on Statistics and Applied Probability. Chapman & Hall Ltd., London (UK). 416p.

Aitchison, J, C. Barcel'o-Vidal, J.J. Egozcue, V. Pawlowsky-Glahn (2002) A concise guide to the algebraic geometric structure of the simplex, the sample space for compositional data analysis, Terra Nostra, Schriften der Alfred Wegener-Stiftung, 03/2003

Malone Brendan, Searle Ross (2021) Updating the Australian digital soil texture mapping (Part 1*): re-calibration of field soil texture class centroids and description of a field soil texture conversion algorithm. Soil Research. https://www.publish.csiro.au/SR/SR20283

Malone Brendan, Searle Ross (2021) Updating the Australian digital soil texture mapping (Part 2*): spatial modelling of merged field and lab measurements. Soil Research. https://doi.org/10.1071/SR20284

Examples

if(
requireNamespace("compositions") &
  requireNamespace("soiltexture")
) {
  
  # sample data, data.frame
  data('sp4')
  
  # filter just Bt horizon data
  ssc <- sp4[grep('^Bt', sp4$name), c('sand', 'silt', 'clay')]
  names(ssc) <- toupper(names(ssc))
  
  # simulate 100 samples
  s <- bootstrapSoilTexture(ssc, n = 100)
  s <- s$samples
  
  # empty soil texture triangle
  TT <- soiltexture::TT.plot(
    class.sys= "USDA-NCSS.TT",
    main= "",
    tri.sum.tst=FALSE,
    cex.lab=0.75,
    cex.axis=0.75,
    frame.bg.col='white',
    class.lab.col='black',
    lwd.axis=1.5,
    arrows.show=TRUE,
    new.mar = c(3, 0, 0, 0)
  )
  
  # add original data points
  soiltexture::TT.points(
    tri.data = s, geo = TT, col='firebrick', 
    pch = 3, cex = 0.5, lwd = 1, 
    tri.sum.tst = FALSE
  )
  
  # add simulated points
  soiltexture::TT.points(
    tri.data = ssc, geo = TT, bg='royalblue', 
    pch = 22, cex = 1, lwd = 1, 
    tri.sum.tst = FALSE
  )
  
  # simple legend
  legend('top', 
         legend = c('Source', 'Simulated'), 
         pch = c(22, 3), 
         col = c('black', 'firebrick'), 
         pt.bg = c('royalblue', NA), 
         horiz = TRUE, bty = 'n'
  )
  
  
}

Description

Compute a multinominal Brier score from predicted class probabilities and observed class label. Lower values are associated with a more accurate classifier.

Usage

brierScore(x, classLabels, actual = "actual")

Arguments

Value

a single Brier score, representative of data in x

Author(s)

D.E. Beaudette

References

Brier, Glenn W. 1950. "Verification of Forecasts Expressed in Terms of Probability." Monthly Weather Review 78 (1): 1-3. doi:10.1175/1520-0493(1950)078<0001:VOFEIT>2.0.CO;2.

Examples

# columns 'a', 'b', 'c' contain predicted probabilities
# column 'actual' contains observed class label

# a good classifier
d.good <- data.frame(
  a = c(0.05, 0.05, 0.10),
  b = c(0.90, 0.85, 0.75),
  c = c(0.05, 0.10, 0.15),
  actual = c('b', 'b', 'b'),
  stringsAsFactors = FALSE
)

# a rather bad classifier
d.bad <- data.frame(
  a = c(0.05, 0.05, 0.10),
  b = c(0.90, 0.85, 0.75),
  c = c(0.05, 0.10, 0.15),
  actual = c('c', 'c', 'c'),
  stringsAsFactors = FALSE
)

# class labels are factors
d.factors <- data.frame(
  a = c(0.05, 0.05, 0.10),
  b = c(0.90, 0.85, 0.75),
  c = c(0.05, 0.10, 0.15),
  actual = c('b', 'b', 'b'),
  stringsAsFactors = TRUE
)

# relatively low value = accurate
brierScore(x = d.good, classLabels = c('a', 'b', 'c'), actual = 'actual')

# high values = not accuate
brierScore(x = d.bad, classLabels = c('a', 'b', 'c'), actual = 'actual')

# message related to conversion of factor -> character
brierScore(x = d.factors, classLabels = c('a', 'b', 'c'), actual = 'actual')

Description

Calculate "Color Development Equivalent" by the method of Buntley & Westin (1965) "A Comparative Study of Developmental Color in a Chestnut-Chernozem-Brunizem Soil Climosequence" DOI: 10.2136/sssaj1965.03615995002900050029x. Originally developed for Mollisols, the Buntley-Westin index has been used as a tool to separate soils based on depth to particular colors.

Usage

buntley.westin.index(hue, chroma)

Arguments

Value

A numeric vector reflecting horizon color development.

Author(s)

Andrew G. Brown

References

Buntley, G.J. and Westin, F.C. (1965), A Comparative Study of Developmental Color in a Chestnut-Chernozem-Brunizem Soil Climosequence. Soil Science Society of America Journal, 29: 579-582. doi:10.2136/sssaj1965.03615995002900050029x

Description

Combine SoilProfileCollection objects or lists of SoilProfileCollection objects. This method provides ... expansion for the pbindlist method.

Usage

## S4 method for signature 'SoilProfileCollection'
c(x, ...)

## S4 method for signature 'SoilProfileCollection'
combine(...)

## S4 method for signature 'list'
combine(...)

Arguments

Value

A SoilProfileCollection

Examples

# example data
spc1 <- random_profile(1, SPC = TRUE)
spc2 <- random_profile(2, SPC = TRUE)
spc3 <- random_profile('A', SPC = TRUE)

# combine into a single SPC, ... interface
spc <- combine(spc1, spc2, spc3)

# combine into a single SPC, list interface
spc <- combine(list(spc1, spc2, spc3))

# input are combined into a single SPC
spc <- c(spc1, spc2, spc3)

# result is a list when a mixture of objects are provided
spc <- c(spc1, bar=spc2, baz="foo")

Description

Site and laboratory data from soils sampled in the central Sierra Nevada Region of California.

Usage

data(ca630)

Format

List containing:

$site : A data frame containing site information.

user_site_id

national user site id

mlra

the MLRA

county

the county

ssa

soil survey area

lon

longitude, WGS84

lat

latitude, WGS84

pedon_key

national soil profile id

user_pedon_id

local soil profile id

cntrl_depth_to_top

control section top depth (cm)

cntrl_depth_to_bot

control section bottom depth (cm)

sampled_taxon_name

soil series name

$lab : A data frame containing horizon information.

pedon_key

national soil profile id

layer_key

national horizon id

layer_sequence

horizon sequence number

hzn_top

horizon top (cm)

hzn_bot

horizon bottom (cm)

hzn_desgn

horizon name

texture_description

USDA soil texture

nh4_sum_bases

sum of bases extracted by ammonium acetate (pH 7)

ex_acid

exchangeable acidity [method ?]

CEC8.2

cation exchange capacity by sum of cations method (pH 8.2)

CEC7

cation exchange capacity by ammonium acetate (pH 7)

bs_8.2

base saturation by sum of cations method (pH 8.2)

bs_7

base saturation by ammonium acetate (pH 7)

Details

These data were extracted from the NSSL database. ca630 is a list composed of site and lab data, each stored as data.frame objects. These data are modeled by a 1:many (site:lab) relation, with the pedon_id acting as the primary key in the site table and as the foreign key in the lab table.

Note

These data are out of date. Pending some new data + documentation. Use with caution

Source

https://ncsslabdatamart.sc.egov.usda.gov/

Examples

## Not run: 
library(tactile)
library(lattice)
library(Hmisc)
library(sp)

# check the data out:
data(ca630)
str(ca630)

# note that pedon_key is the link between the two tables

# make a copy of the horizon data
ca <- ca630$lab

# promote to a SoilProfileCollection class object
depths(ca) <- pedon_key ~ hzn_top + hzn_bot

# add site data, based on pedon_key
site(ca) <- ca630$site

# ID data missing coordinates: '|' is a logical OR
(missing.coords.idx <- which(is.na(ca$lat) | is.na(ca$lon)))

# remove missing coordinates by safely subsetting
if(length(missing.coords.idx) > 0)
	ca <- ca[-missing.coords.idx, ]

# register spatial data
initSpatial(ca) <- ~ lon + lat

# assign a coordinate reference system
prj(ca) <- 'EPSG:4269'

# check the result
print(ca)

# aggregate %BS 7 for all profiles into 1 cm slices
a <- slab(ca, fm= ~ bs_7)

# plot median & IQR by 1 cm slice
xyplot(
top ~ p.q50, 
data = a, 
lower=a$p.q25, 
upper=a$p.q75,
alpha=0.5,
ylim=c(160,-5), 
scales = list(alternating = 1, y = list(tick.num = 7)),
panel = panel.depth_function, 
prepanel = prepanel.depth_function,
ylab='Depth (cm)', xlab='Base Saturation at pH 7',
par.settings = tactile.theme(superpose.line = list(col = 'black', lwd = 2))
)

# aggregate %BS at pH 8.2 for all profiles by MLRA, along 1 cm slices
# note that mlra is stored in @site
a <- slab(ca, mlra ~ bs_8.2)

# keep only MLRA 18 and 22
a <- subset(a, subset=mlra %in% c('18', '22'))

# plot median & IQR by 1 cm slice, using different colors for each MLRA
xyplot(
top ~ p.q50, 
groups = factor(mlra), 
data = a,
lower=a$p.q25, 
upper=a$p.q75,
alpha=0.25,
sync.colors = TRUE,
ylim=c(160,-5), 
scales = list(alternating = 1, y = list(tick.num = 7)),
panel = panel.depth_function, 
prepanel = prepanel.depth_function,
ylab='Depth (cm)', xlab='Base Saturation at pH 7',
par.settings = tactile.theme(superpose.line = list(lwd = 2)),
auto.key = list(lines = TRUE, points = FALSE, columns = 2)
)



# Extract the 2nd horizon from all profiles as SPDF
ca.2 <- ca[, 2]

# subset profiles 1 through 10
ca.1.to.10 <- ca[1:10, ]

# basic plot method: profile plot
par(mar = c(0, 0, 3, 1))
plotSPC(ca.1.to.10, name='hzn_desgn', color = 'CEC7')

## End(Not run)

Description

This function inspects a SoilProfileCollection object, looking for four common errors in horizon depths:

1. bottom depth shallower than top depth

2. equal top and bottom depth

3. missing top or bottom depth (e.g. NA)

4. gap or overlap between adjacent horizons (only if byhz = FALSE)

Usage

checkHzDepthLogic(
  x,
  hzdepths = NULL,
  idname = NULL,
  fast = FALSE,
  byhz = FALSE
)

Arguments

Value

A data.frame containing profile IDs, validity boolean (valid) and test results if fast = FALSE.

The data.frame will have as many rows as profiles in x (length(x)).

• id : Profile IDs, named according to idname(x)

• valid : boolean, profile passes all of the following tests

• depthLogic : boolean, errors related to depth logic

• sameDepth : boolean, errors related to same top/bottom depths

• missingDepth : boolean, NA in top / bottom depths

• overlapOrGap : boolean, gaps or overlap in adjacent horizons (NA when byhz = TRUE)

Author(s)

D.E. Beaudette, A.G. Brown, S.M. Roecker

Examples

## sample data

data(sp3)
depths(sp3) <- id ~ top + bottom

# these data should be clean
res <- checkHzDepthLogic(sp3)

head(res)

# less memory if only concerned about net validity
res <- checkHzDepthLogic(sp3, fast = TRUE)

head(res)

Description

Test for a valid SoilProfileCollection

Usage

checkSPC(x)

Arguments

Details

Test for valid SoilProfileCollection by checking for slots defined in the class prototype. Likely only used between major versions of aqp where internal structure of SoilProfileCollection has changed. Use checkHzDepthLogic to check for common errors in horizon depths.

Value

TRUE or FALSE. Consider using rebuildSPC() if FALSE.

Author(s)

D.E. Beaudette

See Also

rebuildSPC, checkHzDepthLogic

Description

Lookup the n closest Munsell chips from the munsell lookup table from various color notations. This function replaces rgb2munsell().

Usage

col2Munsell(col, space = c("sRGB", "CIELAB"), nClosest = 1)

Arguments

Value

an (NA-padded) data.frame containing hue, value, chroma, and CIE delta-E 2000 color contrast metric between source and nearest matching color(s).

Note

This function is fully vectorized and will pad output with NA-records when NA are present in color.

Author(s)

D.E. Beaudette

References

http://ncss-tech.github.io/AQP/ http://www.brucelindbloom.com/index.html?ColorCalcHelp.html https://www.munsellcolourscienceforpainters.com/MunsellAndKubelkaMunkToolbox/MunsellAndKubelkaMunkToolbox.html http://www.cis.rit.edu/mcsl/online/munsell.php

Examples

# vector of named R colors
col2Munsell(c('red', 'green', 'blue'))

# sRGB matrix in the range of 0-255
col2Munsell(cbind(255, 0, 0))

# sRGB matrix in the range of 0-1
col2Munsell(cbind(1, 0, 0))

# 10YR 5/6 in CIELAB
col2Munsell(
  cbind(51.4337, 9.917916, 38.6889), 
  space = 'CIELAB'
)

# 2.5YR 6/8 in hex notation
col2Munsell("#D18158FF")

# 7.5YR 8/1 in sRGB {0, 1}
col2Munsell(
  cbind(0.8240707, 0.7856834, 0.7541048)
)

# 7.5YR 8/1 in sRGB {0, 255}
col2Munsell(
  cbind(0.8240707, 0.7856834, 0.7541048) * 255
)

# multple colors in CIELAB
col2Munsell(
  parseMunsell(c('10BG 6/6', '2.5YR 4/6'), returnLAB = TRUE),
  space = 'CIELAB'
)

# data.frame input
col2Munsell(
  data.frame(r = 1, g = 0, b = 0),
  space = 'sRGB'
)

# keep examples from using more than 2 cores
data.table::setDTthreads(Sys.getenv("OMP_THREAD_LIMIT", unset = 2))

# Munsell notation to sRGB triplets {0, 1} 
color <- munsell2rgb(
  the_hue = c('10YR', '2.5YR', '5YR'), 
  the_value = c(3, 5, 2.5), 
  the_chroma = c(5, 6, 2), 
  return_triplets = TRUE
)

# result is a data.frame of sRGB {0, 1}
color

# back-transform sRGB -> closest Munsell color
# sigma is the dE00 color contrast metric
col2Munsell(color, space = 'sRGB')

Description

Combines layers and aggregates data by grouping adjacent horizons which match pattern in hzdesgn or, alternately, share a common value in by argument. Numeric properties are combined using the weighted average, and other properties are derived from the dominant condition based on thickness of layers and values in each group.

Usage

collapseHz(
  x,
  pattern = NULL,
  by = NULL,
  hzdesgn = hzdesgnname(x, required = TRUE),
  FUN = function(x, pattern, hzdesgn, ...) grepl(pattern, x[[hzdesgn]], ignore.case =
    FALSE),
  ...,
  AGGFUN = NULL,
  ignore_numerics = NULL,
  na.rm = FALSE
)

Arguments

Details

If a custom matching function (FUN) is used, it should accept arbitrary additional arguments via an ellipsis (...). It is not necessary to do anything with arguments, but the result should match the number of horizons found in the input SoilProfileCollection x.

Custom aggregation functions defined in the AGGFUN argument should either return a single vector value for each group*column combination, or should return a data.frame object with named columns. If the input column name is used as a column name in the result data.frame, then the values of that column name in the result SoilProfileCollection will be replaced by the output of the aggregation function. See examples.

Value

A SoilProfileCollection

Author(s)

Andrew G. Brown

See Also

hz_dissolve()

Examples

data(jacobs2000)

# calculate a new SPC with genhz column based on patterns
new_labels <- c("A", "E", "Bt", "Bh", "C")
patterns <- c("A", "E", "B.*t", "B.*h", "C")
jacobs2000_gen <- generalizeHz(jacobs2000, new = new_labels, pattern = patterns)

# use existing generalized horizon labels
i <- collapseHz(jacobs2000_gen, by = "genhz")

profile_id(i) <- paste0(profile_id(i), "_collapse")

plot(
  c(i, jacobs2000),
  color = "genhz",
  name = "name",
  name.style = "center-center",
  cex.names = 1
)
 
# custom pattern argument  
j <- collapseHz(jacobs2000,
                c(
                  `A` = "^A",
                  `E` = "E",
                  `Bt` = "[ABC]+t",
                  `C` = "^C",
                  `foo` = "bar"
                ))
profile_id(j) <- paste0(profile_id(j), "_collapse")
plot(c(j, jacobs2000), color = "clay")

# custom aggregation function for matrix_color_munsell
k <- collapseHz(jacobs2000,
                pattern = c(
                  `A` = "^A",
                  `E` = "E",
                  `Bt` = "[ABC]+t",
                  `C` = "^C",
                  `foo` = "bar"
                ),
                AGGFUN = list(
                  matrix_color_munsell = function(x, top, bottom) {
                    thk <- bottom - top
                    if (length(x) > 1) {
                      xord <- order(thk, decreasing = TRUE)
                      paste0(paste0(x[xord], " (t=", thk[xord], ")"), collapse = ", ")
                    } else
                      x
                  }
                )
              )
profile_id(k) <- paste0(profile_id(k), "_collapse_custom")

unique(k$matrix_color_munsell)

# custom aggregation function for matrix_color_munsell (returns data.frame)
m <- collapseHz(jacobs2000,
                pattern = c(
                  `A` = "^A",
                  `E` = "E",
                  `Bt` = "[ABC]+t",
                  `C` = "^C",
                  `foo` = "bar"
                ),
                AGGFUN = list(
                  matrix_color_munsell = function(x, top, bottom) {
                    thk <- bottom - top
                    if (length(x) > 1) {
                      xord <- order(thk, decreasing = TRUE)
                      data.frame(matrix_color_munsell = paste0(x, collapse = ";"),
                                 n_matrix_color = length(x))
                    } else {
                      data.frame(matrix_color_munsell = x, 
                                 n_matrix_color = length(x))
                    }
                  }
                )
              )
profile_id(m) <- paste0(profile_id(m), "_collapse_custom")

m$matrix_color_munsell.n_matrix_color

Description

Visualize soil colors in Munsell notation according to within-group frequency.

Usage

colorChart(
  m,
  g = factor("All"),
  size = TRUE,
  annotate = FALSE,
  chip.cex = 3,
  chip.cex.min = 0.1,
  chip.cex.max = 1.5,
  chip.border.col = "black",
  annotate.cex = chip.cex * 0.25,
  annotate.type = c("count", "percentage"),
  threshold = NULL
)

Arguments

Value

a trellis object

Examples

# required for latticeExtra:useOuterStrips
if(!requireNamespace('latticeExtra')) {
  
  # two hue pages
  ric <- expand.grid(
    hue = c('5YR', '7.5YR'),
    value = 2:8,
    chroma = 2:8
  )
  
  # combine hue, value, chroma into standard Munsell notation
  ric <- sprintf("%s %s/%s", ric$hue, ric$value, ric$chroma)
  
  # note that chip frequency-based size is disabled 
  # because all chips have equal frequency
  colorChart(ric, chip.cex = 4, size = TRUE)
  
  # annotation of frequency
  colorChart(ric, chip.cex = 4, annotate = TRUE)
  
  # bootstrap to larger size
  ric.big <- sample(ric, size = 100, replace = TRUE)
  
  # frequency can be encoded in size
  colorChart(ric.big, chip.cex = 3)
  colorChart(ric.big, chip.cex = 5, annotate = TRUE)
     
  # constant size
  colorChart(ric.big, chip.cex = 3, size = FALSE)
  colorChart(ric.big, chip.cex = 3, size = FALSE, chip.border.col = 'NA')
  
  # simulate colors based dE00 thresholding
  p <- list(
    list(m = '10YR 4/4', thresh = 10, hues = c('10YR', '7.5YR'))
  )
  
  # perform 500 simulations
  s <- simulateColor(method = 'dE00', n = 500, parameters = p)
  
  # result is a list, use the first element
  colorChart(s[[1]], chip.cex = 4)
  
  # increase the possible range of color chip sizes
  colorChart(s[[1]], chip.cex = 4, chip.cex.min = 0.01, chip.cex.max = 2)
  
  # slightly funky support for neutral hues
  N <- sprintf('N %s/', 2:8)
  cols <- c(rep(N, times = 5), ric.big)
  
  # note special panel used to show neutral hues
  colorChart(cols, size = FALSE, annotate = TRUE)
  
  # filter proportions below given threshold
  colorChart(cols, size = FALSE, annotate = TRUE, threshold = 0.01, 
  chip.cex = 4, annotate.type = 'percentage')
  
}

Description

Pair-wise comparisons of Munsell color specifications, based on the NCSS color contrast classes (Soil Survey Technical Note 2) and CIE delta-E 2000 metric.

Usage

colorContrast(m1, m2)

Arguments

Details

This function is fully vectorized but expects input to be of the same length. Use expand.grid() to generate suitable input from 1:many or many:1 type comparisons. See this tutorial for an expanded discussion and more examples. Neutral colors are not mentioned in SSTN2: in this function any comparison to a neutral color (e.g. 'N 3/') are assigned a delta-hue of 1. Since SSTN2 expects hues to be counted clock wise from 5R, it possible to get very large delta-hue values for otherwise adjacent colors: '5R' vs. '2.5R'. This will be addressed in an update to the standards.

The most meaningful representation of color contrast is the CIE2000 (dE00) metric.

Value

data.frame with the following columns:

• m1: Munsell color 1

• m2: Munsell color 2

• dH: delta-hue, as computed by huePosition

• dV: delta-value, absolute value of difference in Munsell value (m1 vs. m2)

• dc: delta-chroma, absolute value of difference in Munsell chroma (m1 vs. m2)

• dE00: delta-E00, e.g. the CIE delta-E as refined in 2000

• cc: soil color contrast class, as specified in Soil Survey Technical Note 2.

Note

delta-E00 is computed by the farver package.

Author(s)

D.E. Beaudette

References

1. https://en.wikipedia.org/wiki/Color_difference

See Also

colorContrastPlot, huePosition, huePositionCircle

Examples

# keep examples from using more than 2 cores
data.table::setDTthreads(Sys.getenv("OMP_THREAD_LIMIT", unset = 2))

# two sets of colors to compare
m1 <- c('10YR 6/3', '7.5YR 3/3', '10YR 2/2', '7.5YR 3/4')
m2 <- c('5YR 3/4', '7.5YR 4/4', '2.5YR 2/2', '7.5YR 6/3')

# contrast metrics
colorContrast(m1, m2)

# adjacent chips
colorContrast('10YR 3/3', '10YR 3/4')
colorContrast('10YR 3/3', '7.5YR 3/3')

# highly contrasting colors
# http://colour.granjow.net/fabercastell-polychromos.html
colorContrastPlot('10B 4/13', '10YR 10/15',
labels = c('helioblue-reddish', 'light cadmium yellow')
)


## Note: neutral hues aren't defined in TN2
# approximation / extension of the concept
colorContrast(m1 = 'N 3/', m2 = 'N 6/')

#
colorContrast(m1 = '10YR 3/3', m2 = 'N 3/')

m1 <- c('10YR 6/3', '7.5YR 3/3', '10YR 2/2', 'N 3/')
m2 <- c('5YR 3/4', '7.5YR 4/4', '2.5YR 2/2', '7.5YR 6/3')
colorContrast(m1, m2)

Description

A simple display of two sets of colors, NCSS color contrast class and CIE delta-E00.

Usage

colorContrastPlot(
  m1,
  m2,
  col.cex = 1,
  col.font = 2,
  d.cex = 1,
  cc.font = 3,
  dE00.font = 1,
  labels = c("m1", "m2"),
  label.cex = 1,
  label.font = 1,
  printMetrics = TRUE,
  ...
)

Arguments

Note

This function requires the farver package for calculation of CIE delta-E00.

Author(s)

D.E. Beaudette

See Also

colorContrast()

Examples

# keep examples from using more than 2 cores
data.table::setDTthreads(Sys.getenv("OMP_THREAD_LIMIT", unset = 2))

# two sets of colors to compare
m1 <- c('10YR 6/3', '7.5YR 3/3', '10YR 2/2', '7.5YR 3/4')
m2 <- c('5YR 3/4', '7.5YR 4/4', '2.5YR 2/2', '7.5YR 6/3')

# contrast metrics
colorContrast(m1, m2)

# graphical display
colorContrastPlot(m1, m2)

Description

Estimate central tendency and spread of soil color using marginal quantiles and L1 median of CIELAB coordinates.

Usage

colorQuantiles(soilColors, p = c(0.05, 0.5, 0.95))

Arguments

Details

Colors are converted from sRGB to CIELAB (D65 illuminant), marginal quantiles of (L,A,B) coordinates are estimated, and L1 median (L,A,B) is estimates. The closest Munsell chips (via Munsell/CIELAB lookup table provided by munsell) and R colors are determined by locating chips closest to the marginal quantiles and L1 median.

The results can be conveniently inspected using plotColorQuantiles().

Value

A List containing the following elements:

• marginal: data.frame containing marginal quantiles in CIELAB (D65), closest Munsell chips, and dE00

• L1: L1 median CIELAB (D65) values, closest Munsell chip, and dE00

Author(s)

D.E. Beaudette

Examples

## Not run: 
# example data, see manual page for details
data(sp5)

# slice top 25 cm
# 24-25cm is the last slice
s <- dice(sp5, 0:24 ~ .)

# check some of the data
par(mar=c(0,0,0,0))
plotSPC(sample(s, 25), divide.hz = FALSE, name = '', print.id = FALSE, width = 0.5)

# colors
previewColors(unique(s$soil_color))

# compute marginal quantiles and L1 median
cq <- colorQuantiles(s$soil_color)

# simple graphical display of results
plotColorQuantiles(cq)

## End(Not run)

Description

Compare site level attributes of a SoilProfileCollection object, returning a distance matrix conformal with the output from NCSP(). Values are within the range of 0-1.

Usage

compareSites(x, vars, weights = rep(1, times = length(vars)), ...)

Arguments

Details

This function is typically used in conjunction with the output from NCSP().

Value

dissimilarity / dist class object containing pair-wise distances, row/column names derived from profile_id(x)

See Also

NCSP() cluster::daisy()

Description

compositeSPC() is a convenience function that returns a named list representation of the columns from the @site and @horizons slots.

Usage

compositeSPC(object)

Arguments

Value

A list.

Author(s)

Andrew G. Brown.

Description

Calculate the confusion index of Burrough et al., 1997.

Usage

confusionIndex(x)

Arguments

Value

A single numeric value.

Author(s)

D.E. Beaudette

References

Burrough, P.A., P.F.M. van Gaans, and R. Hootsmans. 1997. "Continuous Classification in Soil Survey: Spatial Correlation, Confusion and Boundaries." Geoderma 77: 115-35. doi:10.1016/S0016-7061(97)00018-9.

Examples

# a very simple example
p <- c(0.25, 0.25, 0.4, 0.05, 0.05)
confusionIndex(p)

# for comparison
shannonEntropy(p)

Description

Compare one or more pages from a simulated Munsell book of soil colors to a reference color.

Usage

contrastChart(
  m,
  hues,
  ccAbbreviate = 1,
  style = "hue",
  gridLines = FALSE,
  de00.cex = 0.6,
  cc.cex = 0.6,
  thresh = NULL,
  returnData = FALSE
)

Arguments

Details

A simulated Munsell color book page or pages are used to demonstrate color contrast between all chips and the reference color m (highlighted in red). NCSS color contrast class and CIE delta-E00 values are printed below all other color chips. Munsell color chips for chroma 5 and 7 are omitted, but axis labels are retained as a reminder of this fact.

Setting style='hue' emphasizes the contrast classes and CIE delta-E00 of chips adjacent to m. Setting style='CC' emphasizes adjacent chips according to respective contrast class via lattice panels.

Two-way panels are used when multiple hues are provided and style='CC'. The default output can be greatly enhanced via:

latticeExtra::useOuterStrips(..., strip = strip.custom(bg=grey(0.85)), strip.left = strip.custom(bg=grey(0.85)) )

Author(s)

D.E. Beaudette

Examples

# keep examples from using more than 2 cores
data.table::setDTthreads(Sys.getenv("OMP_THREAD_LIMIT", unset = 2))

# single hue page
contrastChart(m = '10YR 3/3', hues = '10YR')

# multiple hue pages
contrastChart(m = '10YR 3/3', hues = c('10YR', '2.5Y'))

# contrast class, single hue
contrastChart(m = '10YR 3/3', hues = '10YR', style='CC')

# contrast class, multiple hues
# consider latticeExtra::useOuterStrips()
contrastChart(m = '10YR 5/6', hues = c('10YR', '2.5Y'), style='CC')

Description

Determine soil color contrast class according to methods outlined in the Soil Survey Manual. This function is typically called from colorContrast() which is simpler to use and provides more information.

Usage

contrastClass(v1, c1, v2, c2, dH, dV, dC, verbose = FALSE)

Arguments

Details

This function is fully vectorized but expects all inputs have the same length.

Value

A vector of color contrast classes (ordered factor). A list when verbose is TRUE.

Author(s)

D.E. Beaudette

References

• Soil Survey Technical Note 2 wayback machine URL

See Also

colorContrast

Examples

## standard use, result is an ordered factor
# 10YR 6/3 vs 5YR 3/4
contrastClass(v1=6, c1=3, v2=3, c2=4, dH=2, dV=3, dC=1)

## verbose output, useful for testing rules/cases
# 10YR 6/3 vs 5YR 3/4
contrastClass(v1=6, c1=3, v2=3, c2=4, dH=2, dV=3, dC=1, verbose = TRUE)

Description

Apply rock fragment or salt correction to available water content

Usage

correctAWC(
  awc,
  total_rf = numeric(length(awc)),
  gravel = NULL,
  ec = NULL,
  nullFragsAreZero = TRUE
)

Arguments

Value

A numeric vector (double) containing estimated available water capacities corrected for rock fragments and salts

Examples

# medium organic matter, loam texture 
base.awc <- 0.18 # estimateAWC(texcl = "l", omcl = 2, na.rm = TRUE)

# medium organic matter, loam texture w/ 23% rock fragments by volume 
corrected.awc <- correctAWC(base.awc, total_rf = 23)
corrected.awc

# medium organic matter, loam texture w/ 0% frags by volume and 8 mmhos/cm salts
salty.awc <- correctAWC(base.awc, total_rf = 0, ec = 8)
salty.awc

Description

Given a vector or matrix of "eluvial" horizon clay contents (\ crit.clay.argillic() returns a vector or matrix of minimum clay contents (thresholds) that must be met for an argillic horizon clay increase.

Usage

crit.clay.argillic(eluvial_clay_content)

Arguments

Details

Uses the standard equations for clay contents less than 15 \ and 40 \ the definition of the argillic horizon from 12th Edition Keys to Soil Taxonomy (Soil Survey Staff, 2014).

Value

A vector or matrix (input-dependent) containing minimum "illuvial" horizon clay contents (thresholds) to be met for argillic horizon clay increase.

Note

This function is intended for identifying clay content threshold required for an argillic horizon. These thresholds may not apply depending on the specifics of your soil. E.g. if the upper part of argillic has been plowed (has Ap immediately over upper boundary) the clay increase requirement can be waived (Soil Survey Staff, 2014).

Author(s)

Andrew Gene Brown

References

Soil Survey Staff. 2014. Keys to Soil Taxonomy, 12th ed. USDA-Natural Resources Conservation Service, Washington, DC.

See Also

getArgillicBounds, get.increase.matrix

Examples

# crit.clay.argillic uses different equations for clay content
# less than 15 %, between 15 and 40 %, and >40 %

crit.clay.argillic(eluvial_clay_content=c(5, 20, 45))

Description

Create a (redundant) horizon-level attribute from a site-level attribute. Specify a SoilProfileCollection and a site-level attribute from that SPC (by name) to receive a vector of length equal to the number of horizons containing the site-level values. This vector is directly usable with the SoilProfileCollection horizon setter.

denormalize is the inverse operation for the formula interface that "normalizes" a horizon level variable to site level:

site(object) <- ~ horizonvar

Usage

denormalize(object, attr)

Arguments

Details

"Denormalization" is the process of trying to improve the read performance of a database, at the expense of losing some write performance, by adding redundant copies of data or by grouping data. Sometimes it is beneficial to have site-level attributes denormalized for grouping of horizon-level data in analyses. denormalize achieves this result for SoilProfileCollections.

Value

A vector of values of equal length to the number of rows in the horizon table of the input SPC.

Author(s)

Andrew G. Brown, Dylan Beaudette

Examples

data(sp1)

# create a SoilProfileCollection from horizon data
depths(sp1) <- id ~ top + bottom

# create random site-level attribute `sitevar` with a binary (0/1) outcome
sp1$sitevar <- round(runif(length(sp1)))

# use denormalize() to create a mirror of sitevar in the horizon table
# name the attribute something different (e.g. `hz.sitevar`) to 
# prevent collision with the site attribute
# the attributes can have the same name but you will then need 
# site() or horizons() to access explicitly
sp1$hz.sitevar <- denormalize(sp1, 'sitevar')

# compare number of profiles to number of sitevar assignments
length(sp1)
table(sp1$sitevar)

# compare number of horizons to number of horizon-level copies of sitevar `hz.'sitevar`
nrow(sp1)
table(sp1$hz.sitevar)

Description

Get units of depth measurement from metadata. Default value is centimeters.

Usage

## S4 method for signature 'SoilProfileCollection'
depth_units(object)

## S4 replacement method for signature 'SoilProfileCollection'
depth_units(object) <- value

Arguments

Examples

data(sp5)

## get depth units
du <- depth_units(sp5)

# set alternate units; e.g. inches
depth_units(sp5) <- 'in'

# replace original value (cm)
depth_units(sp5) <- du

Description

The depthOf family of functions calculate depth of occurrence of a horizon designation pattern, or any other value that can be coerced to character and matched with a regular expression.

If you need all depths of occurrence for a particular pattern, depthOf is what you are looking for. minDepthOf and maxDepthOf are wrappers around depthOf that return the minimum and maximum depth. They are all set up to handle missing values and missing "contacts" with the target pattern.

Usage

depthOf(
  p,
  pattern,
  FUN = NULL,
  top = TRUE,
  hzdesgn = hzdesgnname(p, required = TRUE),
  no.contact.depth = NULL,
  no.contact.assigned = NA_real_,
  na.rm = TRUE,
  simplify = TRUE
)

maxDepthOf(
  p,
  pattern,
  top = TRUE,
  hzdesgn = hzdesgnname(p, required = TRUE),
  no.contact.depth = NULL,
  no.contact.assigned = NA,
  na.rm = TRUE,
  simplify = TRUE
)

minDepthOf(
  p,
  pattern,
  top = TRUE,
  hzdesgn = hzdesgnname(p, required = TRUE),
  no.contact.depth = NULL,
  no.contact.assigned = NA,
  na.rm = TRUE,
  simplify = TRUE
)

Arguments

Value

a numeric vector containing specified depth(s) of horizons matching a pattern. If length(p) > 1 then a data.frame containing profile ID, horizon ID, top or bottom depths, horizon designation and pattern.

Author(s)

Andrew G. Brown

Examples

# construct a fake profile
spc <- data.frame(id=1, taxsubgrp = "Lithic Haploxerepts",
                  hzname   = c("A","AB","Bw","BC","R"),
                  hzdept   = c(0,  20, 32, 42,  49),
                  hzdepb   = c(20, 32, 42, 49, 200),
                  clay     = c(19, 22, 22, 21,  NA),
                  texcl    = c("l","l","l", "l","br"),
                  d_value  = c(5,   5,  5,  6,  NA),
                  m_value  = c(2.5, 3,  3,  4,  NA),
                  m_chroma = c(2,   3,  4,  4,  NA))

# promote to SoilProfileCollection
depths(spc) <- id ~ hzdept + hzdepb
hzdesgnname(spc) <- 'hzname'
hztexclname(spc) <- 'texcl'

# multiple horizons contain B
depthOf(spc, "B")

# deepest top depth of horizon containing B
maxDepthOf(spc, "B")

# shallowest top depth
minDepthOf(spc, "B")

# deepest bottom depth
maxDepthOf(spc, "B", top = FALSE)

# deepest bottom depth above 35cm
maxDepthOf(spc, "B", top = FALSE, no.contact.depth = 35)

# assign infinity (Inf) if B horizon does not start within 10cm
minDepthOf(spc, "B", no.contact.depth = 10, no.contact.assigned = Inf)

Description

⁠depths(<data.frame>) <- <formula>⁠: Initialize SoilProfileCollection

⁠depths(<SoilProfileCollection>)⁠: Extract profile ID and horizon depths from SoilProfileCollection

Usage

## S4 method for signature 'SoilProfileCollection'
depths(x, hzID = FALSE, ...)

## S4 replacement method for signature 'SoilProfileCollection'
depths(object) <- value

## S4 replacement method for signature 'data.frame'
depths(object) <- value

Arguments

Details

The input horizon data, and the resulting profile order, is sorted based on unique profile ID and top depth. ID columns are converted to character, depth columns are converted to integer. If NA values exist in all of the top depths, a prototype with 1 horizon per profile ID is returned, with NA in all non-essential columns. If the input object has 0 rows, a prototype with 0 horizons and 0 rows, but same column names as object, is returned.

Value

a data.frame containing profile ID, top depth, and bottom depth

See Also

horizons() idname() hzidname() horizonDepths()

Examples

# load a SoilProfileCollection
data(jacobs2000, package = "aqp")

depths(jacobs2000)
## init SoilProfileCollection objects from data.frame of horizon data

# load demo data
data(sp1)

# promote to SPC
depths(sp1) <- id ~ top + bottom

# plot
plot(sp1)

# number of profiles
length(sp1)

# number of horizons
nrow(sp1)

Description

depthWeights() calculates the contributing fraction for each pair of horizon top and bottom depths, given an upper and lower boundary.

Usage

depthWeights(top, bottom, upper, lower)

Arguments

Value

A named list.

Author(s)

Andrew G. Brown.

Description

Diagnostic horizons describe features of the soil relevant to taxonomic classification. A single profile may have multiple diagnostic features or horizons, each of which may be comprised of multiple horizons.

• diagnostic_hz() (get method): Get diagnostic feature data from a SoilProfileCollection.

• ⁠diagnostic_hz<-⁠ (set method): Set diagnostic feature data for a SoilProfileCollection. The profile ID column from object (idname(object)) must be present in the replacement value object.

Usage

## S4 method for signature 'SoilProfileCollection'
diagnostic_hz(object)

## S4 replacement method for signature 'SoilProfileCollection'
diagnostic_hz(object) <- value

Arguments

Examples

# load test data
data(sp2)

# promote to SPC
depths(sp2) <- id ~ top + bottom

# assign two profiles a zone related to the mollic epipedon
newdata <- data.frame(id = c("hon-1","hon-17"),
                      featkind = "fixed-depth surface sample",
                      featdept = 0,
                      featdepb = 18)

# do left join
diagnostic_hz(sp2) <- newdata

# inspect site table: newvalue TRUE only for horizons
#  with top depth equal to zero
diagnostic_hz(sp2)

Description

Cut ("dice") soil horizons into 1-unit thick slices. This function replaces aqp::slice(), which will be deprecated in aqp 2.0.

Usage

## S4 method for signature 'SoilProfileCollection'
dice(
  x,
  fm = NULL,
  SPC = TRUE,
  pctMissing = FALSE,
  fill = FALSE,
  strict = TRUE,
  byhz = TRUE,
  verbose = FALSE
)

Arguments

Details

For large and potentially messy collections that may include missing horizon depth logic errors, consider using repairMissingHzDepths() before dice(). Consider using accumulateDepths() before invoking dice() on collections that may contain old-style O horizon notation (e.g. 5-0cm).

Value

a SoilProfileCollection object, or data.frame when SPC = FALSE

Author(s)

D.E. Beaudette and A.G. Brown

See Also

repairMissingHzDepths(), accumulateDepths(), fillHzGaps()

Description

A simple function to duplicate the contents of a SoilProfileCollection object. Old profile IDs are saved as a site-level attribute (oldID) and new IDs are generated using a numeric serial number.

Usage

duplicate(x, times = 3, oldID = ".oldID")

Arguments

Value

a SoilProfileCollection object

Author(s)

D.E. Beaudette

Examples

# sample data
data('sp4')

# promote to SPC
depths(sp4) <- id ~ top + bottom

# duplicate each profile 2 times
d <- duplicate(sp4, times = 2)

# graphical check
par(mar = c(0, 0, 3, 1))
plotSPC(d, color = 'Ca', width = 0.25)

Description

This function attempts to move labels along a 1D coordinate system such that overlap (as specified by threshold) is minimized. An electrostatic simulation applies forces of repulsion between labels that are within thresh (e.g. overlapping) and forces of attraction to a uniformly spaced sequence to iteratively perturb affected labels until either no overlap is reported, or a maximum number of iterations (maxIter) has been reached.

Usage

electroStatics_1D(
  x,
  thresh,
  q = 1,
  chargeDecayRate = 0.01,
  QkA_GrowthRate = 0.05,
  maxIter = 100,
  tiny = 1e-04,
  const = 0.001,
  trace = FALSE,
  ...
)

Arguments

Details

Difficult overlap problems can be addressed by reducing thresh and increasing q. Large values of q can lead to chaotic results.

This function will generate unpredictable output when x contains duplicate values.

This function requires input to be pre-sorted, although interesting "artistic" simulations will often result from unsorted x.

Value

When trace = TRUE a list, otherwise numeric vector with converged attribute.

Author(s)

D.E. Beaudette and K.C. Thompson

See Also

fixOverlap(), SANN_1D()

Examples

# vector of object locations, with potential overlap
x <- c(1, 2, 3, 3.3, 3.8, 5, 6, 7, 8, 9, 10)

# full diagnostic output
z <- electroStatics_1D(x, thresh = 0.65, trace = TRUE, q = 1)
txt <- sprintf("Converged %s (%s iterations)", z$converged, length(z$cost))

plot(
seq_along(z$cost),
z$cost, 
las = 1, 
xlab = 'Iteration', 
ylab = 'Overlap Cost', 
type = 'b', 
main = txt
)

abline(h = 0, lty = 2, col = 2)

# final configuration only
xnew <- electroStatics_1D(x, thresh = 0.65, q = 1)

# check for convergence
attr(xnew, 'converged')

# compare original vs. modified
data.frame(orig = x, new = round(xnew, 2))

Description

A pre-calculated lookup list (made with farver::compare_colour) based on pair-wise color contrast (CIE2000 or dE00) evaluated over all "chips" in the aqp::munsell data set.

The intention is to identify Munsell chips that may be "functionally equivalent" to some other given whole chip elsewhere in the Munsell color space – as discretized in the aqp::munsell lookup table.

"Equivalent" chips are based (fairly arbitrarily) on the 0.001 probability level of dE00 (default Type 7 quantile) within the upper triangle of the 8467x8467 contrast matrix. This corresponds to a dE00 threshold of approximately 2.15.

This is a naive (to the subtleties of human color perception, and overall magnitude of contrast between some of the "chips") but computationally consistent approach. Using the lookup list, as opposed to manual contrast via e.g. farver::compare_colour may have some benefits for efficiency in certain applications where the exact contrast value is not as important as the concept of having some threshold that is non-zero, but very small.

Usage

data(equivalent_munsell)

Format

A named list with 8467 elements, each containing a numeric vector of indices corresponding to the munsell data set, which has 8467 rows (unique, whole-number chips). Names have the format HUE VALUE/CHROMA, e.g. "7.5YR 4/4"

References

Gaurav Sharma, Wencheng Wu, Edul N. Dalal. (2005). The CIEDE2000 Color-Difference Formula: Implementation Notes, Supplementary Test Data, and Mathematical Observations. COLOR research and application. 30(1):21-30. http://www2.ece.rochester.edu/~gsharma/ciede2000/ciede2000noteCRNA.pdf

Thomas Lin Pedersen, Berendea Nicolae and Romain Francois (2020). farver: High Performance Colour Space Manipulation. R package version 2.0.3. https://CRAN.R-project.org/package=farver

Dong, C.E., Webb, J.B., Bottrell, M.C., Saginor, I., Lee, B.D. and Stern, L.A. (2020). Strengths, Limitations, and Recommendations for Instrumental Color Measurement in Forensic Soil Characterization. J Forensic Sci, 65: 438-449. https://doi.org/10.1111/1556-4029.14193

See Also

equivalentMunsellChips

Examples

data(equivalent_munsell)

Description

Uses a pre-calculated lookup list (equivalent_munsell) based on pair-wise CIE2000 contrast (dE00) of LAB color with D65 illuminant for all whole value/chroma "chips" in the aqp::munsell data set.

The intention is to identify Munsell chips that may be "functionally equivalent" to some other given whole value/chroma chip elsewhere in the Munsell color space – as discretized in the aqp::munsell data table. This basic assumption needs to be validated against your end goal: probably by visual inspection of some or all of the resulting sets. See colorContrast and colorContrastPlot.

"Equivalent" chips table are based (fairly arbitrarily) on the 0.001 probability level of dE00 (default Type 7 quantile) within the upper triangle of the 8467x8467 contrast matrix. This corresponds to a dE00 contrast threshold of approximately 2.16.

Usage

equivalentMunsellChips(hue = NULL, value = NULL, chroma = NULL)

Arguments

Value

A named list; Each list element contains a data.frame with one or more rows of "equivalent" Munsell, RGB and LAB color coordinates from munsell data set.

References

Gaurav Sharma, Wencheng Wu, Edul N. Dalal. (2005). The CIEDE2000 Color-Difference Formula: Implementation Notes, Supplementary Test Data, and Mathematical Observations. COLOR research and application. 30(1):21-30. http://www2.ece.rochester.edu/~gsharma/ciede2000/ciede2000noteCRNA.pdf

Thomas Lin Pedersen, Berendea Nicolae and Romain François (2020). farver: High Performance Colour Space Manipulation. R package version 2.0.3. https://CRAN.R-project.org/package=farver

Dong, C.E., Webb, J.B., Bottrell, M.C., Saginor, I., Lee, B.D. and Stern, L.A. (2020). Strengths, Limitations, and Recommendations for Instrumental Color Measurement in Forensic Soil Characterization. J Forensic Sci, 65: 438-449. https://doi.org/10.1111/1556-4029.14193

See Also

colorContrast colorContrastPlot equivalent_munsell

Examples

# 7.5YR 4/4 (the one and only)

equivalentMunsellChips("7.5YR", 4, 4)
#>
#> $`7.5YR 4/4`
#>        hue value chroma         r        g         b        L       A       B
#> 8330 7.5YR     4      4 0.4923909 0.352334 0.2313328 41.26403 10.8689 23.5914

# 7.5YR 1/1 (two chips are equivalent; 3 row result)

equivalentMunsellChips("7.5YR", 1, 1)
#>
#> $`7.5YR 1/1`
#>        hue value chroma         r         g          b        L        A        B
#> 1983  10YR     1      1 0.1345633 0.1087014 0.07606787 10.64787 1.621323 6.847629
#> 6189   5YR     1      1 0.1330994 0.1076359 0.09450179 10.63901 2.489012 3.515146
#> 8303 7.5YR     1      1 0.1329483 0.1082380 0.08862581 10.64210 2.065514 4.623922

# 10YR 6/8 (two chips are equivalent; 3 row result)

equivalentMunsellChips("10YR", 6, 8)
#>
#> $`10YR 6/8`
#>       hue value chroma         r         g         b        L        A        B
#> 2039 10YR     6      7 0.7382230 0.5512957 0.2680260 61.76795 10.50886 44.78574
#> 2040 10YR     6      8 0.7519872 0.5472116 0.2157209 61.77496 11.83215 51.15496
#> 2041 10YR     6      9 0.7642826 0.5433189 0.1559069 61.78085 13.09599 57.49773

# compare visually a very red color

veryred <- equivalentMunsellChips("10R", 6, 28)[[1]]

par(mar=c(0,0,1,1))

pie(rep(1, nrow(veryred)), col = with(veryred, munsell2rgb(hue, value, chroma)),
    label = with(veryred, sprintf("%s %s/%s", hue, value, chroma)))

table(veryred$hue) # 2 hues
#> 
#>  10R 7.5R 
#>    8   17

table(veryred$value) # 2 values
#> 
#>  5  6 
#> 11 14

table(veryred$chroma) # 10 chromas
#> 
#> 21 22 23 24 25 26 27 28 29 30 
#>  1  2  2  3  3  4  3  3  2  2

Description

Estimate available water capacity for fine-earth fraction

Usage

estimateAWC(texcl, omcl, precision = 2, FUN = mean, ...)

Arguments

Value

A numeric vector double containing estimated available water capacities for fine-earth fraction.

Examples

# organic matter, loam texture, low medium and high OM
base.awc <- estimateAWC(c("l","l","l"), c(1, 2, 3), na.rm = TRUE)
base.awc

Description

Estimates the upper and lower boundary of the particle size control section for Mineral or Organic soilsby applying a programmatic version of the particle size control section key from the Keys to Soil Taxonomy (13th edition). See details for assumptions and required data elements.

Usage

estimatePSCS(
  p,
  hzdesgn = hzdesgnname(p, required = TRUE),
  clay.attr = hzmetaname(p, "clay", required = TRUE),
  texcl.attr = hztexclname(p, required = TRUE),
  lieutex = hzmetaname(p, "lieutex"),
  tax_order_field = "tax_order",
  bottom.pattern = "Cr|R|Cd|m",
  simplify = TRUE,
  ...
)

Arguments

Details

Requires information to identify argillic horizons (clay contents, horizon designations) with getArgillicBounds() as well as the presence of plow layers and surface organic soil material. Any getArgillicBounds() arguments may be passed to estimatePSCS.

Also, requires information on taxonomic order to handle Andisols.

In aqp 2.1, particle size control sections of organic soils Histosols and Histels are supported. This requires setting the "lieutex" horizon metadata column using ⁠hzmetaname<-()⁠ Horizon designations containing "f" or "W" are recognized as permafrost and water layers, respectively, for application of the organic soils key to control sections. In lieu textures "SPM" and "PEAT" are used to identify low density organic materials used for surface tier thickness. To avoid using the 160 cm surface tier, set the "lieutex" column to any column that does not contain "SPM" or "PEAT" values.

WARNING: Soils in arenic or grossarenic subgroups, with fragipans, or with strongly contrasting PSCs may not be classified correctly. The author would welcome a dataset to develop this functionality for.

Value

A numeric vector (when simplify=TRUE) containing the top and bottom depth of the particle size control section. First value is top, second value is bottom. If p contains more than one profile, the result is a data.frame with profile ID plus PSCS top and bottom depths.

Author(s)

Andrew Gene Brown

References

Soil Survey Staff. 2014. Keys to Soil Taxonomy, 12th ed. USDA-Natural Resources Conservation Service, Washington, DC.

See Also

getArgillicBounds(), getSurfaceHorizonDepth()

Examples

data(sp1, package = 'aqp')
depths(sp1) <- id ~ top + bottom
site(sp1) <- ~ group

# set required metadata
hzdesgnname(sp1) <- 'name'
hztexclname(sp1) <- 'texture'
hzmetaname(sp1, 'clay') <- 'prop'

x <- estimatePSCS(sp1)
x

# change horizon texture and set inlieu texture column to turn
# first profile into an organic soil
sp1$name[1:6] <- c("Oi1", "Oi2", "Oi3", "Oaf", "Cf1", "Cf2")
sp1$texture <- as.character(sp1$texture)
sp1$texture[1:6] <- c("PEAT", "PEAT", "PEAT", "MUCK", "GRVLS", "GRVLS")
sp1$bottom[6] <- 200
hzmetaname(sp1, 'lieutex') <- 'texture'

y <- estimatePSCS(sp1[1, ], simplify = FALSE)

# 74cm lower boundary is 25cm past the upper boundary of permafrost (49cm)
# but minimum depth is 100cm unless there is a root-limiting layer
y