In ncss-tech/soilReports: R Markdown Reports and Convenience Functions for Soil Survey

knitr::opts_chunk$set(
  echo = F,
  message = FALSE,
  warning = FALSE,
  fig.align = 'center',
  out.width = "100%",
  fig.retina = 2,
  dev = 'png',
  antialias = 'cleartype'
)
# LOAD PACKAGES (used in config.R, shiny.Rmd and report.Rmd)
source("packages.R")

Shiny Pedon Summary

r format(Sys.time(), "%Y-%m-%d")

r input$reportname

ranges are (r p.low.rv.high) percentiles

Taxon Names and Pedon Types

Check to make sure that pedons used within this report have been correctly assigned to this component.

# join missing genhz data to site
f <- peds
site(f) <- s$mgz
site.i <- site(f)

kable(site.i[order(site.i$musym, site.i$pedon_id), c('musym', 'pedon_id', 'taxonname', 'taxonkind', 'pedontype', 'bedrckdepth', 'taxpartsize', 'taxsubgrp', 'missing.genhz')], row.names = FALSE, align = 'c')

MUSYM: r input$s.mu; Component: r input$component_name; Modal: r input$modal_pedon

Grouped Profile Plot

#TODO: abstract these with "reactive" and "static" available as a toggle
peds$taxonname <- factor(peds$taxonname)
groupedProfilePlot(peds, groups = 'taxonname', label = 'pedon_id',
                   print.id = TRUE, id.style = 'side', cex.id = 1.2,
                   cex.names = 1,
                   depth.axis = list(cex = 1.25, line = -3.0),
                   y.offset = 7, group.name.cex = 1,
                   group.name.offset = -8, color = input$thematic_field,
                   width = 0.1, shrink = T, shrink.cutoff = 3)

Map

pedon_locations <- sf::st_as_sf(as(peds,'SpatialPointsDataFrame'))
      slot(mapView(pedon_locations, height=2000, width=2000, label=pedon_locations$pedon_id, map.types=c("Stamen.TopOSMRelief",
                                                "Esri.WorldTopoMap",
                                                "Esri.WorldImagery",
                                                "OpenStreetMap.Mapnik")), 'map') 

Slab-wise Profile plot

Blue line shows the median slab value for the selected set of pedons, with the 5th to 95th percentile envelope shown in gray. Thick red line shows the values from selected modal pedon.

if(!input$thematic_field %in% c("moist_soil_color","dry_soil_color")) {
  sl <- slab(peds, fm = as.formula(paste0(" ~ ",input$thematic_field)))
  a <- xyplot(top ~ p.q50, data=sl, ylab='Depth',
              xlab=paste0(input$thematic_field,'\nmedian bounded by 5th and 95th percentiles'),
              lower=sl$p.q5, upper=sl$p.q95, ylim=c(250,-5),
              panel=panel.depth_function, 
              prepanel=prepanel.depth_function,
              cf=sl$contributing_fraction, aspect = 1.5,
              layout=c(1,1), scales=list(x=list(alternating=1)))

  if(!is.na(input$modal_pedon)) {
    modalped <- horizons(peds[which(input$modal_pedon == paste(site(peds)$peiid,site(peds)$pedon_id,sep=":")),]) 
    #TODO: select modal pedon, need to query additional component data? or use dropdown box of available peds
    modalped$middepth <- modalped$hzdept + (modalped$hzdepb - modalped$hzdept)/2

    foo <- 1:(length(modalped$middepth)*3)
    modalpedplot <- data.frame(foo)
    modalpedplot$y <- 1:length(foo)
    modalpedplot$y[which((foo %% 3)==1)] <- modalped$hzdept
    modalpedplot$y[which((foo %% 3)==2)] <- modalped$middepth
    modalpedplot$y[which((foo %% 3)==0)] <- modalped$hzdepb
    modalpedplot$x <- rep(modalped[[input$thematic_field]],each=3)
    b <- xyplot(y~x,data=modalpedplot, type="l", col="RED",lwd=2,ylim=c(250,-5),layout=c(1,1), scales=list(x=list(alternating=1)),par.settings = list(superpose.line = list(lwd=3)))
    (a + as.layer(b))
  } else {
    a
  }
}

Generalized Grouped Profile Plot

