← Back to RTM-Suite

What this manual is

Every function call below actually runs against the real, installed SCOPEinR and ToolsRTM packages — nothing here is illustrative pseudocode. Where the port has real rough edges (a soil energy-balance warning, a non-exported helper, a folder-path quirk), they are shown and explained as they actually occur, not hidden.

library(ToolsRTM)
library(SCOPEinR)
cat("ToolsRTM ", as.character(packageVersion("ToolsRTM")),
    " | SCOPEinR ", as.character(packageVersion("SCOPEinR")), "\n", sep = "")
## ToolsRTM 0.62.5 | SCOPEinR 0.46

SCOPE, and how SCOPEinR relates to ToolsRTM

SCOPE (Soil Canopy Observation, Photochemistry and Energy fluxes) is a coupled radiative-transfer/energy-balance model of a vegetated surface (Van der Tol et al. 2009; Yang et al. 2020, SCOPE 2.0). Given leaf biochemistry, canopy structure, viewing/illumination geometry and meteorology, it simultaneously solves for:

  • Radiative transfer: canopy reflectance and top-of-canopy radiance, in the same sense as PROSAIL/4SAIL.
  • Energy balance: how absorbed radiation partitions into sensible heat, latent heat (transpiration/evaporation) and soil heat flux, for both the canopy and the soil, at leaf/soil temperatures that are solved for rather than prescribed.
  • Photosynthesis / biochemistry: leaf-level CO₂ assimilation (a Farquhar–von Caemmerer–Berry / Collatz-type C3/C4 model) driven by the absorbed light and the solved leaf temperature.
  • Chlorophyll fluorescence: sun-induced fluorescence (SIF) emission, linked mechanistically to the photosynthesis/NPQ (non-photochemical quenching) state of the leaf, and propagated through the canopy to a top-of-canopy radiance.

SCOPEinR is the R port of SCOPE 2.1 (originally MATLAB). It does not reimplement leaf and canopy optics from scratch: the radiative-transfer core (Fluspect leaf optics, 4SAIL canopy geometry) is the same family of models exposed directly in ToolsRTM (see the companion ToolsRTM_PROSAIL_tutorial.Rmd). SCOPEinR calls into that machinery and then adds everything PROSAIL/4SAIL alone does not do: the iterative leaf/soil energy balance (get.ebal()), the biochemical photosynthesis model (get.biochemical()), and the fluorescence radiative transfer (get.RTMf()). ToolsRTM is a hard dependency of SCOPEinR for exactly this reason — you always need both packages loaded to run SCOPE.

Installation and loading

Both packages are distributed from public GitLab repositories and, in this environment, are already installed:

if (!requireNamespace("remotes", quietly = TRUE)) install.packages("remotes")
remotes::install_gitlab("caminoccg/toolsrtm", upgrade = "never")
remotes::install_gitlab("caminoccg/scopeinr", upgrade = "never")

library(SCOPEinR) on its own is not enough — always load ToolsRTM first (or alongside it), since SCOPEinR’s get.SCOPE() dispatches leaf and canopy optics straight into ToolsRTM functions.

Input structure

SCOPE needs two kinds of input beyond the leaf/canopy traits ToolsRTM already required for PROSAIL: simulation options (which sub-models to switch on) and a parameter table with the added biochemical, meteorological and geometric fields the energy balance needs. Both ship inside the installed package under inst/input/, resolved here with system.file() so the paths work regardless of the current working directory:

path_input <- system.file("input", package = "SCOPEinR")
list.files(path_input)
##  [1] "input_borders.csv"      "input_data_default.csv" "inputs_SCOPE.csv"      
##  [4] "inputs_SCOPE.xlsx"      "leafangles"             "LUT_input.csv"         
##  [7] "mSCOPE.csv"             "radiationdata"          "setoptions.csv"        
## [10] "soil_spectra"

setoptions.csv — simulation options

scope_options <- read.table(file.path(path_input, "setoptions.csv"),
                             header = TRUE, sep = ",")
