From bc81ef87a2d9cb567365c2ba78c1da469831ab30 Mon Sep 17 00:00:00 2001 From: Jean-Baptiste Feret Date: Thu, 4 Jan 2024 18:08:41 +0100 Subject: [PATCH] # prospect v1.5.0 version following reviews for publication in JOSS ## addition - added a function download_LeafDB to get leaf databases from gitlab repository ## Changes - modified defaulty merit function as one function instead of two - updated documentation - updated function FitSpectralData --- DESCRIPTION | 9 +- NAMESPACE | 5 +- NEWS.md | 13 ++ R/Lib_PROSPECT.R | 50 ++++++-- R/Lib_PROSPECT_Inversion.R | 85 ++++++------ man/CostVal_RMSE.Rd | 21 --- man/Invert_PROSPECT.Rd | 2 +- man/Invert_PROSPECT_subdomain.Rd | 2 +- ...MSE_PROSPECT.Rd => Merit_PROSPECT_RMSE.Rd} | 6 +- man/download_LeafDB.Rd | 19 +++ man/optimal_features_SFS.Rd | 2 +- paper/paper.bib | 79 ++++++++---- paper/paper.md | 89 +++++-------- tests/testthat/test-Invert_PROSPECT.R | 18 +-- vignettes/prospect2.Rmd | 121 +++++++++++++----- vignettes/prospect3.Rmd | 39 ++---- vignettes/prospect4.Rmd | 20 +-- vignettes/prospect5.Rmd | 31 ++--- 18 files changed, 334 insertions(+), 277 deletions(-) delete mode 100644 man/CostVal_RMSE.Rd rename man/{Merit_RMSE_PROSPECT.Rd => Merit_PROSPECT_RMSE.Rd} (91%) create mode 100644 man/download_LeafDB.Rd diff --git a/DESCRIPTION b/DESCRIPTION index 309fff0..7a56446 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: prospect Title: PROSPECT leaf radiative transfer model and inversion routines -Version: 1.5.0 +Version: 1.5.1 Authors@R: c(person(given = "Jean-Baptiste", family = "Feret", email = "jb.feret@teledetection.fr", @@ -25,12 +25,13 @@ Suggests: rmarkdown, testthat (>= 3.0.0) Imports: + data.table, + dplyr, expint, - pracma, future, future.apply, - progress, - dplyr + pracma, + progress VignetteBuilder: knitr URL: https://gitlab.com/jbferet/prospect BugReports: https://gitlab.com/jbferet/prospect/issues diff --git a/NAMESPACE b/NAMESPACE index 7796044..6a27bb1 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,13 +1,12 @@ # Generated by roxygen2: do not edit by hand export(Complete_Input_PROSPECT) -export(CostVal_RMSE) export(FitSpectralData) export(Get_Nprior) export(Invert_PROSPECT) export(Invert_PROSPECT_OPT) export(Invert_PROSPECT_subdomain) -export(Merit_RMSE_PROSPECT) +export(Merit_PROSPECT_RMSE) export(PROSPECT) export(PROSPECT_LUT) export(SetInitParm) @@ -16,11 +15,13 @@ export(calctav) export(check_prospect_parms) export(check_version_prospect) export(colour_to_ansi) +export(download_LeafDB) export(optimal_features_SFS) export(print_msg) export(reshape_lop4inversion) export(tryInversion) import(dplyr) +importFrom(data.table,fread) importFrom(expint,expint) importFrom(future,multisession) importFrom(future,plan) diff --git a/NEWS.md b/NEWS.md index c4b9588..1482f2b 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,16 @@ +# prospect v1.5.0 +version following reviews for publication in JOSS + +## addition +- added a function download_LeafDB to get leaf databases from gitlab repository + +## Changes +- modified defaulty merit function as one function instead of two +- updated documentation +- updated function FitSpectralData + + + # prospect v1.5.0 version following reviews for publication in JOSS diff --git a/R/Lib_PROSPECT.R b/R/Lib_PROSPECT.R index a196397..2757c3d 100644 --- a/R/Lib_PROSPECT.R +++ b/R/Lib_PROSPECT.R @@ -227,16 +227,15 @@ FitSpectralData <- function(lambda, SpecPROSPECT = NULL, lb <- lambda SpecPROSPECT <- SpecPROSPECT %>% filter(SpecPROSPECT$lambda%in%lb) # if UserDomain is defined - if (!is.null(UserDomain)) { - if (UL_Bounds==TRUE) UserDomain <- seq(min(UserDomain), max(UserDomain)) - if (!is.null(Refl)) Refl <- Refl %>% filter(lambda%in%UserDomain) - if (!is.null(Tran)) Tran <- Tran %>% filter(lambda%in%UserDomain) - lambda <- lambda[lambda%in%UserDomain] - # Adjust PROSPECT - SpecPROSPECT <- SpecPROSPECT %>% filter(SpecPROSPECT$lambda%in%UserDomain) - if (any(!UserDomain%in%lambda)){ - message('leaf optics out of range defined by UserDomain') - } + if (is.null(UserDomain)) UserDomain <- lambda + if (UL_Bounds==TRUE) UserDomain <- seq(min(UserDomain), max(UserDomain)) + if (!is.null(Refl)) Refl <- Refl %>% filter(lambda%in%UserDomain) + if (!is.null(Tran)) Tran <- Tran %>% filter(lambda%in%UserDomain) + lambda <- lambda[lambda%in%UserDomain] + # Adjust PROSPECT + SpecPROSPECT <- SpecPROSPECT %>% filter(SpecPROSPECT$lambda%in%UserDomain) + if (any(!UserDomain%in%lambda)){ + message('leaf optics out of range defined by UserDomain') } RT <- reshape_lop4inversion(Refl = Refl, Tran = Tran, @@ -299,6 +298,36 @@ PROSPECT_LUT <- function(Input_PROSPECT, SpecPROSPECT = NULL) { 'Input_PROSPECT' = Input_PROSPECT)) } +#' Complete the list of PROSPECT parameters with default values +#' +#' @param urldb character. URL for online repository where to download data +#' @param dbName character. name of the database available online +#' +#' @return list. Includes leaf chemistry, refl, tran & number of samples +#' @importFrom data.table fread +#' @export + +download_LeafDB <- function(urldb = NULL, + dbName = 'ANGERS'){ + # repository where data are stored + if (is.null(urldb)) urldb <- 'https://gitlab.com/jbferet/myshareddata/raw/master/LOP/' + # download leaf chemistry and optical properties + DataBioch <- data.table::fread(file.path(urldb,dbName,'DataBioch.txt')) + Refl <- data.table::fread(file.path(urldb,dbName,'ReflectanceData.txt')) + Tran <- data.table::fread(file.path(urldb,dbName,'TransmittanceData.txt')) + # Get wavelengths corresponding to the reflectance & transmittance measurements + lambda <- Refl$wavelength + Refl$wavelength <- Tran$wavelength <- NULL + # Get the number of samples + nbSamples <- ncol(Refl) + return(list('DataBioch' = DataBioch, + 'lambda' = lambda, + 'Refl' = Refl, + 'Tran' = Tran, + 'nbSamples' = nbSamples)) +} + + #' Complete the list of PROSPECT parameters with default values #' #' @param Input_PROSPECT input parameters sent to PROSPECT by user @@ -307,6 +336,7 @@ PROSPECT_LUT <- function(Input_PROSPECT, SpecPROSPECT = NULL) { #' #' @return Input_PROSPECT #' @export + Complete_Input_PROSPECT <- function(Input_PROSPECT, Parm2Add, ExpectedParms) { ii <- 0 nbSamples <- length(Input_PROSPECT[[1]]) diff --git a/R/Lib_PROSPECT_Inversion.R b/R/Lib_PROSPECT_Inversion.R index 3fa18f9..39234f9 100644 --- a/R/Lib_PROSPECT_Inversion.R +++ b/R/Lib_PROSPECT_Inversion.R @@ -66,7 +66,7 @@ Invert_PROSPECT <- function(SpecPROSPECT = NULL, N = 1.5, alpha = 40), Parms2Estimate = "ALL", PROSPECT_version = "D", - MeritFunction = "Merit_RMSE_PROSPECT", + MeritFunction = "Merit_PROSPECT_RMSE", xlub = data.frame( CHL = c(1e-4, 150), CAR = c(1e-4, 25), ANT = c(0, 50), BROWN = c(0, 1), @@ -77,14 +77,6 @@ Invert_PROSPECT <- function(SpecPROSPECT = NULL, verbose = FALSE, progressBar = TRUE) { if (is.null(SpecPROSPECT)) SpecPROSPECT <- prospect::SpecPROSPECT_FullRange - # add default values to xlub in case they were not defined - xlub_default <- data.frame(CHL = c(1e-4, 150), CAR = c(1e-4, 25), - ANT = c(0, 50), BROWN = c(0, 1), - EWT = c(1e-8, 0.1), LMA = c(1e-6, .06), - PROT = c(1e-7, .006), CBC = c(1e-6, .054), - N = c(.5, 4), alpha = c(10, 90)) - AddedParm_LB <- setdiff(names(xlub_default),names(xlub)) - for (ad in AddedParm_LB) xlub[[ad]] <- xlub_default[[ad]] # check if list of parameters applicable to PROSPECT version parms_checked <- check_prospect_parms(PROSPECT_version = PROSPECT_version, Parms2Estimate = Parms2Estimate, @@ -115,7 +107,8 @@ Invert_PROSPECT <- function(SpecPROSPECT = NULL, names(parms_checked$InitValues)) updateInitValues <- parms_checked$InitValues updateInitValues[ModifyInit] <- 1.1*updateInitValues[ModifyInit] - res <- tryInversion(x0 = updateInitValues, MeritFunction = MeritFunction, + res <- tryInversion(x0 = updateInitValues, + MeritFunction = MeritFunction, SpecPROSPECT = SpecPROSPECT, Refl = RT$Refl[[idsample]], Tran =RT$Tran[[idsample]], Parms2Estimate = parms_checked$Parms2Estimate, @@ -211,9 +204,13 @@ tryInversion <- function(x0, MeritFunction, SpecPROSPECT, Refl, Tran, res <- tryCatch( { res <- fmincon( - x0 = as.numeric(x0[Parms2Estimate]), fn = MeritFunction, gr =NULL, - SpecPROSPECT = SpecPROSPECT, Refl = Refl, Tran = Tran, - Input_PROSPECT = x0, Parms2Estimate = Parms2Estimate, + x0 = as.numeric(x0[Parms2Estimate]), + fn = MeritFunction, + gr =NULL, + SpecPROSPECT = SpecPROSPECT, + Refl = Refl, Tran = Tran, + Input_PROSPECT = x0, + Parms2Estimate = Parms2Estimate, method = "SQP", A = NULL, b = NULL, Aeq = NULL, beq = NULL, lb = as.numeric(lb), ub = as.numeric(ub), hin = NULL, heq = NULL, tol = Tolerance, maxfeval = 2000, maxiter = 1000) @@ -256,34 +253,20 @@ tryInversion <- function(x0, MeritFunction, SpecPROSPECT, Refl, Tran, #' #' @return fc estimates of the parameters #' @export -Merit_RMSE_PROSPECT <- function(x, SpecPROSPECT, Refl, Tran, - Input_PROSPECT, Parms2Estimate) { +Merit_PROSPECT_RMSE <- function(x, SpecPROSPECT, + Refl, Tran, + Input_PROSPECT, + Parms2Estimate) { x[x < 0] <- 0 Input_PROSPECT[Parms2Estimate] <- x RT <- do.call("PROSPECT", c(list(SpecPROSPECT = SpecPROSPECT), Input_PROSPECT)) - fc <- CostVal_RMSE(RT, Refl, Tran) + fcr <- fct <- 0 + if (!is.null(Refl)) fcr <- sqrt(sum((Refl - RT$Reflectance)**2) / length(RT$Reflectance)) + if (!is.null(Tran)) fct <- sqrt(sum((Tran - RT$Transmittance)**2) / length(RT$Transmittance)) + fc <- fcr + fct return(fc) } -#' Value of the cost criterion to minimize during PROSPECT inversion -#' @param RT list. Simulated reflectance and transmittance -#' @param Refl numeric. Reflectance on which PROSPECT ins inverted -#' @param Tran numeric. Transmittance on which PROSPECT ins inverted -#' -#' @return fc sum of squared difference bw simulated and measured leaf optics -#' @export -CostVal_RMSE <- function(RT, Refl, Tran) { - if (is.null(Tran)) { - fc <- sqrt(sum((Refl - RT$Reflectance)**2) / length(RT$Reflectance)) - } else if (is.null(Refl)) { - fc <- sqrt(sum((Tran - RT$Transmittance)**2) / length(RT$Transmittance)) - } else { - fc <- sqrt(sum((Refl - RT$Reflectance)**2) / length(RT$Reflectance) + sum((Tran - RT$Transmittance)**2) / length(RT$Transmittance)) - } - return(fc) -} - - #' This function defines a regression model to estimate N from R only or T only #' @param lambda numeric. spectral bands corresponding to experimental data #' @param SpecPROSPECT list. Includes optical constants refractive index, @@ -537,7 +520,7 @@ optimal_features_SFS <- function(Refl = NULL, Tran = NULL, lambda, BiochTruth, N = c(.5, 4), alpha = c(10, 90)), SpecPROSPECT, spectral_domain, spectral_width, number_features,PROSPECT_version = 'D', - MeritFunction = "Merit_RMSE_PROSPECT", + MeritFunction = "Merit_PROSPECT_RMSE", Est_alpha = FALSE, verbose = FALSE, nbCPU = 1,Continue = FALSE, AlreadyDone = NULL){ @@ -652,7 +635,7 @@ Invert_PROSPECT_subdomain <- function(New_Features, Refl, Tran, SpecPROSPECT, PROT = 0.001, CBC = 0.009, N = 1.5, alpha = 40), PROSPECT_version = "D", - MeritFunction = "Merit_RMSE_PROSPECT", + MeritFunction = "Merit_PROSPECT_RMSE", Est_Brown_Pigments = FALSE, Est_alpha = FALSE, xlub = data.frame( @@ -839,9 +822,12 @@ reshape_lop4inversion <- function(Refl, Tran, SpecPROSPECT){ #' @return list of parameters to estimate & corresponding lower/upper boundaries #' @export -check_prospect_parms <- function(PROSPECT_version, Parms2Estimate, - Est_Brown_Pigments, Est_alpha, - xlub, InitValues){ +check_prospect_parms <- function(PROSPECT_version, + Parms2Estimate, + Est_Brown_Pigments, + Est_alpha, + xlub, + InitValues){ # check if version required is available if (!PROSPECT_version %in% c('D', 'PRO')) print_msg(cause = 'WrongVersion') @@ -851,6 +837,24 @@ check_prospect_parms <- function(PROSPECT_version, Parms2Estimate, # add brown pigments and alpha angle if required if (Est_Brown_Pigments==TRUE) allParms <- c(allParms, "BROWN") if (Est_alpha==TRUE) allParms <- c(allParms, "alpha") + + # add default values to xlub in case they were not defined + xlub_default <- data.frame(CHL = c(1e-4, 150), CAR = c(1e-4, 25), + ANT = c(0, 50), BROWN = c(0, 1), + EWT = c(1e-8, 0.1), LMA = c(1e-6, .06), + PROT = c(1e-7, .006), CBC = c(1e-6, .054), + N = c(.5, 4), alpha = c(10, 90)) + AddedParm_LB <- setdiff(names(xlub_default),names(xlub)) + for (ad in AddedParm_LB) xlub[[ad]] <- xlub_default[[ad]] + + InitValues_default = data.frame(CHL = 40, CAR = 10, + ANT = 0.1, BROWN = 0.0, + EWT = 0.01, LMA = 0.01, + PROT = 0.001, CBC = 0.009, + N = 1.5, alpha = 40) + AddedParm_init <- setdiff(names(InitValues_default),names(InitValues)) + for (ad in AddedParm_init) InitValues[[ad]] <- InitValues_default[[ad]] + # if 'ALL' is provided, then assess all parameters available if ("ALL" %in% Parms2Estimate) Parms2Estimate <- allParms # if unknown parameter is provided, then warn @@ -872,6 +876,7 @@ check_prospect_parms <- function(PROSPECT_version, Parms2Estimate, InitValues$PROT <- InitValues$CBC <- 0 if (is.null(InitValues$LMA)) InitValues$LMA <- 0.01 } + if (Est_Brown_Pigments==FALSE) InitValues$BROWN <- 0 xlub <- data.frame(xlub[Parms2Estimate]) lb <- xlub[1,] ub <- xlub[2,] diff --git a/man/CostVal_RMSE.Rd b/man/CostVal_RMSE.Rd deleted file mode 100644 index ca61fac..0000000 --- a/man/CostVal_RMSE.Rd +++ /dev/null @@ -1,21 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/Lib_PROSPECT_Inversion.R -\name{CostVal_RMSE} -\alias{CostVal_RMSE} -\title{Value of the cost criterion to minimize during PROSPECT inversion} -\usage{ -CostVal_RMSE(RT, Refl, Tran) -} -\arguments{ -\item{RT}{list. Simulated reflectance and transmittance} - -\item{Refl}{numeric. Reflectance on which PROSPECT ins inverted} - -\item{Tran}{numeric. Transmittance on which PROSPECT ins inverted} -} -\value{ -fc sum of squared difference bw simulated and measured leaf optics -} -\description{ -Value of the cost criterion to minimize during PROSPECT inversion -} diff --git a/man/Invert_PROSPECT.Rd b/man/Invert_PROSPECT.Rd index 23eed7f..1331180 100644 --- a/man/Invert_PROSPECT.Rd +++ b/man/Invert_PROSPECT.Rd @@ -13,7 +13,7 @@ Invert_PROSPECT( 0.01, PROT = 0.001, CBC = 0.009, N = 1.5, alpha = 40), Parms2Estimate = "ALL", PROSPECT_version = "D", - MeritFunction = "Merit_RMSE_PROSPECT", + MeritFunction = "Merit_PROSPECT_RMSE", xlub = data.frame(CHL = c(1e-04, 150), CAR = c(1e-04, 25), ANT = c(0, 50), BROWN = c(0, 1), EWT = c(1e-08, 0.1), LMA = c(1e-08, 0.06), PROT = c(1e-07, 0.006), CBC = c(1e-06, 0.054), N = c(0.5, 4), alpha = c(10, 90)), diff --git a/man/Invert_PROSPECT_subdomain.Rd b/man/Invert_PROSPECT_subdomain.Rd index cc693e9..189b12a 100644 --- a/man/Invert_PROSPECT_subdomain.Rd +++ b/man/Invert_PROSPECT_subdomain.Rd @@ -18,7 +18,7 @@ Invert_PROSPECT_subdomain( InitValues = data.frame(CHL = 40, CAR = 10, ANT = 0.1, BROWN = 0.01, EWT = 0.01, LMA = 0.01, PROT = 0.001, CBC = 0.009, N = 1.5, alpha = 40), PROSPECT_version = "D", - MeritFunction = "Merit_RMSE_PROSPECT", + MeritFunction = "Merit_PROSPECT_RMSE", Est_Brown_Pigments = FALSE, Est_alpha = FALSE, xlub = data.frame(CHL = c(1e-04, 150), CAR = c(1e-04, 25), ANT = c(0, 50), BROWN = c(0, diff --git a/man/Merit_RMSE_PROSPECT.Rd b/man/Merit_PROSPECT_RMSE.Rd similarity index 91% rename from man/Merit_RMSE_PROSPECT.Rd rename to man/Merit_PROSPECT_RMSE.Rd index 45ef67c..2953842 100644 --- a/man/Merit_RMSE_PROSPECT.Rd +++ b/man/Merit_PROSPECT_RMSE.Rd @@ -1,10 +1,10 @@ % Generated by roxygen2: do not edit by hand % Please edit documentation in R/Lib_PROSPECT_Inversion.R -\name{Merit_RMSE_PROSPECT} -\alias{Merit_RMSE_PROSPECT} +\name{Merit_PROSPECT_RMSE} +\alias{Merit_PROSPECT_RMSE} \title{Merit function for PROSPECT inversion} \usage{ -Merit_RMSE_PROSPECT( +Merit_PROSPECT_RMSE( x, SpecPROSPECT, Refl, diff --git a/man/download_LeafDB.Rd b/man/download_LeafDB.Rd new file mode 100644 index 0000000..3f97028 --- /dev/null +++ b/man/download_LeafDB.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/Lib_PROSPECT.R +\name{download_LeafDB} +\alias{download_LeafDB} +\title{Complete the list of PROSPECT parameters with default values} +\usage{ +download_LeafDB(urldb = NULL, dbName = "ANGERS") +} +\arguments{ +\item{urldb}{character. URL for online repository where to download data} + +\item{dbName}{character. name of the database available online} +} +\value{ +list. Includes leaf chemistry, refl, tran & number of samples +} +\description{ +Complete the list of PROSPECT parameters with default values +} diff --git a/man/optimal_features_SFS.Rd b/man/optimal_features_SFS.Rd index 70c6d42..5976cb8 100644 --- a/man/optimal_features_SFS.Rd +++ b/man/optimal_features_SFS.Rd @@ -21,7 +21,7 @@ optimal_features_SFS( spectral_width, number_features, PROSPECT_version = "D", - MeritFunction = "Merit_RMSE_PROSPECT", + MeritFunction = "Merit_PROSPECT_RMSE", Est_alpha = FALSE, verbose = FALSE, nbCPU = 1, diff --git a/paper/paper.bib b/paper/paper.bib index e472634..b0d5d2f 100644 --- a/paper/paper.bib +++ b/paper/paper.bib @@ -14,6 +14,23 @@ @article{allen1970 pages = {542--547}, } + +@article{colombo_estimation_2008, + title = {Estimation of leaf and canopy water content in poplar plantations by means of hyperspectral indices and inverse modeling}, + volume = {112}, + issn = {00344257}, + url = {http://linkinghub.elsevier.com/retrieve/pii/S0034425707004257}, + doi = {10.1016/j.rse.2007.09.005}, + language = {en}, + number = {4}, + urldate = {2017-07-06}, + journal = {Remote Sensing of Environment}, + author = {Colombo, R and Meroni, M and Marchesi, A and Busetto, L and Rossini, M and Giardino, C and Panigada, C}, + month = apr, + year = {2008}, + pages = {1820--1834}, +} + @article{feret2008, title = {{PROSPECT}-4 and 5: {Advances} in the leaf optical properties model separating photosynthetic pigments}, volume = {112}, @@ -95,22 +112,20 @@ @article{jacquemoud1990 pages = {75--91}, } -@article{spafford2021, - title = {Spectral subdomains and prior estimation of leaf structure improves {PROSPECT} inversion on reflectance or transmittance alone}, - volume = {252}, - issn = {00344257}, - url = {https://linkinghub.elsevier.com/retrieve/pii/S0034425720305496}, - doi = {10.1016/j.rse.2020.112176}, - language = {en}, - urldate = {2020-11-16}, +@article{jacquemoud_prospect+_2009, + title = {{PROSPECT}+ {SAIL} models: {A} review of use for vegetation characterization}, + volume = {113}, + shorttitle = {{PROSPECT}+ {SAIL} models}, + url = {http://www.sciencedirect.com/science/article/pii/S0034425709000765}, + doi = {doi:10.1016/j.rse.2008.01.026}, + urldate = {2016-01-03}, journal = {Remote Sensing of Environment}, - author = {Spafford, Lynsay and le Maire, Guerric and MacDougall, Andrew and de Boissieu, Florian and Féret, Jean-Baptiste}, - month = jan, - year = {2021}, - pages = {112176}, + author = {Jacquemoud, Stéphane and Verhoef, Wout and Baret, Frédéric and Bacour, Cédric and Zarco-Tejada, Pablo J. and Asner, Gregory P. and François, Christophe and Ustin, Susan L.}, + year = {2009}, + pages = {S56--S66}, + file = {[PDF] from csic.es:C\:\\Users\\Feret\\Zotero\\storage\\76ATRJNN\\Jacquemoud et al. - 2009 - PROSPECT+ SAIL models A review of use for vegetat.pdf:application/pdf;Snapshot:C\:\\Users\\Feret\\Zotero\\storage\\UX33MX6F\\S0034425709000765.html:text/html}, } - @article{jay_physically-based_2016, title = {A physically-based model for retrieving foliar biochemistry and leaf orientation using close-range imaging spectroscopy}, volume = {177}, @@ -126,21 +141,35 @@ @article{jay_physically-based_2016 pages = {220--236}, } - -@article{jacquemoud_prospect+_2009, - title = {{PROSPECT}+ {SAIL} models: {A} review of use for vegetation characterization}, - volume = {113}, - shorttitle = {{PROSPECT}+ {SAIL} models}, - url = {http://www.sciencedirect.com/science/article/pii/S0034425709000765}, - doi = {doi:10.1016/j.rse.2008.01.026}, +@article{li_retrieval_2011, + title = {Retrieval of leaf biochemical parameters using {PROSPECT} inversion: {A} new approach for alleviating ill-posed problems}, + volume = {49}, + shorttitle = {Retrieval of leaf biochemical parameters using {PROSPECT} inversion}, + url = {http://ieeexplore.ieee.org/xpls/abs_all.jsp?arnumber=5734840}, + doi = {10.1109/TGRS.2011.2109390}, + number = {7}, urldate = {2016-01-03}, - journal = {Remote Sensing of Environment}, - author = {Jacquemoud, Stéphane and Verhoef, Wout and Baret, Frédéric and Bacour, Cédric and Zarco-Tejada, Pablo J. and Asner, Gregory P. and François, Christophe and Ustin, Susan L.}, - year = {2009}, - pages = {S56--S66}, - file = {[PDF] from csic.es:C\:\\Users\\Feret\\Zotero\\storage\\76ATRJNN\\Jacquemoud et al. - 2009 - PROSPECT+ SAIL models A review of use for vegetat.pdf:application/pdf;Snapshot:C\:\\Users\\Feret\\Zotero\\storage\\UX33MX6F\\S0034425709000765.html:text/html}, + journal = {IEEE Transactions on Geoscience and Remote Sensing,}, + author = {Li, Pingheng and Wang, Quan}, + year = {2011}, + pages = {2499--2506}, + file = {Snapshot:C\:\\Users\\Feret\\Zotero\\storage\\C2VG2NRM\\login.html:text/html}, } +@article{spafford2021, + title = {Spectral subdomains and prior estimation of leaf structure improves {PROSPECT} inversion on reflectance or transmittance alone}, + volume = {252}, + issn = {00344257}, + url = {https://linkinghub.elsevier.com/retrieve/pii/S0034425720305496}, + doi = {10.1016/j.rse.2020.112176}, + language = {en}, + urldate = {2020-11-16}, + journal = {Remote Sensing of Environment}, + author = {Spafford, Lynsay and le Maire, Guerric and MacDougall, Andrew and de Boissieu, Florian and Féret, Jean-Baptiste}, + month = jan, + year = {2021}, + pages = {112176}, +} @article{verhoef_coupled_2007, title = {Coupled soil–leaf-canopy and atmosphere radiative transfer modeling to simulate hyperspectral multi-angular surface reflectance and {TOA} radiance data}, diff --git a/paper/paper.md b/paper/paper.md index c67abd6..358ad29 100644 --- a/paper/paper.md +++ b/paper/paper.md @@ -121,18 +121,22 @@ iterative optimization is the most widespread method to invert PROSPECT and assess leaf chemistry and structure from their optical properties. It usually takes less than 1 second to perform PROSPECT inversion, which is acceptable when processing experimental datasets of leaf optical properties including 100-1000 samples. -Iterative optimization aims at minimizing a cost function comparing +Iterative optimization aims at minimizing a merit function comparing measured and simulated leaf optical properties. This procedure is based on the function `fmincon` included in the package `pracma`. -Various inversion strategies using iterative optimization are described in the literature. -These inversion strategies differ either by the cost function, or by the +Various inversion strategies using iterative optimization are described in the +literature [@colombo_estimation_2008, @li_retrieval_2011, @feret2019]. +These inversion strategies differ either by the merit function, or by the selection of specific spectral domains used to retrieve one or several leaf biophysical properties, or by the introduction of prior information. -The default cost function, `CostVal_RMSE`, corresponds to the root mean square +The default merit function, `Merit_PROSPECT_RMSE`, corresponds to the root mean square of the mean quadratic difference between measured and simulated leaf optical properties (reflectance and/or transmittance). -Users can define their own cost function. +Users can define their own merit function, as long as their function uses the +same input and output variables as `Merit_PROSPECT_RMSE`. +The online documentation of the package provides an example of alternative merit function. + Table \ref{table:2} provides information on the optimal spectral range used to assess leaf chemical constituents from their optical properties, as identified by [@feret2019; @spafford2021]. @@ -155,8 +159,6 @@ LMA: leaf mass per area; PROT: proteins; CBC: carbon based constituents).\label{ # Example 1: running PROSPECT in forward mode -## Individual simulations - PROSPECT simulates leaf directional-hemispherical reflectance and transmittance from the leaf structure parameter and a combination of chemical constituents when running in forward mode. @@ -209,54 +211,17 @@ LRT_VNIR <- PROSPECT(SpecPROSPECT = VNIR$SpecPROSPECT, N = 1.4, CHL = 30, CAR = 6, EWT = 0.02, LMA = 0.01) ``` -## Simulation of a Look-Up-Table - -Look-Up-Tables (LUTs) are widely used in order to infer leaf characteristics from PROSPECT. -The function `PROSPECT_LUT` computes a LUT based on a list of input parameters. -Undefined parameters are set to their default value. Vectors of values are expected -to be the same length. -The output of `PROSPECT_LUT` is a list containing a data frame including the -input parameters, a reflectance data frame and a transmittance data frame. - -```r -# define input parametrs for PROSPECT -Input_PROSPECT <- data.frame('CHL' = 100*runif(1000), - 'CAR' = 25*runif(1000), - 'ANT' = 2*runif(1000), - 'EWT' = 0.04*runif(1000), - 'LMA' = 0.02*runif(1000), - 'N' = 1+2*runif(1000)) -# produce a LUT defined over the VSWIR domain covered by PROSPECT -LUT <- PROSPECT_LUT(Input_PROSPECT = Input_PROSPECT) -# produce a LUT defined over the VNIR domain -LUT_VNIR <- PROSPECT_LUT(SpecPROSPECT = VNIR$SpecPROSPECT, - Input_PROSPECT = Input_PROSPECT) -``` - # Example 2: PROSPECT inversion using iterative optimization -The package `prospect` offers possibilities to adjust these parameters for inversion. -Here, we will illustrate different types of inversion with an experimental database +We will illustrate different types of inversion with an experimental database named __ANGERS__. This database was used to calibrate the model PROSPECT, and is among the most popular public data sets in the domain of leaf spectroscopy. - -## Downloading the ANGERS experimental dataset - -A version of the ANGERS data set can be downloaded from a gitlab repository. +The function `download_LeafDB` allows users to directly download this database +from a gitlab repository. ```r -# use data.table library -library(data.table) -# download ANGERS leaf optics database from gitlab repository -gitlab_Rep <- 'https://gitlab.com/jbferet/myshareddata/raw/master/LOP' -# download leaf biochemical constituents and leaf optical properties -DataBioch <- data.table::fread(file.path(gitlab_Rep,'ANGERS/DataBioch.txt')) -Refl<- data.table::fread(file.path(gitlab_Rep,'ANGERS/ReflectanceData.txt')) -Tran <- data.table::fread(file.path(gitlab_Rep,'ANGERS/TransmittanceData.txt')) -# Get the wavelengths corresponding to reflectance and transmittance measurements -lambda <- Refl$wavelength -Refl$wavelength <- Tran$wavelength <- NULL +LeafDB <- download_LeafDB(dbName = 'ANGERS') ``` ## PROSPECT inversion using the full spectral information @@ -272,23 +237,26 @@ domain covered by the leaf optical properties. ```r -# Adjust spectral domain for SpecPROSPECT to fit leaf optical properties +# Adjust spectral domain for SpecPROSPECT to fit leaf optical properties SubData <- FitSpectralData(SpecPROSPECT = SpecPROSPECT_FullRange, - Refl = Refl, Tran = Tran, - lambda = lambda, UserDomain = lambda) + Refl = LeafDB$Refl, + Tran = LeafDB$Tran, + lambda = LeafDB$lambda, + UserDomain = LeafDB$lambda) ``` The main inversion procedure is called with the function `Invert_PROSPECT`, -which minimizes a cost function. -They can also select the biophysical -properties to assess: +which minimizes a merit function (defined by `Merit_PROSPECT_RMSE` as default option). The full set or a selection of chemical constituents can be assessed from PROSPECT inversion. -This list of parameters is defined with the variable `Parms2Estimate`. -The default parameterization of PROSPECT inversion assesses all parmeters -listed in Table \ref{table:1} except BROWN, which can be assessed by setting -`Est_Brown_Pigments = TRUE` as input for `Invert_PROSPECT`. +This list of PROSPECT parameters assess from PROSPECT inversion is defined with +the variable `Parms2Estimate`. +The default parameterization of PROSPECT inversion assesses all parameters, +including the N structure parameter and all chemical constituents listed in +Table \ref{table:1} except BROWN. +BROWN which can also be assessed by setting `Est_Brown_Pigments = TRUE` as +input for `Invert_PROSPECT`. The value set for parameters which are not assessed is then defined with the `InitValues` input variable. @@ -310,8 +278,9 @@ by user. ```r # Assess a set of parameters using PROSPECT inversion with optimal spectral domains Parms2Estimate <- c('CHL', 'CAR', 'EWT', 'LMA') -res_opt_WL <- Invert_PROSPECT_OPT(lambda = lambda, - Refl = Refl, Tran = Tran +res_opt_WL <- Invert_PROSPECT_OPT(lambda = LeafDB$lambda, + Refl = LeafDB$Refl, + Tran = LeafDB$Tran Parms2Estimate = Parms2Estimate) ``` diff --git a/tests/testthat/test-Invert_PROSPECT.R b/tests/testthat/test-Invert_PROSPECT.R index c919c63..f7bd3d5 100644 --- a/tests/testthat/test-Invert_PROSPECT.R +++ b/tests/testthat/test-Invert_PROSPECT.R @@ -6,7 +6,7 @@ test_that("PROSPECT-D inversion produces accurate biophysical assessment", { leafBP <- Invert_PROSPECT(Refl = lrt$Reflectance, Tran = lrt$Transmittance) expect_true(abs(leafBP$N-BPinit$N)<1e-5) - expect_true(abs(leafBP$CHL-BPinit$CHL)<1e-4) + expect_true(abs(leafBP$CHL-BPinit$CHL)<1e-3) expect_true(abs(leafBP$CAR-BPinit$CAR)<1e-4) expect_true(abs(leafBP$ANT-BPinit$ANT)<1e-4) expect_true(abs(leafBP$EWT-BPinit$EWT)<1e-6) @@ -17,12 +17,12 @@ test_that("PROSPECT-D inversion produces accurate biophysical assessment", { lrt$Transmittance <- lrt$Transmittance*(1+rnorm(length(lrt$Transmittance), 0,0.01)) - leafBP <- Invert_PROSPECT(Refl = lrt$Reflectance, - Tran = lrt$Transmittance) - expect_true(abs(leafBP$N-BPinit$N)<1e-3) - expect_true(abs(leafBP$CHL-BPinit$CHL)<1e-1) - expect_true(abs(leafBP$CAR-BPinit$CAR)<1e-1) - expect_true(abs(leafBP$ANT-BPinit$ANT)<1e-1) - expect_true(abs(leafBP$EWT-BPinit$EWT)<1e-4) - expect_true(abs(leafBP$LMA-BPinit$LMA)<1e-4) + leafBP_noise <- Invert_PROSPECT(Refl = lrt$Reflectance, + Tran = lrt$Transmittance) + expect_true(abs(leafBP_noise$N-BPinit$N)<1e-3) + expect_true(abs(leafBP_noise$CHL-BPinit$CHL)<1e-1) + expect_true(abs(leafBP_noise$CAR-BPinit$CAR)<1e-1) + expect_true(abs(leafBP_noise$ANT-BPinit$ANT)<1e-1) + expect_true(abs(leafBP_noise$EWT-BPinit$EWT)<1e-4) + expect_true(abs(leafBP_noise$LMA-BPinit$LMA)<1e-4) }) diff --git a/vignettes/prospect2.Rmd b/vignettes/prospect2.Rmd index 5670406..250cdb0 100644 --- a/vignettes/prospect2.Rmd +++ b/vignettes/prospect2.Rmd @@ -58,28 +58,35 @@ The list of input variables for inversion is : * `SpecPROSPECT`: data frame including the spectral bands, refractive index and specific absorption coefficients. `SpecPROSPECT = NULL` uses the spectral domain from 400 nm to 2500 nm. Simulation and inversion on different spectral domains can be performed by adjusting `SpecPROSPECT` -* `Refl`: numeric: individual leaf reflectance corresponding to the spectral domain defined in `SpecPROSPECT`. Set to `NULL` if inversion on transmittance only -* `Tran`: numeric: individual leaf transmittance corresponding to the spectral domain defined in `SpecPROSPECT`. Set to `NULL` if inversion on reflectance only +* `Refl`: numeric: individual leaf reflectance corresponding to the spectral +domain defined in `SpecPROSPECT`. Set to `NULL` if inversion on transmittance only +* `Tran`: numeric: individual leaf transmittance corresponding to the spectral +domain defined in `SpecPROSPECT`. Set to `NULL` if inversion on reflectance only * `Parms2Estimate` list. Parameters to estimate. Set to 'ALL' by default. * `InitValues` data frame including initial values of PROSPECT input parameters. -During optimization, they are used either as initialization values for parameters to estimate, or as fix values for other parameters. +During optimization, they are used either as initialization values for parameters +to estimate, or as fix values for other parameters. Parameters are not taken into account when not compatible with PROSPECT_version. * `PROSPECT_version` character. Corresponds to the PROSPECT version. -Should be one of the following versions: '5', '5B', 'D', 'DB', 'PRO', 'PROB'. -Use the version ending with 'B' if you want to estimate brown pigments. -Versions '5' and '5B' are actually based on the specific absorption coefficients of chlorophylls and carotenoids, and the refractive index from PROSPECT-D. -`ANT` is then set to 0 during inversion. -* `MeritFunction` character. name of the function to be used as merit function with given criterion to minimize (default = RMSE) +Should be one of the following versions: 'D', 'PRO'. +* `MeritFunction` character. name of the function to be used as merit function +with given criterion to minimize (default = Merit_PROSPECT_RMSE) * `xlub` data frame. Boundaries of the parameters to estimate. -The data frame must have columns corresponding to \code{Parms2Estimate} first line being the lower boundaries and second line the upper boundaries. -* `alphaEst` boolean. Should `alpha` be estimated or not? Keep in mind that most published results use `alpha` with its default value. +The data frame must have columns corresponding to \code{Parms2Estimate} first +line being the lower boundaries and second line the upper boundaries. +* `Est_Brown_Pigments` boolean. Should `BROWN` be estimated or not? Brown pigments are found during later stages of senescent leaves, but they are usually not found in juvenile, mature or early senescent leaves. +* `Est_alpha` boolean. Should `alpha` be estimated or not? Keep in mind that most published results use `alpha` with its default value. ## Output variables `Invert_PROSPECT` returns a list of estimated PROSPECT parameters defined in `Parms2Estimate`. -## Run PROSPECT-D inversion over the full spectral domain +# Running PROSPECT-D inversion + +## Over the full spectral domain All parameters are estimated, except `alpha` and `BROWN` which are set to their default value. +Brown pigments are characteristic of late senescent stages. +They are absent for all other development stages, hence the assessment of brown pigments is not systematic here. ```{r prospect inverse mode 1} # simulate leaf optical properties @@ -87,19 +94,29 @@ LRT_D <- PROSPECT(CHL = 45, CAR = 10, ANT = 0.2, EWT = 0.012, LMA = 0.01, N = 1.3) # define set of parameters to be estimated Parms2Estimate <- 'ALL' -# define initial values for the inversion (should not impact final results) -InitValues <- data.frame(CHL = 40, CAR = 10, ANT = 0.1, BROWN = 0, - EWT = 0.01, LMA = 0.01, N = 1.5) # invert PROSPECT with simulated leaf optical properties OutPROSPECT <- Invert_PROSPECT(Refl = LRT_D$Reflectance, Tran = LRT_D$Transmittance, - InitValues = InitValues, Parms2Estimate = Parms2Estimate, PROSPECT_version = 'D') + + +# simulate leaf optical properties with brown pigments +LRT_D_brown <- PROSPECT(CHL = 45, CAR = 10, ANT = 0.2, BROWN = 0.1, + EWT = 0.012, LMA = 0.01, N = 1.3) +# define set of parameters to be estimated +Parms2Estimate <- 'ALL' +# invert PROSPECT with simulated leaf optical properties and specify that brown +# pigment should be assessed in addition to all other constituents +OutPROSPECT <- Invert_PROSPECT(Refl = LRT_D_brown$Reflectance, + Tran = LRT_D_brown$Transmittance, + Parms2Estimate = Parms2Estimate, + PROSPECT_version = 'D', + Est_Brown_Pigments = TRUE) ``` -## Run PROSPECT-D inversion over the VNIR domain +## over the VNIR domain All parameters are estimated, except `alpha` and `BROWN` which are set to their default value. `FitSpectralData` directly adjusts reflectance, transmittance and optical constants to user defined spectral bands. These spectral bands should be integer values between 400 and 2500, and can be continuous or not. @@ -130,22 +147,22 @@ OutPROSPECT <- Invert_PROSPECT(SpecPROSPECT = SubData$SpecPROSPECT, PROSPECT_version = 'D') ``` -## Run PROSPECT-D inversion over the VNIR domain with LMA and EWT value set +## over the VNIR domain with LMA and EWT value set Only pigments and `N` are estimated. The same optical properties as previous example are used. ```{r prospect inverse mode 3} # define set of parameters to be estimated -Parms2Estimate <- c("CHL", "CAR", "ANT", "N") -# define initial values for the inversion -InitValues <- data.frame(CHL = 40, CAR = 10, ANT = 0.1, BROWN = 0, - EWT = 0.01, LMA = 0.01, N = 1.5) +Parms2Estimate <- c('CHL', 'CAR', 'ANT', 'N') +# define set of parameters to be estimated +InitValues <- c('EWT' = 0.01, 'LMA' = 0.01) # invert PROSPECT with simulated leaf optical properties OutPROSPECT <- Invert_PROSPECT(SpecPROSPECT = SubData$SpecPROSPECT, Refl = SubData$Refl, Tran = SubData$Tran, Parms2Estimate = Parms2Estimate, - PROSPECT_version = 'D') + PROSPECT_version = 'D', + InitValues = InitValues) ``` ## Run PROSPECT-D inversion over the SWIR domain between 1700 nm and 2400 nm @@ -154,10 +171,7 @@ The same optical properties as previous example are used. ```{r prospect inverse mode 4} # define set of parameters to be estimated -Parms2Estimate <- c("EWT", "LMA", "N") -# define initial values for the inversion -InitValues <- data.frame(CHL = 0, CAR = 0, ANT = 0, BROWN = 0, - EWT = 0.01, LMA = 0.008, N = 1.5) +Parms2Estimate <- c('EWT', 'LMA', 'N') # define spectral subdomain SpectralSubDomain <- seq(1700,2400) # adjust spectral domain @@ -168,14 +182,13 @@ SubData <- FitSpectralData(lambda = LRT_D$wvl, # invert PROSPECT with simulated leaf optical properties OutPROSPECT <- Invert_PROSPECT(SpecPROSPECT = SubData$SpecPROSPECT, - InitValues = InitValues, Refl = SubData$Refl, Tran = SubData$Tran, Parms2Estimate = Parms2Estimate, PROSPECT_version = 'D') ``` -# Invert PROSPECT-D using the optimal configuration for the estimation of leaf constituents +# Invert PROSPECT-D using optimal spectral domains The function `Invert_PROSPECT_OPT` automatically sets the optimal spectral domains during inversion for all constituents to be estimated. Optimal spectral domains and configuration are defined in [Féret et al. (2019)](https://doi.org/10.1016/j.rse.2018.11.002), [Féret et al. (2021)](https://doi.org/10.1016/j.rse.2020.112173), and [Spafford et al. (2021)](https://doi.org/10.1016/j.rse.2020.112176). @@ -186,20 +199,16 @@ Optimal spectral domains and configuration are defined in [Féret et al. (2019)] ```{r prospect inverse mode 5} # define set of parameters to be estimated Parms2Estimate <- c('CHL','CAR','ANT','EWT','LMA') -# define initial values for the inversion -InitValues <- data.frame(CHL = 40, CAR = 8, ANT = 0.1, BROWN = 0, - EWT = 0.01, LMA = 0.01, N = 1.5) # call Invert_PROSPECT_OPT to get optimal estimation of leaf parameters # using optimal spectral domains ParmEst <- Invert_PROSPECT_OPT(lambda = LRT_D$wvl, Refl = LRT_D$Reflectance, Tran = LRT_D$Transmittance, PROSPECT_version = 'D', - Parms2Estimate = Parms2Estimate, - InitValues = InitValues) + Parms2Estimate = Parms2Estimate) ``` -### run PROSPECT-PRO inversion +# run PROSPECT-PRO inversion Such definition of optimal spectral domains can also be set manually. For example, here is how to estimate protein content from leaf optical properties using the optimal spectral domain defined in [Féret et al. (2021)](https://doi.org/10.1016/j.rse.2020.112173). @@ -219,7 +228,7 @@ SubData <- FitSpectralData(lambda = LRT_PRO$wvl, UserDomain = SpectralSubDomain, UL_Bounds = TRUE) -Parms2Estimate <- c("EWT", "PROT", "CBC", "N") +Parms2Estimate <- c('EWT', 'PROT', 'CBC', 'N') # invert PROSPECT with simulated leaf optical properties OutPROSPECT <- Invert_PROSPECT(SpecPROSPECT = SubData$SpecPROSPECT, Refl = SubData$Refl, @@ -228,3 +237,45 @@ OutPROSPECT <- Invert_PROSPECT(SpecPROSPECT = SubData$SpecPROSPECT, PROSPECT_version = 'PRO') ``` + +# Run PROSPECT inversion with user-defined merit function + +A new merit function can be defined, as long as it uses the same input and output +variables as the function [Merit_PROSPECT_RMSE](https://jbferet.gitlab.io/prospect/reference/Merit_PROSPECT_RMSE.html). + +Here is an example of a function using __nMRSE__ instead of __RMSE__ for the +merit function. +This means that for each spectral band, the square difference between measured and +simulated leaf optics is normalized by the value of the measured leaf optics. + +First, the merit function is + +```{r prospect merit function} +Merit_PROSPECT_nRMSE <- function(x, SpecPROSPECT, + Refl, Tran, + Input_PROSPECT, + Parms2Estimate) { + x[x < 0] <- 0 + Input_PROSPECT[Parms2Estimate] <- x + RT <- do.call('PROSPECT', c(list(SpecPROSPECT = SpecPROSPECT), Input_PROSPECT)) + fcr <- fct <- 0 + if (!is.null(Refl)) fcr <- sqrt(sum(((Refl - RT$Reflectance)**2)/Refl) / length(RT$Reflectance)) + if (!is.null(Tran)) fct <- sqrt(sum(((Tran - RT$Transmittance)**2)/Tran) / length(RT$Transmittance)) + fc <- fcr + fct + return(fc) +} +``` + +```{r prospect inverse mode merit function} +# simulate leaf optical properties +LRT <- PROSPECT(CHL = 45, CAR = 10, ANT = 0.2, + EWT = 0.02, LMA = 0.015, N = 1.8) +Parms2Estimate <- c('CHL', 'CAR', 'ANT', 'EWT', 'LMA', 'N') +# invert PROSPECT with simulated leaf optical properties +OutPROSPECT <- Invert_PROSPECT(Refl = LRT$Refl, + Tran = LRT$Tran, + Parms2Estimate = Parms2Estimate, + PROSPECT_version = 'D', + MeritFunction = 'Merit_PROSPECT_nRMSE') +``` + diff --git a/vignettes/prospect3.Rmd b/vignettes/prospect3.Rmd index e58220f..6a9e2a7 100644 --- a/vignettes/prospect3.Rmd +++ b/vignettes/prospect3.Rmd @@ -39,22 +39,7 @@ It also includes a set of measured chemical constituents: `CHL`,`CAR`,`EWT`, and ```{r get ANGERS} # Libraries required library(prospect) -library(data.table) -# repository where data are stored -gitlab_Rep <- 'https://gitlab.com/jbferet/myshareddata/raw/master/LOP/' -# download ANGERS data -dbName <- 'ANGERS' -# files available -fileName <- list('DataBioch.txt','ReflectanceData.txt','TransmittanceData.txt') -DataBioch <- fread(file.path(gitlab_Rep,dbName,fileName[[1]])) -Refl <- fread(file.path(gitlab_Rep,dbName,fileName[[2]])) -Tran <- fread(file.path(gitlab_Rep,dbName,fileName[[3]])) -# Get wavelengths corresponding to the reflectance & transmittance measurements -lambda <- Refl$wavelength -Refl$wavelength <- NULL -Tran$wavelength <- NULL -# Get the number of samples -nbSamples <- ncol(Refl) +LeafDB <- download_LeafDB(dbName = 'ANGERS') ``` # Inversion of PROSPECT-D using full spectral information @@ -67,19 +52,16 @@ Experimental leaf optics and optical constants of PROSPECT need to be adjusted b ```{r Invert PROSPECT-D Full} # Estimate all parameters for PROSPECT-D Parms2Estimate <- 'ALL' -InitValues <- data.frame(CHL = 40, CAR = 10, ANT = 0.1, BROWN = 0, - EWT = 0.01, LMA = 0.01, N = 1.5) # adjust PROSPECT optical constants & experimental leaf optics before inversion -SubData <- FitSpectralData(lambda = lambda, - Refl = Refl, - Tran = Tran) +SubData <- FitSpectralData(lambda = LeafDB$lambda, + Refl = LeafDB$Refl, + Tran = LeafDB$Tran) print('PROSPECT inversion using full spectral range') res <- Invert_PROSPECT(SpecPROSPECT = SubData$SpecPROSPECT, Refl = SubData$Refl, Tran = SubData$Tran, PROSPECT_version = 'D', - Parms2Estimate = Parms2Estimate, - InitValues = InitValues) + Parms2Estimate = Parms2Estimate) ``` ## Results: estimation of `CHL`, `CAR`, `EWT` and `LMA` @@ -108,15 +90,12 @@ CAB, CAR, ANT, EWT, and LMA can be estimated. However, no optimal spectral domai ```{r Invert PROSPECT-D Opt} # Estimate all parameters for PROSPECT-D Parms2Estimate <- c('CHL','CAR','ANT','EWT','LMA') -InitValues <- data.frame(CHL = 40, CAR = 8, ANT = 0.1, BROWN = 0, - EWT = 0.01, LMA = 0.01, N = 1.5) print('PROSPECT inversion using optimal setting') -ParmEst <- Invert_PROSPECT_OPT(lambda = lambda, - Refl = Refl, - Tran = Tran, +ParmEst <- Invert_PROSPECT_OPT(lambda = LeafDB$lambda, + Refl = LeafDB$Refl, + Tran = LeafDB$Tran, PROSPECT_version = 'D', - Parms2Estimate = Parms2Estimate, - InitValues = InitValues) + Parms2Estimate = Parms2Estimate) ``` ## Results: estimation of `CHL`, `CAR`, `EWT` and `LMA` diff --git a/vignettes/prospect4.Rmd b/vignettes/prospect4.Rmd index d299294..c91b887 100644 --- a/vignettes/prospect4.Rmd +++ b/vignettes/prospect4.Rmd @@ -45,22 +45,17 @@ The datasets can be downloaded directly with an R script as described below. ```{r get LOPEX} # Libraries required library(prospect) -library(data.table) -# repository where data are stored -gitlab_Rep <- 'https://gitlab.com/jbferet/myshareddata/raw/master/LOP/' # Datasets dbName <- list('LOPEX_DRY_CAL','LOPEX_FRESH_CAL', 'LOPEX_DRY_VAL','LOPEX_FRESH_VAL') -# files available -fileName <- list('DataBioch.txt','ReflectanceData.txt','TransmittanceData.txt') # download LOPEX data DataBioch <- Refl <- Tran <- lambda <- list() for (db in dbName){ - DataBioch[[db]] <- fread(file.path(gitlab_Rep,db,fileName[[1]])) - Refl[[db]] <- fread(file.path(gitlab_Rep,db,fileName[[2]])) - Tran[[db]] <- fread(file.path(gitlab_Rep,db,fileName[[3]])) - lambda[[db]] <- Refl[[db]]$wavelength - Refl[[db]]$wavelength <- Tran[[db]]$wavelength <- NULL + LeafDB <- download_LeafDB(dbName = db) + DataBioch[[db]] <- LeafDB$DataBioch + Refl[[db]] <- LeafDB$Refl + Tran[[db]] <- LeafDB$Tran + lambda[[db]] <- LeafDB$lambda } ``` @@ -75,16 +70,13 @@ The optimal spectral domains are defined based on the results described in [Fér Parms2Estimate <- c('EWT','PROT','CBC') EWT_mod <- PROT_mod <- CBC_mod <- list() # perform PROSPECT inversion using optimal spectral domains for EWT, PROT & CBC -InitValues <- data.frame(CHL = 40, CAR = 8, ANT = 0.1, BROWN = 0, - EWT = 0.01, CBC = 0.009, PROT = 0.001, N = 1.5) for (db in dbName){ print('PROSPECT inversion using optimal spectral domains') ParmEst <- Invert_PROSPECT_OPT(lambda = lambda[[db]], Refl = Refl[[db]], Tran = Tran[[db]], PROSPECT_version = 'PRO', - Parms2Estimate = Parms2Estimate, - InitValues = InitValues) + Parms2Estimate = Parms2Estimate) EWT_mod[[db]] <- ParmEst$EWT PROT_mod[[db]] <- ParmEst$PROT CBC_mod[[db]] <- ParmEst$CBC diff --git a/vignettes/prospect5.Rmd b/vignettes/prospect5.Rmd index 7e38bdd..9b5a9ce 100644 --- a/vignettes/prospect5.Rmd +++ b/vignettes/prospect5.Rmd @@ -47,26 +47,14 @@ This estimated N value can then be used as prior information when inverting PROS ```{r get prior estimate of N} # Libraries required library(prospect) -library(data.table) -# repository where data are stored -gitlab_Rep <- 'https://gitlab.com/jbferet/myshareddata/raw/master/LOP/' -# download ANGERS data -dbName <- 'ANGERS' -# files available -fileName <- list('DataBioch.txt','ReflectanceData.txt','TransmittanceData.txt') -DataBioch <- fread(file.path(gitlab_Rep,dbName,fileName[[1]])) -Refl <- fread(file.path(gitlab_Rep,dbName,fileName[[2]])) -Tran <- fread(file.path(gitlab_Rep,dbName,fileName[[3]])) -# Get wavelengths corresponding to the reflectance and transmittance measurements -lambda <- Refl$wavelength -Refl$wavelength <- NULL -Tran$wavelength <- NULL -# Get the number of samples -nbSamples <- ncol(Refl) +# download ANGERS dataset +LeafDB <- download_LeafDB(dbName = 'ANGERS') # Prior estimation of N using R only -Nprior_R <- Get_Nprior(lambda = lambda, Refl = Refl) +Nprior_R <- Get_Nprior(lambda = LeafDB$lambda, + Refl = LeafDB$Refl) # Prior estimation of N using T only -Nprior_T <- Get_Nprior(lambda = lambda, Tran = Tran) +Nprior_T <- Get_Nprior(lambda = LeafDB$lambda, + Tran = LeafDB$Tran) ``` ### Prior estimate of `N` based on R or T only vs. `N` estimate from iterative optimization over the full spectral domain @@ -95,8 +83,9 @@ Parms2Estimate <- 'ALL' InitValues <- data.frame(CHL = 40, CAR = 10, ANT = 0.1, BROWN = 0, EWT = 0.01, LMA = 0.01, N = 1.5) # Adjust spectral domain for SpecPROSPECT to fit leaf optical properties -SubData <- FitSpectralData(lambda = lambda, - Refl = Refl, Tran =Tran) +SubData <- FitSpectralData(lambda = LeafDB$lambda, + Refl = LeafDB$Refl, + Tran = LeafDB$Tran) print('PROSPECT inversion using full spectral range') res_Ronly <- Invert_PROSPECT(SpecPROSPECT = SubData$SpecPROSPECT, Refl = SubData$Refl, Tran = NULL, @@ -118,7 +107,7 @@ The performance of PROSPECT inversion when no prior estimation of `N` is provide Parms2Estimate <- c("CHL", "CAR", "ANT", "EWT", "LMA") print('PROSPECT inversion using full spectral range') res_R_Nprior <- res_T_Nprior <- list() -for (i in 1:nbSamples){ +for (i in 1:LeafDB$nbSamples){ print(i) InitValues <- data.frame(CHL = 40, CAR = 10, ANT = 0.1, BROWN = 0, EWT = 0.01, LMA = 0.01, N = Nprior_R$N[i])