peds$taxonname <- factor(peds$taxonname)
groupedProfilePlot(peds, name = 'genhzraw', groups = 'taxonname', label = 'pedon_id',
                   print.id = TRUE, id.style = 'side', cex.id = 1.2,
                   cex.names = 1,  depth.axis = list(cex = 1.25, line = -3.0),
                   y.offset = 7, group.name.cex = 1,
                   group.name.offset = -8, color = input$thematic_field,
                   width = 0.1, shrink = T, shrink.cutoff = 3)

Generalized Horizon Probability

update(s$ml.hz.plot, aspect = 1.5)

Geomorphology

Hillslope position (2D)

df <- categorical.prop.table(peds$hillslopeprof)
kable(df)

if(!is.na(input$modal_pedon)) {
  modal <- peds[which(input$modal_pedon == paste(site(peds)$peiid,site(peds)$pedon_id,sep=":")),] 
  print(paste0("Modal pedon value: ",modal$hillslopeprof))
}

Geomorphic position - Hills (3D)

df <- categorical.prop.table(peds$geomposhill)
kable(df)

if(!is.na(input$modal_pedon)) {
  modal <- peds[which(input$modal_pedon == paste(site(peds)$peiid,site(peds)$pedon_id,sep=":")),] 
  print(paste0("Modal pedon value: ",modal$geomposhill))
}

Geomorphic position - Mountains (3D)

df <- categorical.prop.table(peds$geomposmntn)
kable(df)

if(!is.na(input$modal_pedon)) {
  modal <- peds[which(input$modal_pedon == paste(site(peds)$peiid,site(peds)$pedon_id,sep=":")),] 
  print(paste0("Modal pedon value: ",modal$geomposmntn))
}

Geomorphic position - Flats (3D)

df <- categorical.prop.table(peds$geomposflats)
kable(df)

if(!is.na(input$modal_pedon)) {
  modal <- peds[which(input$modal_pedon == paste(site(peds)$peiid,site(peds)$pedon_id,sep=":")),] 
  print(paste0("Modal pedon value: ",modal$geomposflats))
}

Geomorphons

df <- categorical.prop.table(peds$gis_geomorphons)
kable(df)

if(!is.na(input$modal_pedon)) {
  modal <- peds[which(input$modal_pedon == paste(site(peds)$peiid,site(peds)$pedon_id,sep=":")),] 
  print(paste0("Modal pedon value: ",modal$gis_geomorphons))
}

Drainage class

df <- categorical.prop.table(peds$drainagecl)
kable(df)

if(!is.na(input$modal_pedon)) {
  modal <- peds[which(input$modal_pedon == paste(site(peds)$peiid,site(peds)$pedon_id,sep=":")),] 
  print(paste0("Modal pedon value: ",modal$drainagecl))
}

Surface Shape (DOWN/ACROSS)

down <- factor(as.character(peds$shapedown),labels =  c("V","L","C"), levels=c("convex","linear","concave"))
acro <- factor(as.character(peds$shapeacross),labels =  c("V","L","C"), levels=c("convex","linear","concave"))
shape <- factor(paste(as.character(down),as.character(acro),sep="/"))
shape[grepl(shape,pattern="NA")] <- NA
shape <- factor(shape)
df <- categorical.prop.table(shape)
kable(df)

if(!is.na(input$modal_pedon)) {
  modal <- peds[which(input$modal_pedon == paste(site(peds)$peiid,site(peds)$pedon_id,sep=":")),] 
  print(paste0("Modal pedon value: ",modal$shapedown,modal$shapeacross))
}

Aspect

if(sum(!is.na(peds$aspect_field)) > 2)
  aspect.plot(peds$aspect_field, q=p.low.rv.high, plot.title=input$pedon_pattern, pch=21, bg='RoyalBlue', col='black', arrow.col=c('grey', 'red', 'grey'))
if(!is.na(input$modal_pedon)) {
  modal <- peds[which(input$modal_pedon == paste(site(peds)$peiid,site(peds)$pedon_id,sep=":")),] 
  print(paste0("Modal pedon value: ",modal$aspect_field))
}

Ecology

Ecological site

df <-  categorical.prop.table(peds$ecositeid)
kable(df)