scope_options[, c("Options", "Value")]
##                 Options Value
## 1                  lite     1
## 2     calc_fluorescence     1
## 3  calc_spectrum_planck     1
## 4   calc_xanthophyllabs     1
## 5          soilspectrum     0
## 6    Fluorescence_model     0
## 7             applTcorr     1
## 8                verify     1
## 9                mSCOPE     0
## 10           simulation     0
## 11     calc_directional     0
## 12   calc_vert_profiles     0
## 13     soil_heat_method     2
## 14         calc_rss_rbs     0
## 15         MoninObukhov     1
## 16                 LIDF     0
## 17           irradiance     0

Each row switches a SCOPE sub-model or numerical behaviour on/off (0/1) or picks a variant (e.g. Fluorescence_model). The ones that matter most for what follows: calc_fluorescence = 1 turns on the Fluspect/RTMf fluorescence chain (without it, data.rad$reflapp and the fluorescence_* outputs are not produced); calc_spectrum_planck = 1 turns on the thermal-emission spectrum; soil_heat_method = 2 picks how the soil heat flux G is estimated (0.35 * Rn when not running a time series); LIDF = 0 means the leaf-inclination distribution is built from LIDFa/LIDFb in the LUT rather than from a measured-angles file.

Parameter tables: LUT_input.csv and inputs_SCOPE.csv

LUT_input.csv is a ready-to-run, single-row example parameter set (one SCOPE simulation) with every column SCOPE expects — leaf biochemistry (N, Cab, Car, …), photosynthetic capacity (Vcmax25, BallBerrySlope, …), canopy structure (LAI, hc, LIDFa, LIDFb, …), soil/meteorology (Ta, Rin, Rli, p, ea, RH, u, …) and viewing/illumination geometry (tts, tto, psi):

LUT_default <- read.table(file.path(path_input, "LUT_input.csv"),
                           header = TRUE, sep = ",")
dim(LUT_default)
## [1]  1 77
LUT_default[1, c("N","Cab","Car","LAI","Vcmax25","tts","tto","Ta","Rin")]
##     N Cab Car LAI Vcmax25 tts tto Ta Rin
## 1 1.5  60  10   3      60  30   0 20 600

inputs_SCOPE.csv is the companion ranges table used to generate new, randomized parameter sets (Section 6): one row per variable, with a lower/upper bound, a sampling Distribution (Uniform, Gaussian or Fixed), and a default value used when the variable is held fixed.

inputLUT <- read.table(file.path(path_input, "inputs_SCOPE.csv"),
                        header = TRUE, sep = ",")
inputLUT[inputLUT$variable %in% c("Cab","LAI","Vcmax25","EWT","LIDFa","tts"),
         c("variable","lower","upper","Distribution","default")]
##    variable lower upper Distribution default
## 2       Cab  5.00  90.0     Gaussian  40.000
## 6       EWT  0.00   0.2      Uniform   0.009
## 15  Vcmax25  0.75 250.0      Uniform  70.000
## 40      LAI  0.10   7.0      Uniform   3.000
## 42    LIDFa -1.00   1.0      Uniform  -0.350
## 74      tts  0.00  15.0      Uniform  30.000

A single simulation with get.SCOPE()

The default LUT_input.csv row plus the default options are enough to run one full simulation. optipar selects the wavelength-dependent leaf optical coefficients for the chosen leaf.model — optipar2021.Pro.CX pairs with leaf.model = "fluspect-CX" (PROSPECT-PRO pigment set, plus the Cx xanthophyll de-epoxidation state for PRI-type simulations).

invisible(capture.output(
  db.sim <- SCOPEinR::get.SCOPE(
    LUT           = LUT_default[1, ],
    options.SCOPE = scope_options,
    optipar       = SCOPEinR::optipar2021.Pro.CX,
    leaf.model    = "fluspect-CX",
    canopy.model  = "fourSAIL",
    get.outputs   = "ALL",
    get.plots     = FALSE
  )
))

get.SCOPE() always returns a list of simulations (one element per LUT row); for a single row that is a list of length 1, and db.sim[[1]] is itself a list with 19 named elements:

length(db.sim)
## [1] 1
length(db.sim[[1]])
## [1] 19
names(db.sim[[1]])
##  [1] "data.spectral"    "data.angles"      "data.rad"         "atmo"            
##  [5] "resist_out"       "data.fluxes"      "data.soil"        "data.leafopt"    
##  [9] "data.leafbio"     "data.canopy"      "data.profiles"    "data.gap"        
## [13] "data.meteo"       "data.thermal"     "data.bcu"         "data.bch"        
## [17] "data.directional" "iter.ebal"        "data.opts"

The elements that matter for everyday use:

Element Contents
data.rad Reflectance/radiance/irradiance spectra (rdd, rdo, rsd, rso, refl, reflapp, Lo_, Esun_, Esky_, fluorescence spectra, absorbed-PAR terms, …) — 87 fields
data.fluxes Scalar energy/carbon fluxes: Rnctot, lEctot, Hctot, Actot, Tcave, and the soil/total equivalents
data.thermal Solved sunlit/shaded canopy and soil temperatures (Tcu, Tch, Tsu, Tsh)
data.bcu / data.bch Biochemical model state for sunlit (bcu) / shaded (bch) leaves: A (assimilation), Ja (electron transport), Kn (NPQ), eta (fluorescence yield factor)
data.canopy, data.leafbio, data.soil, data.meteo, data.angles Echo back the resolved canopy/leaf/soil/meteo/geometry inputs actually used
iter.ebal Energy-balance solver diagnostics: iteration counter, maxit, and the residual maxEBercu/maxEBerch/maxEBers errors (canopy sunlit/shaded/soil)
data.opts, data.spectral, data.gap, data.profiles, data.directional, atmo, resist_out Options actually used, wavelength grid, gap probabilities, vertical profiles, directional (BRDF) output, atmosphere, aerodynamic/surface resistances 
db.sim[[1]]$data.fluxes
## $Rnctot
## [1] 381.7305
## 
## $lEctot
## [1] 131.1323
## 
## $Hctot
## [1] 250.724
## 
## $Actot
## [1] 19.8942
## 
## $Tcave
## [1] 22.04278
## 
## $Rnstot
## [1] 114.9653
## 
## $lEstot
## [1] 50.53012
## 
## $Hstot
## [1] 24.25024
## 
## $Gtot
## [1] 40.23785
## 
## $Tsave
## [1] 26.7952
## 
## $Rntot
## [1] 496.6957
## 
## $lEtot
## [1] 181.6625
## 
## $Htot
## [1] 274.9742
db.sim[[1]]$iter.ebal
## $counter
## [1] 7
## 
## $maxit
## [1] 100
## 
## $maxEBercu
## [1] 0.0365275
## 
## $maxEBerch
## [1] 0.08928034
## 
## $maxEBers
## [1] 0.1082583

iter.ebal$counter is how many iterations the energy-balance solver used to converge; maxEBercu/maxEBerch/maxEBers are its final canopy sunlit/shaded/soil residual errors.

Reflectance components: rdd, rdo, rsd, rso, refl, reflapp

SCOPE separates canopy reflectance into four bidirectional components before summing them, plus the sensor-observed “apparent” reflectance:

  • rdd — diffuse-in / diffuse-out reflectance. Sky light entering the canopy and leaving it also as diffuse radiation; describes multiple scattering inside the canopy.
  • rdo — diffuse-in / directional-out reflectance. Diffuse incoming radiation (sky light) leaving the canopy in the sensor’s specific viewing direction; relevant to BRDF effects under overcast conditions.
  • rsd — direct-in (sun beam) / diffuse-out reflectance. Direct sunlight entering the canopy and emerging as diffuse scattered light; shows sun–canopy interactions and internal scattering.
  • rso — direct-in (sun beam) / directional-out reflectance. Direct sun illumination leaving the canopy exactly toward the sensor; contains hot-spot effects and strong geometry interactions.
  • refl — total canopy reflectance, rdd + rdo + rsd + rso; the “true” physical reflectance of the canopy.
  • reflapp — apparent reflectance, derived from the outgoing radiance at sensor level (accounts for the actual illumination geometry and, when fluorescence is on, differs from refl in the fluorescence emission band because it is radiance-derived). This is what a real sensor (e.g. Sentinel-2, PRISMA) would observe.