if(!is.na(input$modal_pedon)) {
  modal <- peds[which(input$modal_pedon == paste(site(peds)$peiid,site(peds)$pedon_id,sep=":")),] 
  print(paste0("Modal pedon value: ",modal$ecositeid))
}

Plant list

TODO: Print component plant data

Horizon

Field-described versus Generalized Horizonation

    if(length(peds)) {
      m <- as.matrix(table(peds$hzname,peds$genhz))
      colz <- colnames(m)
      rowz <- rownames(m)
      x <- length(colz)
      y <- length(rowz)
      buf <- data.frame(field.hz=rowz)
      for(i in 1:x) {
        buf <- cbind(buf, data.frame(m[,i]))
      }
      colnames(buf) <- c("Field HZ", colz)
      rownames(buf) <- NULL
      DT::datatable(buf, options = list(
        pageLength = 250
      ))
    }
#kable(table(peds$genhz,peds$hzname))

Horizonation Network Graph

m <- genhzTableToAdjMat(table(peds$genhz,peds$hzname))
plotSoilRelationGraph(m, graph.mode = 'directed', edge.arrow.size=0.5, vertex.label.family='sans')

Texture

kable(s$tt)
#texture.triangle.low.rv.high(data.frame(sand=peds$sand,silt=peds$silt,clay=peds$clay), p=c(0.05, 0.5, 0.95))  

Color

aggregateColorPlot(aggregateColor(peds, groups = 'genhz', col = 'soil_color'), label.font = 1, label.cex = 0.95, print.n.hz = TRUE)

Morphology

kable(s$rt)

Surface Fragments

kable(s$sf)

Diagnostics

kable(s$dt)

Diagnostics plot

 if(length(peds) & nrow(diagnostic_hz(peds))) {
      #return(diagnostic_hz(peds))
      v.possible <- c("anthropic epipedon","abrupt textural change",
              "andic soil properties","cambic horizon","calcic horizon","durinodes",
              "petrocalcic horizon",
             "fragipan","gypsic horizon","gilgai","glossic horizon","aquic conditions",
             "histic epipedon","sombric horizon","petrogypsic horizon","placic horizon",
             "kandic horizon","lithic contact","secondary carbonates","mollic epipedon",
             "melanic epipedon","natric horizon","ochric epipedon","plaggen epipedon",
             "petroferric contact","permafrost","plinthite","albic horizon",
             "agric horizon","spodic horizon","slickensides","argillic horizon",
             "umbric epipedon","sulfuric horizon","paralithic contact","oxic horizon",
             "salic horizon","duripan","sulfidic materials",
             "interfingering of albic materials","fibric soil materials",
             "hemic soil materials","sapric soil materials","humilluvic material",
             "limnic materials","coprogenous earth","diatomaceous earth","marl",
             "albic materials","mottles with chroma 2 or less","lamellae",
             "ortstein","fragic soil properties","densic contact","densic materials",
             "paralithic materials","anhydrous conditions","cryoturbation",
             "gelic materials","glacic layer","folistic epipedon",
             "redox concentrations","redox depletions with chroma 2 or less",
             "reduced matrix","endosaturation","episaturation","anthric saturation",
             "lithologic discontinuity","strongly contrasting particle size class",
             "gypsum accumulations","salt accumulations","volcanic glass",
             "n value > 0.7","spodic materials","artifacts","fibers",
             "free carbonates","resistant minerals","weatherable minerals",
             "anhydritic horizon","human-altered material","human-transported material",
             "manufactured layer","manufactured layer contact")
      v.use <- gsub(v.possible, pattern = "\\s{1}", replacement=".")
      v <- v.use[v.use %in% siteNames(peds)]
      diagnosticPropertyPlot(peds, v = v, k = 2)
    } else {
      print("No pedons matching criteria or diagnostic features populated.")
    }

Mapunit Summary

print(s$pmg)

---

This document is based on aqp version r utils::packageDescription("aqp", field="Version") and soilDB version r utils::packageDescription("soilDB", field="Version"). Interface based on shiny version r utils::packageDescription("shiny", field="Version") and rmarkdown version r utils::packageDescription("rmarkdown", field="Version").

ncss-tech/soilReports documentation built on Feb. 4, 2024, 11:47 p.m.