bands <- SCOPEinR::define.bands()
wlS <- bands$wlS[1:2001]   # 400-2500 nm, shortwave-only part of data.rad's 2162-length vectors

rad <- db.sim[[1]]$data.rad
refl_components <- data.frame(
  Wavelength = rep(wlS, 6),
  Value = c(rad$rdd[1:2001], rad$rdo[1:2001], rad$rsd[1:2001],
            rad$rso[1:2001], rad$refl[1:2001], rad$reflapp[1:2001]),
  Component = rep(c("rdd","rdo","rsd","rso","refl","reflapp"), each = 2001)
)

library(ggplot2)
ggplot(refl_components, aes(x = Wavelength, y = Value, color = Component)) +
  geom_line(linewidth = 0.5) +
  theme_bw() + xlim(400, 2500) +
  xlab("Wavelength (nm)") + ylab("Reflectance") +
  ggtitle("Reflectance components for the default LUT_input.csv simulation")

reflapp tracks refl closely outside the fluorescence band but departs from it around 640-850nm — the signature of SCOPE’s fluorescence emission riding on top of the reflected signal, exactly what a real sensor would pick up in that region.

Building a randomized LUT with getLUT.SCOPE()

getLUT.SCOPE(inputLUT, nLUT, setseed) samples nLUT parameter sets from the ranges table read above: Uniform rows are sampled with stats::runif(), Fixed rows are repeated at their default, and anything else (Gaussian) is sampled with ToolsRTM::gauss_byMin_Max() (truncated Gaussian within [lower, upper], centered at Mean_D with spread Std_D).

n_samples <- 10
LUT <- getLUT.SCOPE(inputLUT = inputLUT, nLUT = n_samples, setseed = 42)
dim(LUT)
## [1] 10 76
summary(LUT[, c("Cab","LAI","Vcmax25")])
##       Cab             LAI           Vcmax25      
##  Min.   :22.22   Min.   :1.020   Min.   : 29.71  
##  1st Qu.:44.34   1st Qu.:1.938   1st Qu.: 47.81  
##  Median :47.61   Median :3.754   Median : 96.58  
##  Mean   :51.67   Mean   :3.726   Mean   : 97.08  
##  3rd Qu.:60.26   3rd Qu.:5.494   3rd Qu.:137.01  
##  Max.   :80.23   Max.   :6.434   Max.   :182.64

Real traits co-vary. ToolsRTM::getCor() (covered in depth in the ToolsRTM manual) generates correlated pairs instead of independent draws — this is exactly the pattern used in Carlos’s production scripts (Scripts/R/ForSCOPE/1-getSCOPE-v3_withChunck.R) to keep the leaf-inclination parameters LIDFa/LIDFb jointly plausible instead of independently uniform:

lidf_cor <- suppressMessages(ToolsRTM::getCor(
  n_inputs = 2, setseed = 7, distribution = "Uniform",
  nLUT = n_samples, rho = 0.20, Varnames = c("LIDFa", "LIDFb"),
  MinRange = c(-0.5, -0.5), MaxRange = c(0.2, 0.2)
))
LUT$LIDFa <- lidf_cor$LUT$LIDFa
LUT$LIDFb <- lidf_cor$LUT$LIDFb
cat("Achieved correlation:", round(cor(LUT$LIDFa, LUT$LIDFb), 2), "\n")
## Achieved correlation: -0.06
plot(LUT$LIDFa, LUT$LIDFb, xlab = "LIDFa", ylab = "LIDFb",
     main = "Correlated leaf-angle-distribution parameters", pch = 19)

The scatter leans along a diagonal instead of filling a square — rho = 0.20 is a modest but real correlation, not independent sampling.

Running many simulations: serial vs. parallel

get.SCOPE() accepts a multi-row LUT directly and loops over it serially (n.LUT must match nrow(LUT)). get.SCOPE.parallel() splits the same LUT into n.cores chunks, runs them on a parallel/doParallel PSOCK cluster, and (optionally) writes CSV outputs itself. Timed back-to-back on the same 10-row LUT:

t0 <- Sys.time()
sims_serial <- SCOPEinR::get.SCOPE(
  LUT = LUT, n.LUT = n_samples, options.SCOPE = scope_options,
  optipar = SCOPEinR::optipar2021.Pro.CX, leaf.model = "fluspect-CX",
  canopy.model = "fourSAIL", get.outputs = "ALL", get.plots = FALSE
)
t_serial <- as.numeric(Sys.time() - t0, units = "secs")

t0 <- Sys.time()
sims_parallel <- SCOPEinR::get.SCOPE.parallel(
  LUT = LUT, options.SCOPE = scope_options, optipar = SCOPEinR::optipar2021.Pro.CX,
  leaf.model = "fluspect-CX", canopy.model = "fourSAIL",
  parallel = TRUE, n.cores = 4, get.outputs = "ALL",
  get.plots = FALSE, get.csv = FALSE
)
## Total simulations: 10 
## Total execution time: 8.91504
t_parallel <- as.numeric(Sys.time() - t0, units = "secs")

cat("Serial (1 core):    ", round(t_serial, 1), "s for", length(sims_serial), "sims\n")
## Serial (1 core):     12.8 s for 10 sims
cat("Parallel (4 cores):  ", round(t_parallel, 1), "s for", length(sims_parallel), "sims\n")
## Parallel (4 cores):   10.7 s for 10 sims

For a LUT this small the parallel path still pays a fixed cost — spinning up a PSOCK cluster and exporting SCOPEinR/ToolsRTM to every worker — so the speed-up is modest here; it becomes worthwhile once nLUT is in the hundreds/thousands, which is the realistic use case for get.LUT() + inversion training data.

The chunked production pattern

For large LUTs, Carlos’s real workflow (Scripts/R/ForSCOPE/1-getSCOPE-v3_withChunck.R) splits the LUT into chunks and calls get.SCOPE.parallel() once per chunk with get.csv = TRUE, so intermediate results are written to disk instead of being held in memory all at once:

# A dedicated, larger LUT for this section -- get.SCOPE.plots() (Section 8)
# needs at least 10 simulations per chunk to have something to bin, and the
# inversion demo below (Section 11) reuses this same data: get.inversion()
# uses 10-fold x 10-repeat cross-validation internally, which needs enough
# rows per fold to fit PLSR/RF reliably -- 20 total (14 in the 70% training
# split) is too few and makes caret::train() fail.
n_prod <- 60
LUT_prod <- getLUT.SCOPE(inputLUT = inputLUT, nLUT = n_prod, setseed = 43)
lidf_cor_prod <- suppressMessages(ToolsRTM::getCor(
  n_inputs = 2, setseed = 8, distribution = "Uniform",
  nLUT = n_prod, rho = 0.20, Varnames = c("LIDFa", "LIDFb"),
  MinRange = c(-0.5, -0.5), MaxRange = c(0.2, 0.2)
))
LUT_prod$LIDFa <- lidf_cor_prod$LUT$LIDFa
LUT_prod$LIDFb <- lidf_cor_prod$LUT$LIDFb

out_root <- paste0(tempdir(), "/scope_outs/")   # trailing slash required by path.out below
dir.create(out_root, showWarnings = FALSE, recursive = TRUE)

n_chunks   <- 2
n_rows     <- nrow(LUT_prod)
chunk_size <- n_rows / n_chunks

for (i in seq_len(n_chunks)) {
  start_row <- ((i - 1) * chunk_size) + 1
  end_row   <- min(i * chunk_size, n_rows)
  LUT_chunk <- LUT_prod[start_row:end_row, ]

  SCOPEinR::get.SCOPE.parallel(
    LUT           = LUT_chunk,
    options.SCOPE = scope_options,
    optipar       = SCOPEinR::optipar2021.Pro.CX,
    leaf.model    = "fluspect-CX",
    canopy.model  = "fourSAIL",
    parallel      = TRUE,
    n.cores       = 4,
    get.outputs   = "ALL",
    get.plots     = FALSE,
    get.csv       = TRUE,
    path.out      = out_root
  )
}
subdirs <- list.dirs(out_root, recursive = FALSE)
basename(subdirs)
## [1] "2026-08-20_17-19-34" "2026-08-20_17-19-58"

get.SCOPE.parallel()/get.SCOPE.outputs() build the timestamped output folder with paste(path.out, current_time, sep = "") — plain string concatenation, not file.path(). path.out must end in /, or the “subfolder” ends up as a sibling directory with the timestamp glued onto its name instead of nested inside it (an easy trap; tempdir() never ends in a slash on its own).

Exploring SCOPE outputs

Each chunk folder holds one CSV per output variable, plus a Parameters/ subfolder with the exact LUT used:

last_subdir <- subdirs[length(subdirs)]
list.files(last_subdir)
##  [1] "aPAR.csv"                      "Eout_spectrum.csv"            
##  [3] "Esky.csv"                      "Esun.csv"                     
##  [5] "fluorescence.csv"              "fluorescence_All_Leaves.csv"  
##  [7] "fluorescence_hemis.csv"        "fluorescence_ReabsCorr.csv"   
##  [9] "fluorescence_scalar.csv"       "fluorescence_scalar_units.csv"
## [11] "fluorescence_scattered.csv"    "fluorescence_shaded.csv"      
## [13] "fluorescence_soil.csv"         "fluorescence_sunlit.csv"      
## [15] "fluxes.csv"                    "fluxes_units.csv"             
## [17] "Lo_spectrum.csv"               "Lo_spectrum_includingF.csv"   
## [19] "Parameters"                    "radiation.csv"                
## [21] "radiation_units.csv"           "rdd.csv"                      
## [23] "rdo.csv"                       "refl.csv"                     
## [25] "reflapp.csv"                   "resistance.csv"               
## [27] "resistance_units.csv"          "rsd.csv"                      
## [29] "rso.csv"                       "SigmaF.csv"                   
## [31] "vegatation.csv"                "vegatation_units.csv"
list.files(file.path(last_subdir, "Parameters"))
## [1] "inputLUT.csv"     "SCOPEversion.txt"

One CSV per output variable (fluxes.csv, reflapp.csv, …), plus the exact input LUT that produced them saved alongside for reproducibility.

Merging chunks

Reading each chunk’s CSV directly and row-binding across chunk folders merges chunked output back into one table:

fluxes_all <- do.call(rbind, lapply(subdirs, function(d) read.csv(file.path(d, "fluxes.csv"))))
dim(fluxes_all)
## [1] 60 14
head(fluxes_all[, c("Rnctot","Actot","Tcave","Rnstot","Tsave")])
##     Rnctot      Actot    Tcave     Rnstot    Tsave
## 1 530.7934 22.7952724 22.14358  39.489187 23.42798
## 2 582.9905 20.2861382 22.46935  10.602271 25.42557
## 3 616.1376 26.2266501 22.12216   9.564002 25.13658
## 4 201.9446  0.9033706 21.84580 272.589532 26.79398
## 5 425.9027 27.2858073 22.87809 104.218383 25.02557
## 6 596.2558 26.0440659 21.19897  17.022793 23.27270

Plotting: get.SCOPE.plots()

get.SCOPE.plots(path.files, plant.trait, get.plots) groups a chunk’s spectra by an input trait (into 8 bins) and saves averaged spectra plots to path.files/Plots.checks/. It needs at least 10 simulations per file to have something to bin. Called with ::: since it’s an internal helper, not part of the exported user-facing API:

SCOPEinR:::get.SCOPE.plots(path.files = last_subdir, plant.trait = c("Vcmax25", "Cab"),
                            get.plots = "reflectance")
SCOPEinR:::get.SCOPE.plots(path.files = last_subdir, plant.trait = c("Vcmax25", "Cab"),
                            get.plots = "fluorescence")
SCOPEinR:::get.SCOPE.plots(path.files = last_subdir, plant.trait = c("Vcmax25", "Cab"),
                            get.plots = "radiance")
png_files <- list.files(file.path(last_subdir, "Plots.checks"), full.names = TRUE)
basename(png_files)
refl_png <- png_files[grepl("^1-reflectance_refl_by_Vcmax25", basename(png_files))]
knitr::include_graphics(refl_png)
refl grouped by Vcmax25, produced by get.SCOPE.plots()

refl grouped by Vcmax25, produced by get.SCOPE.plots()

Reflectance averaged within each of the 8 Vcmax25 bins – a quick way to see whether a trait leaves a visible signature on the spectrum before running a full inversion.

Inverting traits from simulated reflectance

The apparent reflectance (reflapp) SCOPE produces is exactly the kind of “sensor-observed” spectrum inversion works from — same idea as the PROSAIL inversion in the ToolsRTM manual, just generated with a coupled energy-balance model instead of pure radiative transfer. This section runs the short version of that workflow; see ToolsRTM_PROSAIL_tutorial.Rmd for the full 12-algorithm reference and the deep-learning (CNN/Hidden-layers) path.

A 10-row LUT is too small to train anything meaningful. Rather than running a second, separate batch of simulations, this section reuses the 20 already produced by the chunked-parallel run above (out_root, two 10-row chunks) — same effective LUT size, without simulating twice:

inv_subdirs <- list.dirs(out_root, recursive = FALSE)
# A chunk can come back with no reflapp.csv if every simulation in it hit an
# unrecoverable error (e.g. a randomly-drawn LIDFa/LIDFb combination outside
# the physically valid range for fourSAIL) -- drop incomplete chunks rather
# than fail the whole comparison.
inv_subdirs <- inv_subdirs[file.exists(file.path(inv_subdirs, "reflapp.csv"))]
length(inv_subdirs)
## [1] 2
lut_used <- do.call(rbind, lapply(inv_subdirs, function(d) {
  read.csv(file.path(d, "Parameters", "inputLUT.csv"))
}))
reflapp <- do.call(rbind, lapply(inv_subdirs, function(d) {
  read.csv(file.path(d, "reflapp.csv"))
}))

wlS <- SCOPEinR::define.bands()$wlS[1:2001]
colnames(reflapp) <- paste0("R.", wlS)

df.indices <- ToolsRTM::getIndices(reflapp, pattern.rfl = "R.", spectral.domain = "VNIR")
## [1] "Estimating indices using VNIR domain..."
## 
df.indices <- df.indices[, colSums(is.na(df.indices)) == 0]

target_wavelengths <- paste0("R.", seq(400, 800, by = 5))
dataset <- cbind(Cab = lut_used$Cab, reflapp[, target_wavelengths], df.indices)
dataset <- dataset[, !duplicated(names(dataset))]
dim(dataset)
## [1]   60 1968

Two caret-backed algorithms via ToolsRTM::get.inversion(), using the same interface the ToolsRTM manual documents in full for all 12 options:

inputs_inv <- colnames(dataset)[-1]

inv_plsr <- get.inversion(data = dataset, depVar = "Cab", inputs = inputs_inv,
                           algorithm = "PLSR", n.samples = nrow(dataset), seed = 123)

inv_rf   <- get.inversion(data = dataset, depVar = "Cab", inputs = inputs_inv,
                           algorithm = "RF", n.samples = nrow(dataset), seed = 123)

cat("PLSR trained. Object class:", class(inv_plsr)[1], "\n")
## PLSR trained. Object class: list
cat("RF trained. Object class:  ", class(inv_rf)[1], "\n")
## RF trained. Object class:   list

Summary

Want to… Use
Run one simulation get.SCOPE(LUT[1, ], options.SCOPE, optipar, leaf.model, canopy.model, get.outputs, get.plots)
Sample a random parameter set getLUT.SCOPE(inputLUT, nLUT, setseed)
Correlate two traits (e.g. LIDFa/LIDFb) ToolsRTM::getCor(...)
Run many simulations in parallel get.SCOPE.parallel(LUT, ..., parallel = TRUE, n.cores = , get.csv = TRUE, path.out = "somewhere/")
Merge chunked CSV outputs read.csv() + rbind() per chunk folder
Inspect reflectance/fluorescence/radiance by trait SCOPEinR:::get.SCOPE.plots(path.files, plant.trait, get.plots)
Invert a trait from simulated reflapp ToolsRTM::getIndices() + ToolsRTM::get.inversion() — full detail in ToolsRTM_PROSAIL_tutorial.Rmd

Every code chunk above executed for real while writing this manual.