tests/testthat/test_partGL.R
In bgctw/REddyProc: Post Processing of (Half-)Hourly Eddy-Covariance Measurements
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#+++ Unit tests for fConvertTimeToPosix functions +++
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# Author: TW
#require(testthat)
context("partGL")
if (!exists(".partGPAssociateSpecialRows")) .partGPAssociateSpecialRows <-
REddyProc:::.partGPAssociateSpecialRows
if (!exists(".binUstar")) .binUstar <- REddyProc:::.binUstar
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# 10 days from June from Example_DETha98.txt shipped with REddyProc
Example_DETha98_Filled <- getFilledExampleDETha98Data()
tzEx <- REddyProc:::getTZone(Example_DETha98_Filled$sDateTime)
test_that("example dataset starts at midngiht",{
expect_true( Example_DETha98_Filled$sDateTime[1] ==
as.POSIXct("1998-01-01 00:15:00",tz = tzEx) )
})
dsNEE <- subset(Example_DETha98_Filled
, sDateTime >= as.POSIXct("1998-06-01 00:15:00",tz = tzEx) &
sDateTime <= as.POSIXct("1998-06-09 21:45:00",tz = tzEx))
dsNEE$Temp <- dsNEE$Tair_f
dsNEE$isNight <- (dsNEE$Rg_f <= 4 & dsNEE$PotRad_NEW == 0)
dsNEE$isDay = (dsNEE$Rg_f > 4 & dsNEE$PotRad_NEW != 0)
# 8 first days of June from IT-MBo.2005.txt
# 10 days from June from Example_DETha98.txt shipped with REddyProc
.tmp.f <- function(){
#save(dsNEE, file = "tmp/dsNEE_Tharandt.RData")
load("tmp/dsNEE_Tharandt.RData") # dsNEE
dsNEE$Temp <- dsNEE$Tair_f
dsNEE$Rg_f <- dsNEE$Rg
dsNEE$isNight <- (dsNEE$Rg_f <= 4 & dsNEE$PotRad_NEW == 0)
dsNEE$isDay = (dsNEE$Rg_f > 4 & dsNEE$PotRad_NEW != 0)
}
.tmp.f2 <- function(){
# else stop in partitionNEEGL, and grap ds:
attr(ds$sDateTime, "tzone") <- "UTC"
dsJune <- ds
dsNEE <- dsJune[ dsJune$sDateTime >= as.POSIXct("1998-06-01", tz = "UTC") &
dsJune$sDateTime < as.POSIXct("1998-06-10", tz = "UTC")
, c("sDateTime","NEE_f","NEE_fqc","NEE_fsd","Tair_f","Tair_fqc","VPD_f"
,"VPD_fqc","Rg_f","PotRad_NEW")]
save(dsNEE, file = "tmp/dsNEE_Tharandt.RData")
dput(dsNEE)
}
# regression result from dput in
resLRCEx1 <- structure(list(
iWindow = 1:4, dayStart = c(1L, 3L, 5L, 7L), dayEnd = c(4L, 6L, 8L, 9L)
, iRecStart = c(1L, 97L, 193L, 289L), iRecEnd = c(192L, 288L, 384L, 428L)
, iCentralRec = c(97L, 193L, 289L, 385L), nValidRec = c(98L,98L, 75L, 56L)
, iMeanRec = c(97L, 194L, 283L, 385L)
, convergence = c(0L,0L, 0L, 0L), parms_out_range = c(0L, 0L, 0L, 0L)
, k = c(
0.15556543911127, 0.100336071150363, 0.0292272454449355, 0.110909241852741)
, beta = c(
34.2363908612952, 28.7053014210465, 11.965007963525, 35.6061046130658)
, alpha = c(
0.0408287988682114, 0.0317164138158062, 0.12184364495971, 0.0565989757582912)
, RRef = c(
1.30588274435898, 2.30342378833466, 4.39970416578445, 3.02347421699605)
, E0 = structure(
c(62.633343296311, 62.633343296311, 62.633343296311, 62.633343296311)
, .Dim = c(4L, 1L), .Dimnames = list(NULL, NULL))
, k_sd = c(
0.0162347703672834, 0.00819228863397153, 0.0036484758449877, 0.0143126637303779)
, beta_sd = c(
2.6633120785356, 2.53043346467805, 0.480085981083196, 3.86498333963868)
, alpha_sd = c(
0.00282359981111383, 0.00245807510249875, 0.0189581631026536, 0.00624526974304941)
, RRef_sd = c(
0.243651749533536, 0.200219643557999, 0.379555581976541, 0.408032180867498)
, E0_sd = structure(
c(4.25106879296635, 4.22511397480583, 4.22511362599932, 4.46593783063789)
, .Dim = c(4L, 1L), .Dimnames = list(NULL, NULL))
, GPP2000 = c(
24.1225750077742, 19.7622684102263, 11.4050231061514, 27.0862112802389)
, isValid = c(TRUE, TRUE, TRUE, TRUE)
, resOpt = list(structure(list(thetaOpt = structure(
c(0.15556543911127, 34.2363908612952, 0.0408287988682114
, 1.30588274435898, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, iOpt = 1:5, thetaInitialGuess = structure(c(
0.05, 24.3542, 0.1, 5.17329116845493, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, covParms = structure(c(
0.000263567768878423, 0.0299004485274175, -2.01722085100335e-05
, -0.00170017321010028, 0, 0.0299004485274174, 7.09323122767363
, -0.00442973195793744, -0.184157427431281, 0, -2.01722085100335e-05
, -0.00442973195793744, 7.97271589332207e-06, 0.000574595628371077, 0
, -0.00170017321010028, -0.18415742743128, 0.000574595628371077
, 0.0593661750507529, 0, 0, 0, 0, 0, 18.0715858825324)
, .Dim = c(5L, 5L), .Dimnames = list(
structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))))
, convergence = 0L)
, .Names = c(
"thetaOpt", "iOpt", "thetaInitialGuess", "covParms", "convergence"))
, structure(list(thetaOpt = structure(
c(0.100336071150363, 28.7053014210465, 0.0317164138158062
, 2.30342378833466, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, iOpt = 1:5
, thetaInitialGuess = structure(
c(0.05, 24.8403, 0.1, 5.14770049230269, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, covParms = structure(
c(6.71135930622992e-05, 0.0178801034106948, -1.74280275277554e-05
, -0.00134394852514135, 0, 0.0178801034106948, 6.40309351916256
, -0.00493560933704443, -0.271156219332448, 0, -1.74280275277554e-05
, -0.00493560933704443, 6.04213320952425e-06, 0.000418391213333105
, 0, -0.00134394852514135, -0.271156219332449, 0.000418391213333105
, 0.040087905666492, 0, 0, 0, 0, 0, 17.8515881000995)
, .Dim = c(5L, 5L), .Dimnames = list(
structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))))
, convergence = 0L)
, .Names = c(
"thetaOpt", "iOpt", "thetaInitialGuess", "covParms", "convergence"))
, structure(list(thetaOpt = structure(
c(0.0292272454449355, 11.965007963525, 0.12184364495971
, 4.39970416578445, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, iOpt = 1:5
, thetaInitialGuess = structure(
c(0.05, 24.1194, 0.1, 5.14770049230269, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, covParms = structure(
c(1.33113759914587e-05, 0.000439566354933103, -5.94748188903169e-05
, -0.00111894361237909, 0, 0.000439566354933103, 0.230482549232615
, -0.000116871653668728, 0.053260722726946, 0, -5.94748188903169e-05
, -0.000116871653668727, 0.000359411948226816, 0.00657062463616557
, 0, -0.00111894361237909, 0.0532607227269459, 0.00657062463616558
, 0.144062439809551, 0, 0, 0, 0, 0, 17.8515851526051)
, .Dim = c(5L, 5L), .Dimnames = list(
structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))))
, convergence = 0L)
, .Names = c(
"thetaOpt", "iOpt", "thetaInitialGuess", "covParms", "convergence"))
, structure(list(thetaOpt = structure(
c(0.110909241852741, 35.6061046130658, 0.0565989757582912
, 3.02347421699605, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, iOpt = 1:5
, thetaInitialGuess = structure(
c(0.05, 21.665, 0.1, 5.21501814931109, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, covParms = structure(c(
0.000204852343058874, 0.0495123146513288, -6.84478020597351e-05
, -0.00386015942484132, 0, 0.0495123146513288, 14.9380962156846
, -0.0178713249608098, -0.786540464781837, 0, -6.8447802059735e-05
, -0.0178713249608097, 3.90033941634485e-05, 0.00226964418347153
, 0, -0.00386015942484131, -0.786540464781834, 0.00226964418347153
, 0.166490260623486, 0, 0, 0, 0, 0, 19.9446007071226)
, .Dim = c(5L, 5L), .Dimnames = list(structure(
c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, structure(
c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))))
, convergence = 0L)
, .Names = c(
"thetaOpt", "iOpt", "thetaInitialGuess", "covParms", "convergence")))
, E0_bootstrap_sd = c(
4.25106879296635, 4.22511397480583, 4.22511362599932, 4.46593783063789)
, RRef_night = c(
5.17329116845493, 5.14770049230269, 5.14770049230269, 5.21501814931109))
, .Names = c(
"iWindow", "dayStart", "dayEnd", "iRecStart", "iRecEnd", "iCentralRec"
, "nValidRec", "iMeanRec", "convergence", "parms_out_range", "k", "beta"
, "alpha", "RRef", "E0", "k_sd", "beta_sd", "alpha_sd", "RRef_sd", "E0_sd"
, "GPP2000", "isValid", "resOpt", "E0_bootstrap_sd", "RRef_night")
, row.names = c(NA, -4L)
, class = c("tbl_df", "tbl", "data.frame")
)
resLRCEx1Nonrectangular <- structure(list(
iWindow = 1:4, dayStart = c(1L, 3L, 5L, 7L), dayEnd = c(4L, 6L, 8L, 9L)
, iRecStart = c(1L, 97L, 193L, 289L), iRecEnd = c(192L, 288L, 384L, 428L)
, iCentralRec = c(97L, 193L, 289L, 385L), nValidRec = c(98L, 98L, 75L, 56L)
, iMeanRec = c(97L, 194L, 283L, 385L), convergence = c(1005L, 1005L, 1005L, 0L)
, E0 = structure(c(
62.633343296311, 62.633343296311, 62.633343296311, 62.633343296311)
, .Dim = c(4L, 1L), .Dimnames = list( NULL, NULL))
, E0_sd = structure(
c(4.25106879296635, 4.22511397480583, 4.22511362599932, 4.46593783063789)
, .Dim = c(4L, 1L), .Dimnames = list( NULL, NULL))
, RRefNight = c(5.17329116845493, 5.14770049230269, 5.14770049230269, NA)
, parms_out_range = c(NA, NA, NA, 0L), k = c(NA, NA, NA, 0.15301818456226)
, beta = c(NA, NA, NA, 49.3407957035085 )
, alpha = c(NA, NA, NA, 0.0444237102808321)
, RRef = c(NA, NA, NA, 2.35368171175251)
, logitconv = c(NA, NA, NA, -3.10163338239329 )
, k_sd = c(NA, NA, NA, 0.0197088986436212)
, beta_sd = c(NA, NA, NA, 7.58853578454936)
, alpha_sd = c(NA, NA, NA, 0.00432079403235484 )
, RRef_sd = c(NA, NA, NA, 0.350880271587216)
, logitconv_sd = c(NA, NA, NA, 1.92633463362266)
, GPP2000 = c(NA, NA, NA, 32.0432128516724 )
, isValid = c(NA, NA, NA, TRUE)
, resOpt = list(NULL, NULL, NULL, structure(
list(thetaOpt = structure(c(
0.15301818456226, 49.3407957035085, 0.0444237102808321, 2.35368171175251
, 62.633343296311, -3.10163338239329 )
, .Names = c("k", "beta", "alpha", "RRef", "E0", "logitconv" ))
, iOpt = c(1, 2, 3, 4, 6, 5)
, thetaInitialGuess = structure(c(
0.05, 21.665, 0.1, 5.21501814931109, 62.633343296311, 1.09861228866811 )
, .Names = c("k", "beta", "alpha", "RRef", "E0", "logitconv" ))
, covParms = structure(c(
0.000388440685744535, 0.13746635747882, -6.51906051810543e-05
, -0.00449519721692684, 0, 0.00207560732637497, 0.13746635747882
, 57.5858753533862, -0.0238832979670161, -1.36443384766771, 0
, -1.52788892519458, -6.51906051810543e-05, -0.023883297967016
, 1.86692610700332e-05, 0.00131289852541519, 0, -0.00181474371141437
, -0.00449519721692684, -1.36443384766771, 0.00131289852541519
, 0.123116964989119, 0, -0.0864363703508266, 0, 0, 0, 0, 19.9446007071226
, 0, 0.00207560732637499, -1.52788892519458, -0.00181474371141437
, -0.0864363703508268, 0, 3.71076512069413 )
, .Dim = c(6L, 6L), .Dimnames = list(structure(c(
"k", "beta", "alpha", "RRef", "E0", "logitconv")
, .Names = c("k", "beta", "alpha", "RRef", "E0", "logitconf"))
, structure(c(
"k", "beta", "alpha", "RRef", "E0", "logitconv")
, .Names = c("k", "beta", "alpha", "RRef", "E0", "logitconf"))))
, convergence = 0L)
, .Names = c(
"thetaOpt", "iOpt", "thetaInitialGuess", "covParms", "convergence")))
, E0_bootstrap_sd = c(
4.25106879296635, 4.22511397480583, 4.22511362599932, 4.46593783063789)
, RRef_night = c(
5.17329116845493, 5.14770049230269, 5.14770049230269, 5.21501814931109))
, .Names = c(
"iWindow", "dayStart", "dayEnd", "iRecStart", "iRecEnd", "iCentralRec"
, "nValidRec", "iMeanRec", "convergence", "E0", "E0_sd", "RRefNight"
, "parms_out_range", "k", "beta", "alpha", "RRef", "logitconv", "k_sd"
, "beta_sd", "alpha_sd", "RRef_sd", "logitconv_sd", "GPP2000", "isValid"
, "resOpt", "E0_bootstrap_sd", "RRef_night")
, row.names = c(NA, -4L)
, class = c("tbl_df", "tbl", "data.frame")
)
test_that("replaceMissingSdByPercentage",{
n <- 100
x <- rnorm(n,10)
sdX <- sdX0 <- pmax(0.7, rnorm(n, 10*0.1 ))
iMissing <- sample(n, 20)
sdX[iMissing] <- NA
sdXf <- REddyProc:::replaceMissingSdByPercentage(sdX,x, 0.2, 1.7)
#plot( sdXf ~ x, col = "blue" );points( sdX ~ x )
expect_equal( sdXf[-iMissing], sdX[-iMissing])
expect_true( all(is.finite(sdXf)))
expect_true( all(sdXf[iMissing] >= 1.7))
expect_true( all(sdXf[iMissing] >= 0.2*x[iMissing]))
#
sdX <- sdX0
iMissing <- TRUE
sdX[iMissing] <- NA
sdXf <- REddyProc:::replaceMissingSdByPercentage(sdX,x, 0.2, 1.7)
#plot( sdXf ~ x, col = "blue" );points( sdX ~ x )
expect_true( all(is.finite(sdXf)))
expect_true( all(sdXf[iMissing] >= 1.7))
expect_true( all(sdXf[iMissing] >= 0.2*x[iMissing]))
#
sdX <- sdX0
iMissing <- sample(n, 20)
sdX[iMissing] <- NA
sdXf <- REddyProc:::replaceMissingSdByPercentage(sdX,x, NA, 1.7)
expect_true( all(sdXf[iMissing] == 1.7))
#plot( sdXf ~ x, col = "blue" );points( sdX ~ x )
sdXf <- REddyProc:::replaceMissingSdByPercentage(sdX,x, 0.2, NA)
expect_true( all(sdXf[iMissing] == 0.2*x[iMissing]))
sdXf <- REddyProc:::replaceMissingSdByPercentage(sdX,x, NA, NA)
expect_true( all(is.na(sdXf[iMissing])))
})
test_that("estimating temperature sensitivity oneWindow are in accepted range",{
skip_if_not_installed("mlegp")
dss <- dsNEE[ dsNEE$Rg_f <= 0 & dsNEE$PotRad_NEW <= 0
& as.POSIXlt(dsNEE$sDateTime)$mday %in% 1:8, ]
dss <- dss[ order(dss$Temp), ]
dss$NEE <- dss$NEE_f
resE0 <- REddyProc:::partGLEstimateTempSensInBoundsE0Only(
dss$NEE_f, dss$Temp + 273.15)
expect_true( resE0$E0 >= 50 && resE0$E0 < 400 )
medianResp <- median(dss$NEE_f,na.rm = TRUE)
expect_true( abs(resE0$RRefFit - medianResp)/medianResp < 0.2 )
E0Win <- as.data.frame(resE0)
res <- REddyProc:::partGLFitNightRespRefOneWindow(
dss, data.frame(iWindow = 1L), E0Win = E0Win)
RRef <- res[1]
expect_true( RRef >= 0)
.tmp.plot <- function(){
plot( NEE_f ~ Temp, dss) # FP_VARnight negative?
lines( fLloydTaylor(RRef, resE0$E0, dss$Temp+273.15, TRef = 273.15+15) ~
dss$Temp)
}
})
test_that("estimating temperature sensitivity on record with some freezing temperatures",{
skip_if_not_installed("mlegp")
dss <- dsNEE[ dsNEE$Rg_f <= 0 & dsNEE$PotRad_NEW <= 0 &
as.POSIXlt(dsNEE$sDateTime)$mday %in% 1:8, ]
dss$NEE <- dss$NEE_f
dss <- dss[ order(dss$Temp), ]
# only the last 13 records will have temperature above -1degC
dss$Temp <- dss$Temp - dss$Temp[ nrow(dss) - 14L ] - 1
isValid <- REddyProc:::isValidNightRecord(dss)
expect_true(sum(isValid) >= 13 &&
sum(REddyProc:::isValidNightRecord(dss)) < nrow(dss) )
#
resE0 <- REddyProc:::partGLEstimateTempSensInBoundsE0Only(
dss$NEE_f[isValid], dss$Temp[isValid] + 273.15)
expect_true( is.na(resE0$E0) )
.tmp.f <- function(){
expect_true( resE0$E0 >= 50 && resE0$E0 <= 400 )
medianResp <- median(dss$NEE_f[isValid],na.rm = TRUE)
expect_true( abs(resE0$RRefFit - medianResp)/medianResp < 0.2 )
E0Win <- as.data.frame(resE0)
res <- REddyProc:::partGLFitNightRespRefOneWindow(
dss, data.frame(iWindow = 1L), E0Win = E0Win)
RRef <- res[[2]]$RRef[1]
expect_true( RRef >= 0 )
}
.tmp.plot <- function(){
plot( NEE_f ~ Temp, dss) # FP_VARnight negative?
points( NEE_f ~ Temp, dss[isValid,], col = "red")
lines( fLloydTaylor(RRef, resE0$E0, dss$Temp+273.15, TRef = 273.15+15) ~
dss$Temp)#
}
})
test_that("applyWindows",{
nRec <- nrow(dsNEE)
nRecInDay <- 10L
ds <- within(dsNEE, {
iRec <- 1:nRec
# specifying the day for each record assuming equidistand records
iDayOfRec <- ((c(1:nRec) - 1L) %/% nRecInDay) + 1L
})
fReportTime <- function(dss, winInfo, prevRes){
nRecS <- nrow(dss)
list( res2 = data.frame(
startRec = dss$iRec[1]
,endRec = dss$iRec[nRecS]
,startDay = dss$iDayOfRec[1]
,endDay = dss$iDayOfRec[nRecS]
)
,sumRes = prevRes$sumRes + 1L
)
}
prevRes <- list(sumRes = 0)
fReportTime(ds,0,prevRes)
# larger than reference window of 4 days
resApply <- REddyProc:::applyWindows(
ds, fReportTime, prevRes, winSizeInDays = 6L, nRecInDay = nRecInDay )
res <- cbind( resApply$winInfo, tmp <- do.call(
rbind, lapply(resApply$resFUN, "[[", 1L)))
nRecRes <- nrow(res)
expect_equal( res$dayStart, res$startDay )
expect_equal( res$dayEnd, res$endDay )
expect_equal( res$iRecStart, res$startRec )
expect_equal( res$iRecEnd, res$endRec )
#
expect_true( all(diff(res$startDay[-1]) == 2L)) # shifted starting day
# day boundary before startRec
expect_true( all((ds$iDayOfRec[res$startRec[-1]] -
ds$iDayOfRec[res$startRec[-1] - 1]) == 1L))
# day boundary after endRec
expect_true( all((ds$iDayOfRec[res$startRec[-nRecRes]] + 1 -
ds$iDayOfRec[res$startRec[-nRecRes]]) == 1L))
# prevRes accumulated
expect_equal( resApply$resFUN[[nRecRes]]$sumRes, nrow(res) )
})
test_that("simplifyApplyWindows",{
nRec <- nrow(dsNEE)
nRecInDay <- 10L
ds <- within(dsNEE, {
iRec <- 1:nRec
# specifying the day for each record assuming equidistand records
iDayOfRec <- ((c(1:nRec) - 1L) %/% nRecInDay) + 1L
})
fReportTimeSimple <- function(dss, winInfo, prevRes = list()){
nRecS <- nrow(dss)
c(
startRec = dss$iRec[1]
,endRec = dss$iRec[nRecS]
,startDay = dss$iDayOfRec[1]
,endDay = dss$iDayOfRec[nRecS]
)
}
fReportTimeSimple(ds,0)
# larger than reference window of 4 days
resApply <- REddyProc:::applyWindows(
ds, fReportTimeSimple, winSizeInDays = 6L, nRecInDay = nRecInDay )
res <- REddyProc:::simplifyApplyWindows(resApply)
nRecRes <- nrow(res)
expect_equal( res$dayStart, res$startDay )
expect_equal( res$dayEnd, res$endDay )
expect_equal( res$iRecStart, res$startRec )
expect_equal( res$iRecEnd, res$endRec )
#
# repeat with function returning a single-row data.frame
fReportTimeSimpleDs <- function(dss, winInfo, prevRes = list()){
nRecS <- nrow(dss)
data.frame(
startRec = dss$iRec[1]
,endRec = dss$iRec[nRecS]
,startDay = dss$iDayOfRec[1]
,endDay = dss$iDayOfRec[nRecS]
)
}
fReportTimeSimpleDs(ds,0)
# larger than reference window of 4 days
resApply <- REddyProc:::applyWindows(
ds, fReportTimeSimpleDs, winSizeInDays = 6L, nRecInDay = nRecInDay )
resDs <- REddyProc:::simplifyApplyWindows(resApply)
expect_equal(res, resDs )
})
test_that("estimating temperature sensitivity windows outputs are in accepted range",{
skip_if_not_installed("mlegp")
dss <- dsNEE[ dsNEE$Rg_f <= 0 & dsNEE$PotRad_NEW <= 0 &
as.POSIXlt(dsNEE$sDateTime)$mday %in% 1:12, ]
dss$NEE <- dss$NEE_f
dss <- dss[ order(dss$Temp), ]
res <- REddyProc:::simplifyApplyWindows(
REddyProc:::applyWindows(
dss, REddyProc:::partGLFitNightTempSensOneWindow
, prevRes = data.frame(E0 = NA)
, winSizeInDays = 12L
#,controlGLPart = controlGLPart
))
#res <- partGLEstimateTempSensInBounds(dss$NEE_f, dss$Temp+273.15)
expect_true( res$E0 >= 50 && res$E0 <= 400 )
expect_true( res$RRefFit > 0 )
})
test_that("partGLFitLRCWindows outputs are in accepted range",{
skip_if_not_installed("mlegp")
ds <- partGLExtractStandardData(dsNEE)
#
#yday <- as.POSIXlt(dsNEE$sDateTime)$yday
dsTempSens <- dsTempSens0 <- REddyProc:::partGLFitNightTimeTRespSens(
ds
, nRecInDay = 48L
, controlGLPart = partGLControl()
)
lrcFitter <- RectangularLRCFitter()
#lrcFitter <- NonrectangularLRCFitter()
resFits <- REddyProc:::partGLFitLRCWindows(
ds, nRecInDay = 48L, dsTempSens = dsTempSens
, controlGLPart = partGLControl(nBootUncertainty = 10L)
, lrcFitter = lrcFitter)
expect_true( var(resFits$k) > 0)
expect_true( all( sapply(resFits$resOpt, length) != 0 ))
.tmp.f <- function(){
# in order to replicate, use nBoot = 0L
resParms0 <- REddyProc:::partGLFitLRCWindows(
ds, nRecInDay = 48L, dsTempSens = dsTempSens
, controlGLPart = partGLControl(nBootUncertainty = 0L)
, lrcFitter = lrcFitter)
iTooMany <- which(resParms0$iRecStart >= nrow(ds) - 2*48)
if (length(iTooMany)) {
resParms0 <- slice(resParms0,-iTooMany)
}
# regression result for resLRCEx1
# resLRCEx1 <- resParms0
# resLRCEx1Nonrectangular <- resParms0
dput(resParms0)
}
# check the conditions of Lasslop10 Table A1
expect_true( all(resFits$E0 >= 50 & resFits$E0 <= 400) )
expect_true( all(resFits$RRef > 0) )
expect_true( all(resFits$alpha >= 0 ) )
# first value may be greater, due to previous estimate
expect_true( all(resFits$alpha[-1] < 0.22) )
expect_true( all(resFits$beta >= 0 & resFits$beta < 250) )
expect_true( all(ifelse(resFits$beta > 100, resFits$beta_sd < resFits$beta, TRUE) ))
expect_true( all(resFits$k >= 0) )
expect_true( !all(is.na(resFits$RRef_sd)))
expect_true( all(resFits$iMeanRec < nrow(ds)) )
expect_true( all(resFits$iCentralRec < nrow(ds)) )
})
test_that("partGLFitLRCWindows error on low uncertainty",{
skip_if_not_installed("mlegp")
# see LightResponseCurveFitter_optimLRCOnAdjustedPrior
ds <- partGLExtractStandardData(dsNEE)
ds$sdNEE[100:200] <- 0
dsTempSens <- dsTempSens0 <- REddyProc:::partGLFitNightTimeTRespSens(
ds
, nRecInDay = 48L
, controlGLPart = partGLControl()
)
lrcFitter <- RectangularLRCFitter()
#lrcFitter <- NonrectangularLRCFitter()
if (exists("resFits")) rm(resFits)
expect_error(
resFits <- REddyProc:::partGLFitLRCWindows(
ds, nRecInDay = 48L, dsTempSens = dsTempSens
, controlGLPart = partGLControl(nBootUncertainty = 10L)
, lrcFitter = lrcFitter)
,"zeros in uncertainty"
)
expect_true(!exists("resFits"))
})
test_that(".partGPAssociateSpecialRows correct next lines",{
skip_if_not_installed("mlegp")
expect_error(
res <- REddyProc:::.partGPAssociateSpecialRows(integer(0),9)
)
res <- REddyProc:::.partGPAssociateSpecialRows(c(3,6,7,9),12)
expect_true( all(res$wBefore + res$wAfter == 1))
# special rows
expect_equal( c(iRec = 3,iSpecialBefore = 1L,iSpecialAfter = 1L
, iBefore = 3, iAfter = 3, wBefore = 0.5, wAfter = 0.5)
, unlist(res[3,]))
expect_equal( c(iRec = 6,iSpecialBefore = 2L,iSpecialAfter = 2L
, iBefore = 6, iAfter = 6, wBefore = 0.5, wAfter = 0.5)
, unlist(res[6,]))
# first rows and last rows
expect_equal( c(iRec = 1,iSpecialBefore = 1L,iSpecialAfter = 1L
, iBefore = 3, iAfter = 3, wBefore = 0.5, wAfter = 0.5)
, unlist(res[1,]))
expect_equal( c(iRec = 12,iSpecialBefore = 4L,iSpecialAfter = 4L
, iBefore = 9, iAfter = 9, wBefore = 0.5, wAfter = 0.5)
, unlist(res[12,]))
# weights after and before special row 6
expect_equal( rep(3,2), unlist(res[4:5,"iBefore"]))
expect_equal( rep(6,2), unlist(res[4:5,"iAfter"]))
expect_equal( rep(1,2), unlist(res[4:5,"iSpecialBefore"]))
expect_equal( rep(2,2), unlist(res[4:5,"iSpecialAfter"]))
expect_equal( (2:1)/3, unlist(res[4:5,"wBefore"]))
expect_equal( (1:2)/3, unlist(res[4:5,"wAfter"]))
# test last row is special
res <- REddyProc:::.partGPAssociateSpecialRows(c(3,6,7,9),9)
})
testGLInterpolateFluxes <- function(lrcFitter, resEx ){
tmp <- REddyProc:::partGLInterpolateFluxes(
dsNEE$Rg_f, dsNEE$VPD_f, dsNEE$Temp, resEx, lrcFitter = lrcFitter)
expect_equal( nrow(dsNEE), nrow(tmp) )
.tmp.plot <- function(){
tmp$time <- dsNEE$sDateTime
plot( Reco ~ time, tmp)
plot( GPP ~ time, tmp)
}
expect_true( all(c("GPP","Reco") %in% names(tmp) ))
tmp <- REddyProc:::partGLInterpolateFluxes(
dsNEE$Rg_f, dsNEE$VPD_f, dsNEE$Temp, resEx
, controlGLPart = partGLControl(isSdPredComputed = TRUE)
, lrcFitter = lrcFitter)
expect_equal( nrow(dsNEE), nrow(tmp) )
expect_true( all(c("GPP","Reco") %in% names(tmp) ))
expect_true( all(c("sdGPP","sdReco") %in% names(tmp) ))
}
test_that("partGLInterpolateFluxes runs with rectangular LRCFitter",{
lrcFitter = RectangularLRCFitter(); resEx <- resLRCEx1
testGLInterpolateFluxes( lrcFitter, resEx)
})
test_that("partGLInterpolateFluxes runs with rectangular LRCFitter",{
lrcFitter = NonrectangularLRCFitter(); resEx <- resLRCEx1Nonrectangular
testGLInterpolateFluxes( lrcFitter, resEx)
})
#resLRCEx1
test_that("partitionNEEGL",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
resEx <- resLRCEx1
#DoY.V.n <- as.POSIXlt(dsNEE1$sDateTime)$yday + 1L
#Hour.V.n <- as.POSIXlt(dsNEE1$sDateTime)$hour + as.POSIXlt(dsNEE1$sDateTime)$min/60
#dsNEE1$PotRad_NEW <- fCalcPotRadiation(DoY.V.n, Hour.V.n, LatDeg = 45.0, LongDeg = 1, TimeZoneHour = 0 )
tmp <- partitionNEEGL( dsNEE1,RadVar = 'Rg_f')
#tmp <- partitionNEEGL( dsNEE1,RadVar = 'Rg_f', isVerbose = FALSE) #no messages
#tmp <- partitionNEEGL( dsNEE1,RadVar = 'Rg_f', controlGLPart = partGLControl(nBootUncertainty = 0L, isAssociateParmsToMeanOfValids = FALSE))
expect_equal( nrow(dsNEE1), nrow(tmp) )
#tmp[ is.finite(tmp$FP_beta), ] # note FP_dRecPar is not zero, because iCentralRec != iMeanRec
#iPar <- which(is.finite(tmp$FP_E0))
#plot( tmp$FP_beta[iPar] ~ dsNEE1$sDateTime[iPar],type = "l", ylim = range(c(tmp$FP_beta[iPar],tmp$FP_GPP2000[iPar]))); lines( tmp$FP_GPP2000[iPar] ~ dsNEE1$sDateTime[iPar], col = "red")
#
# now test with different suffix: u50
dsNEE2 <- dsNEE
names(dsNEE2)[ match(c("NEE_f", "NEE_fqc", "NEE_fsd"),names(dsNEE2))] <-
c("NEE_u50_f", "NEE_u50_fqc", "NEE_u50_fsd")
tmp <- partitionNEEGL(
dsNEE2, RadVar = 'Rg_f', suffix = "u50"
, controlGLPart = partGLControl(
nBootUncertainty = 0L, isAssociateParmsToMeanOfValids = FALSE) )
expect_equal( nrow(dsNEE1), nrow(tmp) )
expect_true( all(is.finite(tmp$GPP_DT_u50)))
expect_true( all(tmp$GPP_DT_u50 >= 0))
expect_true( all(tmp$GPP_DT_u50 < 250))
expect_true( all(tmp$Reco_DT_u50 < 10))
expect_true( all(tmp$Reco_DT_u50 > 0))
expect_true( all(tmp$Reco_DT_u50_SD > 0))
expect_true( all(tmp$GPP_DT_u50_SD >= 0))
expect_true( all(abs(diff(tmp$Reco_DT_u50)) < 0.6)) #smooth
# reporting good values at central records
# tmp[resEx$iCentralRec,]
iRowsOpt <- which(is.finite(tmp$FP_E0))
expect_true( length(iRowsOpt) == nrow(resEx) )
#expect_true( all((tmp$FP_alpha[resEx$iCentralRec] - resEx$alpha)[
#resEx$parms_out_range == 0L] < 1e-2) )
.tmp.plot <- function(){
tmp$time <- dsNEE1$sDateTime
plot( Reco_DT_u50 ~ time, tmp)
#plot( diff(Reco_DT_u50) ~ time[-1], tmp)
plot( GPP_DT_u50 ~ time, tmp)
plot( FP_RRef_Night ~ time, tmp)
plot( FP_RRef ~ FP_RRef_Night, tmp)
}
})
test_that("partitionNEEGL with Lasslop options",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
#ds <- data.frame( NEE = dsNEE1$NEE_f, sdNEE = dsNEE1$NEE_fsd, Rg = dsNEE1$Rg_f, VPD = dsNEE1$VPD_f, Temp = dsNEE1$Temp, isDay = dsNEE1$isDay, isNight = dsNEE$isNight )
resEx <- resLRCEx1
dsNEE2 <- dsNEE
partGLControl()
names(dsNEE2)[ match(c("NEE_f", "NEE_fqc", "NEE_fsd"),names(dsNEE2))] <-
c("NEE_u50_f", "NEE_u50_fqc", "NEE_u50_fsd")
dsNEE2$NEE_u50_fsd[24] <- NA
tmp <- partitionNEEGL( dsNEE2, RadVar = 'Rg_f', suffix = "u50"
, controlGLPart = partGLControlLasslopCompatible() )
tmp$sDateTime <- dsNEE2$sDateTime
tmp$iRec <- 1:nrow(tmp)
#subset(tmp, is.finite(FP_errorcode))
expect_equal( nrow(dsNEE1), nrow(tmp) )
expect_true( all(is.finite(tmp$GPP_DT_u50)))
expect_true( all(tmp$GPP_DT_u50 >= 0))
expect_true( all(tmp$GPP_DT_u50 < 250))
expect_true( all(tmp$Reco_DT_u50 < 10))
expect_true( all(tmp$Reco_DT_u50 > 0))
expect_true( all(tmp$Reco_DT_u50_SD > 0))
expect_true( all(tmp$GPP_DT_u50_SD >= 0))
expect_true( all(abs(diff(tmp$Reco_DT_u50)) < 0.6)) #smooth
# reporting good values at central records
# tmp[resEx$iCentralRec,]
#expect_true( sum( is.finite(tmp$FP_alpha) ) == nrow(resEx) )
# options yield large differences, hence do not compare to resEx1
#expect_true( all(abs(tmp$FP_alpha[resEx$iCentralRec] - resEx$alpha)[
#resEx$parms_out_range == 0L & is.finite(tmp$FP_alpha[resEx$iCentralRec])] < 0.05) )
#expect_true( all((is.na(tmp$FP_alpha[resLRCEx1$iFirstRec] - resLRCEx1$a)[
#resLRCEx1$parms_out_range != 0L])) )
expect_true( length(is.finite(tmp$FP_RRef_Night)) > 0 )
.tmp.plot <- function(){
tmp$time <- dsNEE1$sDateTime
plot( Reco_DT_u50 ~ time, tmp)
#plot( diff(Reco_DT_u50) ~ time[-1], tmp)
plot( GPP_DT_u50 ~ time, tmp)
plot( FP_RRef_Night ~ time, tmp)
plot( FP_RRef ~ FP_RRef_Night, tmp)
}
})
isTimeInTestPeriod <- function(sDateTime){
(sDateTime >= as.POSIXct("1998-06-03 00:00:00",tz = tzEx)) &
(sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx) |
sDateTime >= as.POSIXct("1998-06-05 22:00:00",tz = tzEx))
}
.tmp.f <- function(){
isTimeInTestPeriodNoTemp <- function(sDateTime){
# strange: 4 more valid NEE yiels in temperature estimate failing TODO inspect
(sDateTime >= as.POSIXct("1998-06-03 00:00:00",tz = tzEx)) &
(sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx) |
sDateTime >= as.POSIXct("1998-06-06 00:00:00",tz = tzEx))
}
dsNEE2 <- dsNEE
dsNEE2$NEE_fqc[ isTimeInTestPeriodNoTemp(dsNEE$sDateTime) ] <- 2L
#table(dsNEE1$NEE_fqc)
#plot( NEE_f ~ sDateTime, dsNEE1 )
ds2 <- partGLExtractStandardData(dsNEE2)
dsTempSens2 <- REddyProc:::partGLFitNightTimeTRespSens( ds
, nRecInDay = 48L
, controlGLPart = partGLControl()
)
}
test_that("partitionNEEGL sparse data",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
#flag all data except one day bad in order to associate the same
#data with several windows
dsNEE1$NEE_fqc[ isTimeInTestPeriod(dsNEE$sDateTime) ] <- 2L
#table(dsNEE1$NEE_fqc)
#plot( NEE_f ~ sDateTime, dsNEE1 )
ds <- partGLExtractStandardData(dsNEE1)
#
dsTempSens <- REddyProc:::partGLFitNightTimeTRespSens(
ds
, nRecInDay = 48L
, controlGLPart = partGLControl()
)
resLRC <- REddyProc:::partGLFitLRCWindows(
ds, dsTempSens = dsTempSens, lrcFitter = RectangularLRCFitter() )
expect_true( resLRC$iMeanRec[2] == resLRC$iMeanRec[3])
#
tmp <- partitionNEEGL( dsNEE1)
expect_equal( nrow(dsNEE1), nrow(tmp) )
#
expect_true( all(is.finite(tmp$GPP_DT)))
expect_true( all(tmp$GPP_DT >= 0))
expect_true( all(tmp$GPP_DT < 250))
expect_true( all(tmp$Reco_DT < 10))
expect_true( all(tmp$Reco_DT > 0))
expect_true( all(tmp$Reco_DT_SD > 0))
expect_true( all(tmp$GPP_DT_SD >= 0))
#TODO expect_true( all(abs(diff(tmp$Reco_DT)) < 0.6)) #smooth
# reporting good values at first row
expect_true( sum( is.finite(tmp$FP_alpha) ) == sum(is.finite(resLRC$alpha)) )
expect_true( all((is.na(tmp$FP_alpha[resLRC$iCentralRec] - resLRC$alpha)[
resLRC$parms_out_range != 0L & is.finite(tmp$FP_alpha[resLRC$iCentralRec])])) )
.tmp.plot <- function(){
tmp$time <- dsNEE1$sDateTime
plot( Reco_DT ~ time, tmp)
#plot( diff(Reco_DT_u50) ~ time[-1], tmp)
plot( GPP_DT ~ time, tmp)
}
})
test_that("partitionNEEGL isNeglectVPDEffect",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
#flag all VPD data except one day as bad to associate the same data with several windows
dsNEE1$VPD_fqc[ (dsNEE$sDateTime >= as.POSIXct("1998-06-03 00:00:00",tz = tzEx)) &
(dsNEE$sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx) |
dsNEE$sDateTime >= as.POSIXct("1998-06-06 00:00:00",tz = tzEx)) ] <- 2L
#plot(VPD_fqc ~ sDateTime, dsNEE1)
ctrl <- partGLControl(isFilterMeteoQualityFlag = TRUE, isNeglectVPDEffect = TRUE)
ds <- partGLExtractStandardData(dsNEE1, controlGLPart = ctrl)
#plot(VPD ~ sDateTime, ds)
#
dsTempSens <- REddyProc:::partGLFitNightTimeTRespSens( ds
, nRecInDay = 48L
, controlGLPart = ctrl
)
resLRC <- REddyProc:::partGLFitLRCWindows(
ds, dsTempSens = dsTempSens, lrcFitter = RectangularLRCFitter() )
expect_true( resLRC$iMeanRec[2] == resLRC$iMeanRec[3])
#resLRC$nValidRec
#
resLRC2 <- REddyProc:::partGLFitLRCWindows(
ds, dsTempSens = dsTempSens, lrcFitter = RectangularLRCFitter()
,controlGLPart = ctrl )
expect_true( resLRC2$iMeanRec[2] != resLRC2$iMeanRec[3])
# used all records despite missing VPD
expect_true( all( resLRC2$nValidRec >= 50L) )
#resLRC2
#
resLRC2C <- REddyProc:::partGLFitLRCWindows(
ds, dsTempSens = dsTempSens, lrcFitter = RectangularLRCFitterCVersion()
,controlGLPart = ctrl )
# used all records despite missing VPD
expect_true( all( resLRC2C$nValidRec >= 50L) )
# equal unless bootstrapped uncertainty and optimization result
iSd <- grep("_sd$|resOpt$",names(resLRC2))
expect_equal( as.data.frame(resLRC2[,-iSd]), as.data.frame(resLRC2C[,-iSd]) )
#
tmp <- partitionNEEGL( dsNEE1, controlGLPart = ctrl )
#tmp[ resLRC2$iCentralRec,]
expect_equal( tmp[ resLRC2$iCentralRec,"FP_alpha"], resLRC2$alpha )
expect_equal( nrow(dsNEE1), nrow(tmp) )
#
# finite prediction despite missing VPD
expect_true( all(is.finite(tmp$GPP_DT)))
expect_true( all(tmp$GPP_DT >= 0))
expect_true( all(tmp$GPP_DT < 250))
expect_true( all(tmp$Reco_DT < 10))
expect_true( all(tmp$Reco_DT > 0))
expect_true( all(tmp$Reco_DT_SD > 0))
expect_true( all(tmp$GPP_DT_SD >= 0))
})
test_that("partGLPartitionFluxes missing night time data",{
skip_if_not_installed("mlegp")
dsNEE1 <- dsNEE
#set all data except one day to NA to associate the same data
#with several windows
dsNEE1$hourOfDay <- as.POSIXlt(dsNEE1$sDateTime)$hour
dsNEE1$NEE_f[ (dsNEE1$sDateTime >= as.POSIXct("1998-06-01 00:00:00",tz = tzEx)) &
(dsNEE1$sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx)) &
((dsNEE1$hourOfDay >= 18) | (dsNEE1$hourOfDay <= 6))] <- NA
ds <- dsNEE1
ds$NEE <- ds$NEE_f
ds$sdNEE <- ds$NEE_fsd
ds$Temp <- ds$Tair_f
ds$VPD <- ds$VPD_f
ds$Rg <- ds$Rg_f
#
dsTempSens <- REddyProc:::partGLFitNightTimeTRespSens( ds
, nRecInDay = 48L
, controlGLPart = partGLControl()
, winSizeNight = 4L, winExtendSizes = c()
)
expect_true( all( is.finite(dsTempSens$RRef)) )
resLRC <- REddyProc:::partGLFitLRCWindows(
ds, lrcFitter = RectangularLRCFitter(), dsTempSens = dsTempSens )
expect_true( all( is.finite(resLRC$RRef_night)) )
})
test_that("partGLPartitionFluxes filter Meteo flag not enough without VPD",{
skip_if_not_installed("mlegp")
dsNEE1 <- dsNEE
# test setting VPD_fqc to other than zero, for omitting those rows
dsNEE1$VPD_fqc[ isTimeInTestPeriod(dsNEE1$sDateTime) ] <- 1L
#plot( VPD_fqc ~ sDateTime, dsNEE1 )
#plot( NEE_f ~ sDateTime, dsNEE1 )
ds <- dsNEE1
#ds$NEE <- ds$NEE_f
#ds$sdNEE <- ds$NEE_fsd
#ds$Temp <- ds$Tair_f
#ds$VPD <- ds$VPD_f
ds$Rg <- ds$Rg_f
#
.tmp.f <- function(){
tmp0 <- partitionNEEGL(
dsNEE1, controlGLPart = partGLControl(isFilterMeteoQualityFlag = TRUE) )
tmpPar0 <- tmp0[is.finite(tmp0$FP_beta),]
}
tmp <- partitionNEEGL(
ds, controlGLPart = partGLControl(isFilterMeteoQualityFlag = TRUE) )
expect_equal( nrow(ds), nrow(tmp) )
tmpPar <- tmp[is.finite(tmp$FP_RRef_Night),]
# estimate despite missing VPD
expect_equal( sum(is.finite(tmpPar[,"FP_beta"])), 4L )
expect_equal( tmpPar[4L,"FP_k"], 0 ) # k = 0 for missing VPD
#
# now set temperature to bad quality flag,
dsNEE1$Tair_fqc[ isTimeInTestPeriod(dsNEE1$sDateTime) ] <- 1L
#dsNEE1$VPD_fqc[ (dsNEE1$sDateTime > as.POSIXct("1998-06-03 00:00:00",tz = tzEx)) &
#(dsNEE1$sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx) |
#dsNEE1$sDateTime >= as.POSIXct("1998-06-06 00:00:00",tz = tzEx)) ] <- 1L
ds <- dsNEE1
ds$Rg <- ds$Rg_f
tmp <- partitionNEEGL(
ds, controlGLPart = partGLControl(isFilterMeteoQualityFlag = TRUE) )
expect_equal( nrow(ds), nrow(tmp) )
tmpPar <- tmp[is.finite(tmp$FP_RRef_Night),]
# one row less with parameter estimates
expect_equal( sum(is.finite(tmpPar[,"FP_beta"])), 3L )
})
test_that("partGLPartitionFluxes missing prediction VPD",{
skip_if_not_installed("mlegp")
dsNEE1 <- dsNEE
# test setting VPD_fqc to other than zero, for omitting those rows
dsNEE1$VPD_fqc[ (dsNEE1$sDateTime > as.POSIXct("1998-06-03 00:00:00",tz = tzEx)) &
(dsNEE1$sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx) |
dsNEE1$sDateTime >= as.POSIXct("1998-06-06 00:00:00",tz = tzEx)) ] <- 1L
#plot( VPD_fqc ~ sDateTime, dsNEE1 )
#plot( NEE_f ~ sDateTime, dsNEE1 )
iMissingVPD <- c(10,400)
dsNEE1$VPD_f[iMissingVPD] <- NA
ds <- dsNEE1
#ds$NEE <- ds$NEE_f
#ds$sdNEE <- ds$NEE_fsd
#ds$Temp <- ds$Tair_f
#ds$VPD <- ds$VPD_f
ds$Rg <- ds$Rg_f
#
expect_warning(tmp <- partitionNEEGL(
ds, controlGLPart = partGLControl(
isFilterMeteoQualityFlag = TRUE
, isRefitMissingVPDWithNeglectVPDEffect = FALSE) ))
expect_equal( nrow(ds), nrow(tmp) )
tmpPar <- tmp[is.finite(tmp$FP_RRef_Night),]
# estimate despite missing VPD
expect_equal( sum(is.finite(tmpPar[,"FP_beta"])), 4L )
expect_equal( tmpPar[4L,"FP_k"], 0 ) # k = 0 for missing VPD
expect_true( is.na(tmp$GPP_DT[iMissingVPD[1]]) ) # default could not predict
# second: both parameter sets with k = 0 -> VPD not evaluated for prediction
expect_true( is.finite(tmp$GPP_DT[iMissingVPD[2]]) )
#
# repeat with (default) refitting with neglecting VPD
tmp2 <- partitionNEEGL(
ds, controlGLPart = partGLControl(
isFilterMeteoQualityFlag = TRUE
, isRefitMissingVPDWithNeglectVPDEffect = TRUE) )
expect_equal( tmp2[-iMissingVPD,"GPP_DT"], tmp[-iMissingVPD,"GPP_DT"] )
expect_equal( tmp2[-iMissingVPD,"Reco_DT"], tmp[-iMissingVPD,"Reco_DT"] )
expect_equal( tmp2[iMissingVPD[2],c("GPP_DT","Reco_DT")],
tmp2[iMissingVPD[2],c("GPP_DT","Reco_DT")])
# now also estimate for first missing VPD
expect_true( is.finite(tmp2$GPP_DT[iMissingVPD[1]]) )
})
test_that("partitionNEEGL fixed tempSens",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
#
ds <- dsNEE1
ds$NEE <- ds$NEE_f
ds$sdNEE <- ds$NEE_fsd
ds$Temp <- ds$Tair_f
ds$VPD <- ds$VPD_f
ds$Rg <- ds$Rg_f
#
dsTempSens0 <- REddyProc:::partGLFitNightTimeTRespSens( ds )
startRecs <- REddyProc:::getStartRecsOfWindows( nrow(ds) )
fixedTempSens <- data.frame(E0 = 80, sdE0 = 10 )
tmp <- partitionNEEGL(
dsNEE1, RadVar = "Rg_f", controlGLPart = partGLControl(
fixedTempSens = fixedTempSens) )
expect_equal( nrow(dsNEE1), nrow(tmp) )
#
expect_true( all(is.finite(tmp$GPP_DT)))
expect_true( all(tmp$GPP_DT >= 0))
expect_true( all(tmp$GPP_DT < 250))
expect_true( all(tmp$Reco_DT < 10))
expect_true( all(tmp$Reco_DT > 0))
expect_true( all(tmp$Reco_DT_SD > 0))
expect_true( all(tmp$GPP_DT_SD >= 0))
.tmp.plot <- function(){
tmp$time <- dsNEE1$sDateTime
plot( Reco_DT ~ time, tmp)
#plot( diff(Reco_DT_u50) ~ time[-1], tmp)
plot( GPP_DT ~ time, tmp)
}
})
test_that("partitionNEEGL: missing sdNEE",{
skip_if_not_installed("mlegp")
dsNEE1 <- dsNEE
#summary(dsNEE1)
dsNEE1$NEE_fsd <- NA_real_
#
tmp <- partitionNEEGL(
dsNEE1, RadVar = "Rg_f", controlGLPart = partGLControl() )
tmpWithPar <- tmp[ is.finite(tmp$FP_alpha),]
expect_equal( nrow(tmpWithPar), 4) # fitted with missing fsd
#
rm(tmp)
expect_error(
tmp <- partitionNEEGL(
dsNEE1, RadVar = "Rg_f"
, controlGLPart = partGLControl(replaceMissingSdNEEParms = c(NA,NA)) )
)
})
# see .profilePartGL in develop/profile/profile.R
test_that("partitionNEEGL long gap",{
skip_if_not_installed("mlegp")
skip("cannot repoduce SmoothTempSens failing")
dsNEE1 <- Example_DETha98_Filled
# introduce a long four month gap
bo <- dsNEE1$DoY >= 80 & dsNEE1$DoY <= 360
dsNEE1$NEE_f[bo] <- NA
dsNEE1$NEE_fqc[bo] <- 3
tmp <- partitionNEEGL( dsNEE1)
tmp2 <- cbind(select(dsNEE1, sDateTime), tmp)
plot( FP_qc ~ sDateTime, tmp2)
plot( GPP_DT ~ sDateTime, tmp2)
})
test_that("partitionNEETK",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
resEx <- resLRCEx1
#DoY.V.n <- as.POSIXlt(dsNEE1$sDateTime)$yday + 1L
#Hour.V.n <- as.POSIXlt(dsNEE1$sDateTime)$hour + as.POSIXlt(dsNEE1$sDateTime)$min/60
#dsNEE1$PotRad_NEW <- fCalcPotRadiation(DoY.V.n, Hour.V.n, LatDeg = 45.0, LongDeg = 1, TimeZoneHour = 0 )
tmp <- partitionNEEGL(
dsNEE1,RadVar = 'Rg_f'
, controlGLPart = partGLControl(useNightimeBasalRespiration = TRUE))
#tmp <- partitionNEEGL( dsNEE1,RadVar = 'Rg_f', controlGLPart = partGLControl(nBootUncertainty = 0L, isAssociateParmsToMeanOfValids = FALSE))
expect_equal( nrow(dsNEE1), nrow(tmp) )
#tmp[ is.finite(tmp$FP_beta), ] # note FP_dRecPar is not zero, because iCentralRec != iMeanRec
#iPar <- which(is.finite(tmp$FP_E0))
#plot( tmp$FP_beta[iPar] ~ dsNEE1$sDateTime[iPar],type = "l", ylim = range(c(tmp$FP_beta[iPar],tmp$FP_GPP2000[iPar]))); lines( tmp$FP_GPP2000[iPar] ~ dsNEE1$sDateTime[iPar], col = "red")
#
# now test with different suffix: u50
dsNEE2 <- dsNEE
names(dsNEE2)[ match(c("NEE_f", "NEE_fqc", "NEE_fsd"),names(dsNEE2))] <-
c("NEE_u50_f", "NEE_u50_fqc", "NEE_u50_fsd")
tmp <- partitionNEEGL(
dsNEE2, RadVar = 'Rg_f', suffix = "u50"
, controlGLPart = partGLControl(
nBootUncertainty = 0L, isAssociateParmsToMeanOfValids = FALSE
,useNightimeBasalRespiration = TRUE )
)
expect_equal( nrow(dsNEE1), nrow(tmp) )
expect_true( all(is.finite(tmp$GPP_DT_u50)))
expect_true( all(tmp$GPP_DT_u50 >= 0))
expect_true( all(tmp$GPP_DT_u50 < 250))
expect_true( all(tmp$Reco_DT_u50 < 10))
expect_true( all(tmp$Reco_DT_u50 > 0))
expect_true( all(tmp$Reco_DT_u50_SD >= 0))
expect_true( all(tmp$GPP_DT_u50_SD >= 0))
# jumps between daytime and nighttime
#expect_true( all(abs(diff(tmp$Reco_DT_u50)) < 0.6)) #smooth
# reporting good values at central records
# tmp[resEx$iCentralRec,]
iRowsOpt <- which(is.finite(tmp$FP_E0))
expect_true( length(iRowsOpt) == nrow(resEx) )
#expect_true( all((tmp$FP_alpha[resEx$iCentralRec] - resEx$alpha)[
#resEx$parms_out_range == 0L] < 1e-2) )
.tmp.plot <- function(){
tmp$time <- dsNEE1$sDateTime
plot( Reco_DT_u50 ~ time, tmp)
#plot( diff(Reco_DT_u50) ~ time[-1], tmp)
plot( GPP_DT_u50 ~ time, tmp)
plot( FP_RRef_Night ~ time, tmp)
plot( FP_RRef ~ FP_RRef_Night, tmp)
plot( Reco_DT_u50_SD ~ time, tmp)
}
})
test_that("report missing VPD_f column in error",{
skip("only interactively test issue #34")
DETha98 <- fConvertTimeToPosix(Example_DETha98, 'YDH', Year = 'Year',
Day = 'DoY', Hour = 'Hour')[-(2:4)]
EProc <- sEddyProc$new('DE-Tha', DETha98,
c('NEE', 'Rg', 'Tair', 'VPD', 'Ustar'))
EProc$sMDSGapFillAfterUstar('NEE', uStarTh = 0.3, FillAll = TRUE)
# missing 'VPD'...
for (i in c('Tair', 'Rg')) EProc$sMDSGapFill(i, FillAll = TRUE)
EProc$sSetLocationInfo(LatDeg = 51.0, LongDeg = 13.6, TimeZoneHour = 1)
# ...produces Error here
# Error should report missing column VPD_f
expect_error(
EProc$sGLFluxPartition(suffix = "uStar")
,"VPD_f"
)
})
test_that("no nighttime data",{
skip_if_not_installed("mlegp")
skip_on_cran()
cols <- c('NEE','Rg','Tair','VPD', 'Ustar'
,"Tair_f","VPD_f","Rg_f","NEE_f"
,"NEE_fqc", "Tair_fqc", "Rg_fqc", "NEE_fsd")
#ds <- Example_DETha98_Filled[,c('DateTime',cols)]
ds <- subset(Example_DETha98_Filled
, sDateTime >= as.POSIXct("1998-05-01 00:15:00",tz = tzEx) &
sDateTime <= as.POSIXct("1998-09-01 21:45:00",tz = tzEx))
ds <- ds[,c('DateTime',cols)]
ds$NEE_f[ds$Rg_f <= 10] <- NA # simulate having no night-time records
EProc <- suppressWarnings(sEddyProc$new('DE-ThaSummer', ds, cols))
EProc$sSetLocationInfo(LatDeg = 51.0, LongDeg = 13.6, TimeZoneHour = 1)
EProc$sGLFluxPartition(controlGLPart = partGLControl(fixedTempSens = data.frame(
E0 = 150, sdE0 = 50, RRef = 20)))
})
.test_sEddyProc_sGLFluxPartitionUStarScens <- function(){}
test_that("sEddyProc_sGLFluxPartitionUStarScens wrong suffix",{
skip_if_not_installed("mlegp")
dsTest <- Example_DETha98_Filled %>% mutate(
NEE_uStar_f = NEE, NEE_uStar_fqc = NEE_fqc, NEE_uStar_fsd = NEE_fsd,
NEE_U50_f = NEE, NEE_U50_fqc = NEE_fqc, NEE_U50_fsd = NEE_fsd
) %>%
select(!"sDateTime")
EProc <- sEddyProc$new(
'DE-Tha', dsTest, setdiff(names(dsTest), c("NEE_fnum","NEE_fwin")))
EProc$sSetUStarSeasons(1L)
EProc$sSetUstarScenarios(data.frame(season = 1L, uStar = 0.28, U50 = 0.38))
EProc$sSetLocationInfo(LatDeg = 51.0, LongDeg = 13.6, TimeZoneHour = 1)
expect_error(
EProc$sGLFluxPartitionUStarScens(uStarScenKeep = "nonexisting")
,"nonexisting")
})
test_that("sEddyProc_sGLFluxPartitionUStarScens",{
skip_if_not_installed("mlegp")
dsTest <- Example_DETha98_Filled %>% mutate(
NEE_uStar_f = NEE, NEE_uStar_fqc = NEE_fqc, NEE_uStar_fsd = NEE_fsd,
NEE_U50_f = NEE, NEE_U50_fqc = NEE_fqc, NEE_U50_fsd = NEE_fsd
) %>%
select(!"sDateTime")
EProc <- sEddyProc$new(
'DE-Tha', dsTest, setdiff(names(dsTest), c("NEE_fnum","NEE_fwin")))
EProc$sSetUStarSeasons(1L)
EProc$sSetUstarScenarios(data.frame(season = 1L, uStar = 0.28, U50 = 0.38))
EProc$sSetLocationInfo(LatDeg = 51.0, LongDeg = 13.6, TimeZoneHour = 1)
#EProc$sMDSGapFill('Tair', FillAll = FALSE, minNWarnRunLength = NA)
#EProc$sMDSGapFill('VPD', FillAll = FALSE, minNWarnRunLength = NA)
#EProc$trace(sGLFluxPartitionUStarScens, recover); # EProc$untrace(sGLFluxPartitionUStarScens)
#EProc$sGLFluxPartitionUStarScens(uStarScenKeep = "U50")
EProc$sGetUstarScenarios()
EProc$sGetUstarSuffixes()
EProc$sGLFluxPartitionUStarScens(warnOnOtherErrors=TRUE, uStarScenKeep="U50")
expect_true(all(c("GPP_DT_uStar","GPP_DT_U50", "GPP_DT_U50_SD") %in%
names(EProc$sTEMP)))
})
bgctw/REddyProc documentation built on March 26, 2024, 11:35 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#+++ Unit tests for fConvertTimeToPosix functions +++
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# Author: TW
#require(testthat)
context("partGL")
if (!exists(".partGPAssociateSpecialRows")) .partGPAssociateSpecialRows <-
REddyProc:::.partGPAssociateSpecialRows
if (!exists(".binUstar")) .binUstar <- REddyProc:::.binUstar
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# 10 days from June from Example_DETha98.txt shipped with REddyProc
Example_DETha98_Filled <- getFilledExampleDETha98Data()
tzEx <- REddyProc:::getTZone(Example_DETha98_Filled$sDateTime)
test_that("example dataset starts at midngiht",{
expect_true( Example_DETha98_Filled$sDateTime[1] ==
as.POSIXct("1998-01-01 00:15:00",tz = tzEx) )
})
dsNEE <- subset(Example_DETha98_Filled
, sDateTime >= as.POSIXct("1998-06-01 00:15:00",tz = tzEx) &
sDateTime <= as.POSIXct("1998-06-09 21:45:00",tz = tzEx))
dsNEE$Temp <- dsNEE$Tair_f
dsNEE$isNight <- (dsNEE$Rg_f <= 4 & dsNEE$PotRad_NEW == 0)
dsNEE$isDay = (dsNEE$Rg_f > 4 & dsNEE$PotRad_NEW != 0)
# 8 first days of June from IT-MBo.2005.txt
# 10 days from June from Example_DETha98.txt shipped with REddyProc
.tmp.f <- function(){
#save(dsNEE, file = "tmp/dsNEE_Tharandt.RData")
load("tmp/dsNEE_Tharandt.RData") # dsNEE
dsNEE$Temp <- dsNEE$Tair_f
dsNEE$Rg_f <- dsNEE$Rg
dsNEE$isNight <- (dsNEE$Rg_f <= 4 & dsNEE$PotRad_NEW == 0)
dsNEE$isDay = (dsNEE$Rg_f > 4 & dsNEE$PotRad_NEW != 0)
}
.tmp.f2 <- function(){
# else stop in partitionNEEGL, and grap ds:
attr(ds$sDateTime, "tzone") <- "UTC"
dsJune <- ds
dsNEE <- dsJune[ dsJune$sDateTime >= as.POSIXct("1998-06-01", tz = "UTC") &
dsJune$sDateTime < as.POSIXct("1998-06-10", tz = "UTC")
, c("sDateTime","NEE_f","NEE_fqc","NEE_fsd","Tair_f","Tair_fqc","VPD_f"
,"VPD_fqc","Rg_f","PotRad_NEW")]
save(dsNEE, file = "tmp/dsNEE_Tharandt.RData")
dput(dsNEE)
}
# regression result from dput in
resLRCEx1 <- structure(list(
iWindow = 1:4, dayStart = c(1L, 3L, 5L, 7L), dayEnd = c(4L, 6L, 8L, 9L)
, iRecStart = c(1L, 97L, 193L, 289L), iRecEnd = c(192L, 288L, 384L, 428L)
, iCentralRec = c(97L, 193L, 289L, 385L), nValidRec = c(98L,98L, 75L, 56L)
, iMeanRec = c(97L, 194L, 283L, 385L)
, convergence = c(0L,0L, 0L, 0L), parms_out_range = c(0L, 0L, 0L, 0L)
, k = c(
0.15556543911127, 0.100336071150363, 0.0292272454449355, 0.110909241852741)
, beta = c(
34.2363908612952, 28.7053014210465, 11.965007963525, 35.6061046130658)
, alpha = c(
0.0408287988682114, 0.0317164138158062, 0.12184364495971, 0.0565989757582912)
, RRef = c(
1.30588274435898, 2.30342378833466, 4.39970416578445, 3.02347421699605)
, E0 = structure(
c(62.633343296311, 62.633343296311, 62.633343296311, 62.633343296311)
, .Dim = c(4L, 1L), .Dimnames = list(NULL, NULL))
, k_sd = c(
0.0162347703672834, 0.00819228863397153, 0.0036484758449877, 0.0143126637303779)
, beta_sd = c(
2.6633120785356, 2.53043346467805, 0.480085981083196, 3.86498333963868)
, alpha_sd = c(
0.00282359981111383, 0.00245807510249875, 0.0189581631026536, 0.00624526974304941)
, RRef_sd = c(
0.243651749533536, 0.200219643557999, 0.379555581976541, 0.408032180867498)
, E0_sd = structure(
c(4.25106879296635, 4.22511397480583, 4.22511362599932, 4.46593783063789)
, .Dim = c(4L, 1L), .Dimnames = list(NULL, NULL))
, GPP2000 = c(
24.1225750077742, 19.7622684102263, 11.4050231061514, 27.0862112802389)
, isValid = c(TRUE, TRUE, TRUE, TRUE)
, resOpt = list(structure(list(thetaOpt = structure(
c(0.15556543911127, 34.2363908612952, 0.0408287988682114
, 1.30588274435898, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, iOpt = 1:5, thetaInitialGuess = structure(c(
0.05, 24.3542, 0.1, 5.17329116845493, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, covParms = structure(c(
0.000263567768878423, 0.0299004485274175, -2.01722085100335e-05
, -0.00170017321010028, 0, 0.0299004485274174, 7.09323122767363
, -0.00442973195793744, -0.184157427431281, 0, -2.01722085100335e-05
, -0.00442973195793744, 7.97271589332207e-06, 0.000574595628371077, 0
, -0.00170017321010028, -0.18415742743128, 0.000574595628371077
, 0.0593661750507529, 0, 0, 0, 0, 0, 18.0715858825324)
, .Dim = c(5L, 5L), .Dimnames = list(
structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))))
, convergence = 0L)
, .Names = c(
"thetaOpt", "iOpt", "thetaInitialGuess", "covParms", "convergence"))
, structure(list(thetaOpt = structure(
c(0.100336071150363, 28.7053014210465, 0.0317164138158062
, 2.30342378833466, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, iOpt = 1:5
, thetaInitialGuess = structure(
c(0.05, 24.8403, 0.1, 5.14770049230269, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, covParms = structure(
c(6.71135930622992e-05, 0.0178801034106948, -1.74280275277554e-05
, -0.00134394852514135, 0, 0.0178801034106948, 6.40309351916256
, -0.00493560933704443, -0.271156219332448, 0, -1.74280275277554e-05
, -0.00493560933704443, 6.04213320952425e-06, 0.000418391213333105
, 0, -0.00134394852514135, -0.271156219332449, 0.000418391213333105
, 0.040087905666492, 0, 0, 0, 0, 0, 17.8515881000995)
, .Dim = c(5L, 5L), .Dimnames = list(
structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))))
, convergence = 0L)
, .Names = c(
"thetaOpt", "iOpt", "thetaInitialGuess", "covParms", "convergence"))
, structure(list(thetaOpt = structure(
c(0.0292272454449355, 11.965007963525, 0.12184364495971
, 4.39970416578445, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, iOpt = 1:5
, thetaInitialGuess = structure(
c(0.05, 24.1194, 0.1, 5.14770049230269, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, covParms = structure(
c(1.33113759914587e-05, 0.000439566354933103, -5.94748188903169e-05
, -0.00111894361237909, 0, 0.000439566354933103, 0.230482549232615
, -0.000116871653668728, 0.053260722726946, 0, -5.94748188903169e-05
, -0.000116871653668727, 0.000359411948226816, 0.00657062463616557
, 0, -0.00111894361237909, 0.0532607227269459, 0.00657062463616558
, 0.144062439809551, 0, 0, 0, 0, 0, 17.8515851526051)
, .Dim = c(5L, 5L), .Dimnames = list(
structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, structure(c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))))
, convergence = 0L)
, .Names = c(
"thetaOpt", "iOpt", "thetaInitialGuess", "covParms", "convergence"))
, structure(list(thetaOpt = structure(
c(0.110909241852741, 35.6061046130658, 0.0565989757582912
, 3.02347421699605, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, iOpt = 1:5
, thetaInitialGuess = structure(
c(0.05, 21.665, 0.1, 5.21501814931109, 62.633343296311)
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, covParms = structure(c(
0.000204852343058874, 0.0495123146513288, -6.84478020597351e-05
, -0.00386015942484132, 0, 0.0495123146513288, 14.9380962156846
, -0.0178713249608098, -0.786540464781837, 0, -6.8447802059735e-05
, -0.0178713249608097, 3.90033941634485e-05, 0.00226964418347153
, 0, -0.00386015942484131, -0.786540464781834, 0.00226964418347153
, 0.166490260623486, 0, 0, 0, 0, 0, 19.9446007071226)
, .Dim = c(5L, 5L), .Dimnames = list(structure(
c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))
, structure(
c("k", "beta", "alpha", "RRef", "E0")
, .Names = c("k", "beta", "alpha", "RRef", "E0"))))
, convergence = 0L)
, .Names = c(
"thetaOpt", "iOpt", "thetaInitialGuess", "covParms", "convergence")))
, E0_bootstrap_sd = c(
4.25106879296635, 4.22511397480583, 4.22511362599932, 4.46593783063789)
, RRef_night = c(
5.17329116845493, 5.14770049230269, 5.14770049230269, 5.21501814931109))
, .Names = c(
"iWindow", "dayStart", "dayEnd", "iRecStart", "iRecEnd", "iCentralRec"
, "nValidRec", "iMeanRec", "convergence", "parms_out_range", "k", "beta"
, "alpha", "RRef", "E0", "k_sd", "beta_sd", "alpha_sd", "RRef_sd", "E0_sd"
, "GPP2000", "isValid", "resOpt", "E0_bootstrap_sd", "RRef_night")
, row.names = c(NA, -4L)
, class = c("tbl_df", "tbl", "data.frame")
)
resLRCEx1Nonrectangular <- structure(list(
iWindow = 1:4, dayStart = c(1L, 3L, 5L, 7L), dayEnd = c(4L, 6L, 8L, 9L)
, iRecStart = c(1L, 97L, 193L, 289L), iRecEnd = c(192L, 288L, 384L, 428L)
, iCentralRec = c(97L, 193L, 289L, 385L), nValidRec = c(98L, 98L, 75L, 56L)
, iMeanRec = c(97L, 194L, 283L, 385L), convergence = c(1005L, 1005L, 1005L, 0L)
, E0 = structure(c(
62.633343296311, 62.633343296311, 62.633343296311, 62.633343296311)
, .Dim = c(4L, 1L), .Dimnames = list( NULL, NULL))
, E0_sd = structure(
c(4.25106879296635, 4.22511397480583, 4.22511362599932, 4.46593783063789)
, .Dim = c(4L, 1L), .Dimnames = list( NULL, NULL))
, RRefNight = c(5.17329116845493, 5.14770049230269, 5.14770049230269, NA)
, parms_out_range = c(NA, NA, NA, 0L), k = c(NA, NA, NA, 0.15301818456226)
, beta = c(NA, NA, NA, 49.3407957035085 )
, alpha = c(NA, NA, NA, 0.0444237102808321)
, RRef = c(NA, NA, NA, 2.35368171175251)
, logitconv = c(NA, NA, NA, -3.10163338239329 )
, k_sd = c(NA, NA, NA, 0.0197088986436212)
, beta_sd = c(NA, NA, NA, 7.58853578454936)
, alpha_sd = c(NA, NA, NA, 0.00432079403235484 )
, RRef_sd = c(NA, NA, NA, 0.350880271587216)
, logitconv_sd = c(NA, NA, NA, 1.92633463362266)
, GPP2000 = c(NA, NA, NA, 32.0432128516724 )
, isValid = c(NA, NA, NA, TRUE)
, resOpt = list(NULL, NULL, NULL, structure(
list(thetaOpt = structure(c(
0.15301818456226, 49.3407957035085, 0.0444237102808321, 2.35368171175251
, 62.633343296311, -3.10163338239329 )
, .Names = c("k", "beta", "alpha", "RRef", "E0", "logitconv" ))
, iOpt = c(1, 2, 3, 4, 6, 5)
, thetaInitialGuess = structure(c(
0.05, 21.665, 0.1, 5.21501814931109, 62.633343296311, 1.09861228866811 )
, .Names = c("k", "beta", "alpha", "RRef", "E0", "logitconv" ))
, covParms = structure(c(
0.000388440685744535, 0.13746635747882, -6.51906051810543e-05
, -0.00449519721692684, 0, 0.00207560732637497, 0.13746635747882
, 57.5858753533862, -0.0238832979670161, -1.36443384766771, 0
, -1.52788892519458, -6.51906051810543e-05, -0.023883297967016
, 1.86692610700332e-05, 0.00131289852541519, 0, -0.00181474371141437
, -0.00449519721692684, -1.36443384766771, 0.00131289852541519
, 0.123116964989119, 0, -0.0864363703508266, 0, 0, 0, 0, 19.9446007071226
, 0, 0.00207560732637499, -1.52788892519458, -0.00181474371141437
, -0.0864363703508268, 0, 3.71076512069413 )
, .Dim = c(6L, 6L), .Dimnames = list(structure(c(
"k", "beta", "alpha", "RRef", "E0", "logitconv")
, .Names = c("k", "beta", "alpha", "RRef", "E0", "logitconf"))
, structure(c(
"k", "beta", "alpha", "RRef", "E0", "logitconv")
, .Names = c("k", "beta", "alpha", "RRef", "E0", "logitconf"))))
, convergence = 0L)
, .Names = c(
"thetaOpt", "iOpt", "thetaInitialGuess", "covParms", "convergence")))
, E0_bootstrap_sd = c(
4.25106879296635, 4.22511397480583, 4.22511362599932, 4.46593783063789)
, RRef_night = c(
5.17329116845493, 5.14770049230269, 5.14770049230269, 5.21501814931109))
, .Names = c(
"iWindow", "dayStart", "dayEnd", "iRecStart", "iRecEnd", "iCentralRec"
, "nValidRec", "iMeanRec", "convergence", "E0", "E0_sd", "RRefNight"
, "parms_out_range", "k", "beta", "alpha", "RRef", "logitconv", "k_sd"
, "beta_sd", "alpha_sd", "RRef_sd", "logitconv_sd", "GPP2000", "isValid"
, "resOpt", "E0_bootstrap_sd", "RRef_night")
, row.names = c(NA, -4L)
, class = c("tbl_df", "tbl", "data.frame")
)
test_that("replaceMissingSdByPercentage",{
n <- 100
x <- rnorm(n,10)
sdX <- sdX0 <- pmax(0.7, rnorm(n, 10*0.1 ))
iMissing <- sample(n, 20)
sdX[iMissing] <- NA
sdXf <- REddyProc:::replaceMissingSdByPercentage(sdX,x, 0.2, 1.7)
#plot( sdXf ~ x, col = "blue" );points( sdX ~ x )
expect_equal( sdXf[-iMissing], sdX[-iMissing])
expect_true( all(is.finite(sdXf)))
expect_true( all(sdXf[iMissing] >= 1.7))
expect_true( all(sdXf[iMissing] >= 0.2*x[iMissing]))
#
sdX <- sdX0
iMissing <- TRUE
sdX[iMissing] <- NA
sdXf <- REddyProc:::replaceMissingSdByPercentage(sdX,x, 0.2, 1.7)
#plot( sdXf ~ x, col = "blue" );points( sdX ~ x )
expect_true( all(is.finite(sdXf)))
expect_true( all(sdXf[iMissing] >= 1.7))
expect_true( all(sdXf[iMissing] >= 0.2*x[iMissing]))
#
sdX <- sdX0
iMissing <- sample(n, 20)
sdX[iMissing] <- NA
sdXf <- REddyProc:::replaceMissingSdByPercentage(sdX,x, NA, 1.7)
expect_true( all(sdXf[iMissing] == 1.7))
#plot( sdXf ~ x, col = "blue" );points( sdX ~ x )
sdXf <- REddyProc:::replaceMissingSdByPercentage(sdX,x, 0.2, NA)
expect_true( all(sdXf[iMissing] == 0.2*x[iMissing]))
sdXf <- REddyProc:::replaceMissingSdByPercentage(sdX,x, NA, NA)
expect_true( all(is.na(sdXf[iMissing])))
})
test_that("estimating temperature sensitivity oneWindow are in accepted range",{
skip_if_not_installed("mlegp")
dss <- dsNEE[ dsNEE$Rg_f <= 0 & dsNEE$PotRad_NEW <= 0
& as.POSIXlt(dsNEE$sDateTime)$mday %in% 1:8, ]
dss <- dss[ order(dss$Temp), ]
dss$NEE <- dss$NEE_f
resE0 <- REddyProc:::partGLEstimateTempSensInBoundsE0Only(
dss$NEE_f, dss$Temp + 273.15)
expect_true( resE0$E0 >= 50 && resE0$E0 < 400 )
medianResp <- median(dss$NEE_f,na.rm = TRUE)
expect_true( abs(resE0$RRefFit - medianResp)/medianResp < 0.2 )
E0Win <- as.data.frame(resE0)
res <- REddyProc:::partGLFitNightRespRefOneWindow(
dss, data.frame(iWindow = 1L), E0Win = E0Win)
RRef <- res[1]
expect_true( RRef >= 0)
.tmp.plot <- function(){
plot( NEE_f ~ Temp, dss) # FP_VARnight negative?
lines( fLloydTaylor(RRef, resE0$E0, dss$Temp+273.15, TRef = 273.15+15) ~
dss$Temp)
}
})
test_that("estimating temperature sensitivity on record with some freezing temperatures",{
skip_if_not_installed("mlegp")
dss <- dsNEE[ dsNEE$Rg_f <= 0 & dsNEE$PotRad_NEW <= 0 &
as.POSIXlt(dsNEE$sDateTime)$mday %in% 1:8, ]
dss$NEE <- dss$NEE_f
dss <- dss[ order(dss$Temp), ]
# only the last 13 records will have temperature above -1degC
dss$Temp <- dss$Temp - dss$Temp[ nrow(dss) - 14L ] - 1
isValid <- REddyProc:::isValidNightRecord(dss)
expect_true(sum(isValid) >= 13 &&
sum(REddyProc:::isValidNightRecord(dss)) < nrow(dss) )
#
resE0 <- REddyProc:::partGLEstimateTempSensInBoundsE0Only(
dss$NEE_f[isValid], dss$Temp[isValid] + 273.15)
expect_true( is.na(resE0$E0) )
.tmp.f <- function(){
expect_true( resE0$E0 >= 50 && resE0$E0 <= 400 )
medianResp <- median(dss$NEE_f[isValid],na.rm = TRUE)
expect_true( abs(resE0$RRefFit - medianResp)/medianResp < 0.2 )
E0Win <- as.data.frame(resE0)
res <- REddyProc:::partGLFitNightRespRefOneWindow(
dss, data.frame(iWindow = 1L), E0Win = E0Win)
RRef <- res[[2]]$RRef[1]
expect_true( RRef >= 0 )
}
.tmp.plot <- function(){
plot( NEE_f ~ Temp, dss) # FP_VARnight negative?
points( NEE_f ~ Temp, dss[isValid,], col = "red")
lines( fLloydTaylor(RRef, resE0$E0, dss$Temp+273.15, TRef = 273.15+15) ~
dss$Temp)#
}
})
test_that("applyWindows",{
nRec <- nrow(dsNEE)
nRecInDay <- 10L
ds <- within(dsNEE, {
iRec <- 1:nRec
# specifying the day for each record assuming equidistand records
iDayOfRec <- ((c(1:nRec) - 1L) %/% nRecInDay) + 1L
})
fReportTime <- function(dss, winInfo, prevRes){
nRecS <- nrow(dss)
list( res2 = data.frame(
startRec = dss$iRec[1]
,endRec = dss$iRec[nRecS]
,startDay = dss$iDayOfRec[1]
,endDay = dss$iDayOfRec[nRecS]
)
,sumRes = prevRes$sumRes + 1L
)
}
prevRes <- list(sumRes = 0)
fReportTime(ds,0,prevRes)
# larger than reference window of 4 days
resApply <- REddyProc:::applyWindows(
ds, fReportTime, prevRes, winSizeInDays = 6L, nRecInDay = nRecInDay )
res <- cbind( resApply$winInfo, tmp <- do.call(
rbind, lapply(resApply$resFUN, "[[", 1L)))
nRecRes <- nrow(res)
expect_equal( res$dayStart, res$startDay )
expect_equal( res$dayEnd, res$endDay )
expect_equal( res$iRecStart, res$startRec )
expect_equal( res$iRecEnd, res$endRec )
#
expect_true( all(diff(res$startDay[-1]) == 2L)) # shifted starting day
# day boundary before startRec
expect_true( all((ds$iDayOfRec[res$startRec[-1]] -
ds$iDayOfRec[res$startRec[-1] - 1]) == 1L))
# day boundary after endRec
expect_true( all((ds$iDayOfRec[res$startRec[-nRecRes]] + 1 -
ds$iDayOfRec[res$startRec[-nRecRes]]) == 1L))
# prevRes accumulated
expect_equal( resApply$resFUN[[nRecRes]]$sumRes, nrow(res) )
})
test_that("simplifyApplyWindows",{
nRec <- nrow(dsNEE)
nRecInDay <- 10L
ds <- within(dsNEE, {
iRec <- 1:nRec
# specifying the day for each record assuming equidistand records
iDayOfRec <- ((c(1:nRec) - 1L) %/% nRecInDay) + 1L
})
fReportTimeSimple <- function(dss, winInfo, prevRes = list()){
nRecS <- nrow(dss)
c(
startRec = dss$iRec[1]
,endRec = dss$iRec[nRecS]
,startDay = dss$iDayOfRec[1]
,endDay = dss$iDayOfRec[nRecS]
)
}
fReportTimeSimple(ds,0)
# larger than reference window of 4 days
resApply <- REddyProc:::applyWindows(
ds, fReportTimeSimple, winSizeInDays = 6L, nRecInDay = nRecInDay )
res <- REddyProc:::simplifyApplyWindows(resApply)
nRecRes <- nrow(res)
expect_equal( res$dayStart, res$startDay )
expect_equal( res$dayEnd, res$endDay )
expect_equal( res$iRecStart, res$startRec )
expect_equal( res$iRecEnd, res$endRec )
#
# repeat with function returning a single-row data.frame
fReportTimeSimpleDs <- function(dss, winInfo, prevRes = list()){
nRecS <- nrow(dss)
data.frame(
startRec = dss$iRec[1]
,endRec = dss$iRec[nRecS]
,startDay = dss$iDayOfRec[1]
,endDay = dss$iDayOfRec[nRecS]
)
}
fReportTimeSimpleDs(ds,0)
# larger than reference window of 4 days
resApply <- REddyProc:::applyWindows(
ds, fReportTimeSimpleDs, winSizeInDays = 6L, nRecInDay = nRecInDay )
resDs <- REddyProc:::simplifyApplyWindows(resApply)
expect_equal(res, resDs )
})
test_that("estimating temperature sensitivity windows outputs are in accepted range",{
skip_if_not_installed("mlegp")
dss <- dsNEE[ dsNEE$Rg_f <= 0 & dsNEE$PotRad_NEW <= 0 &
as.POSIXlt(dsNEE$sDateTime)$mday %in% 1:12, ]
dss$NEE <- dss$NEE_f
dss <- dss[ order(dss$Temp), ]
res <- REddyProc:::simplifyApplyWindows(
REddyProc:::applyWindows(
dss, REddyProc:::partGLFitNightTempSensOneWindow
, prevRes = data.frame(E0 = NA)
, winSizeInDays = 12L
#,controlGLPart = controlGLPart
))
#res <- partGLEstimateTempSensInBounds(dss$NEE_f, dss$Temp+273.15)
expect_true( res$E0 >= 50 && res$E0 <= 400 )
expect_true( res$RRefFit > 0 )
})
test_that("partGLFitLRCWindows outputs are in accepted range",{
skip_if_not_installed("mlegp")
ds <- partGLExtractStandardData(dsNEE)
#
#yday <- as.POSIXlt(dsNEE$sDateTime)$yday
dsTempSens <- dsTempSens0 <- REddyProc:::partGLFitNightTimeTRespSens(
ds
, nRecInDay = 48L
, controlGLPart = partGLControl()
)
lrcFitter <- RectangularLRCFitter()
#lrcFitter <- NonrectangularLRCFitter()
resFits <- REddyProc:::partGLFitLRCWindows(
ds, nRecInDay = 48L, dsTempSens = dsTempSens
, controlGLPart = partGLControl(nBootUncertainty = 10L)
, lrcFitter = lrcFitter)
expect_true( var(resFits$k) > 0)
expect_true( all( sapply(resFits$resOpt, length) != 0 ))
.tmp.f <- function(){
# in order to replicate, use nBoot = 0L
resParms0 <- REddyProc:::partGLFitLRCWindows(
ds, nRecInDay = 48L, dsTempSens = dsTempSens
, controlGLPart = partGLControl(nBootUncertainty = 0L)
, lrcFitter = lrcFitter)
iTooMany <- which(resParms0$iRecStart >= nrow(ds) - 2*48)
if (length(iTooMany)) {
resParms0 <- slice(resParms0,-iTooMany)
}
# regression result for resLRCEx1
# resLRCEx1 <- resParms0
# resLRCEx1Nonrectangular <- resParms0
dput(resParms0)
}
# check the conditions of Lasslop10 Table A1
expect_true( all(resFits$E0 >= 50 & resFits$E0 <= 400) )
expect_true( all(resFits$RRef > 0) )
expect_true( all(resFits$alpha >= 0 ) )
# first value may be greater, due to previous estimate
expect_true( all(resFits$alpha[-1] < 0.22) )
expect_true( all(resFits$beta >= 0 & resFits$beta < 250) )
expect_true( all(ifelse(resFits$beta > 100, resFits$beta_sd < resFits$beta, TRUE) ))
expect_true( all(resFits$k >= 0) )
expect_true( !all(is.na(resFits$RRef_sd)))
expect_true( all(resFits$iMeanRec < nrow(ds)) )
expect_true( all(resFits$iCentralRec < nrow(ds)) )
})
test_that("partGLFitLRCWindows error on low uncertainty",{
skip_if_not_installed("mlegp")
# see LightResponseCurveFitter_optimLRCOnAdjustedPrior
ds <- partGLExtractStandardData(dsNEE)
ds$sdNEE[100:200] <- 0
dsTempSens <- dsTempSens0 <- REddyProc:::partGLFitNightTimeTRespSens(
ds
, nRecInDay = 48L
, controlGLPart = partGLControl()
)
lrcFitter <- RectangularLRCFitter()
#lrcFitter <- NonrectangularLRCFitter()
if (exists("resFits")) rm(resFits)
expect_error(
resFits <- REddyProc:::partGLFitLRCWindows(
ds, nRecInDay = 48L, dsTempSens = dsTempSens
, controlGLPart = partGLControl(nBootUncertainty = 10L)
, lrcFitter = lrcFitter)
,"zeros in uncertainty"
)
expect_true(!exists("resFits"))
})
test_that(".partGPAssociateSpecialRows correct next lines",{
skip_if_not_installed("mlegp")
expect_error(
res <- REddyProc:::.partGPAssociateSpecialRows(integer(0),9)
)
res <- REddyProc:::.partGPAssociateSpecialRows(c(3,6,7,9),12)
expect_true( all(res$wBefore + res$wAfter == 1))
# special rows
expect_equal( c(iRec = 3,iSpecialBefore = 1L,iSpecialAfter = 1L
, iBefore = 3, iAfter = 3, wBefore = 0.5, wAfter = 0.5)
, unlist(res[3,]))
expect_equal( c(iRec = 6,iSpecialBefore = 2L,iSpecialAfter = 2L
, iBefore = 6, iAfter = 6, wBefore = 0.5, wAfter = 0.5)
, unlist(res[6,]))
# first rows and last rows
expect_equal( c(iRec = 1,iSpecialBefore = 1L,iSpecialAfter = 1L
, iBefore = 3, iAfter = 3, wBefore = 0.5, wAfter = 0.5)
, unlist(res[1,]))
expect_equal( c(iRec = 12,iSpecialBefore = 4L,iSpecialAfter = 4L
, iBefore = 9, iAfter = 9, wBefore = 0.5, wAfter = 0.5)
, unlist(res[12,]))
# weights after and before special row 6
expect_equal( rep(3,2), unlist(res[4:5,"iBefore"]))
expect_equal( rep(6,2), unlist(res[4:5,"iAfter"]))
expect_equal( rep(1,2), unlist(res[4:5,"iSpecialBefore"]))
expect_equal( rep(2,2), unlist(res[4:5,"iSpecialAfter"]))
expect_equal( (2:1)/3, unlist(res[4:5,"wBefore"]))
expect_equal( (1:2)/3, unlist(res[4:5,"wAfter"]))
# test last row is special
res <- REddyProc:::.partGPAssociateSpecialRows(c(3,6,7,9),9)
})
testGLInterpolateFluxes <- function(lrcFitter, resEx ){
tmp <- REddyProc:::partGLInterpolateFluxes(
dsNEE$Rg_f, dsNEE$VPD_f, dsNEE$Temp, resEx, lrcFitter = lrcFitter)
expect_equal( nrow(dsNEE), nrow(tmp) )
.tmp.plot <- function(){
tmp$time <- dsNEE$sDateTime
plot( Reco ~ time, tmp)
plot( GPP ~ time, tmp)
}
expect_true( all(c("GPP","Reco") %in% names(tmp) ))
tmp <- REddyProc:::partGLInterpolateFluxes(
dsNEE$Rg_f, dsNEE$VPD_f, dsNEE$Temp, resEx
, controlGLPart = partGLControl(isSdPredComputed = TRUE)
, lrcFitter = lrcFitter)
expect_equal( nrow(dsNEE), nrow(tmp) )
expect_true( all(c("GPP","Reco") %in% names(tmp) ))
expect_true( all(c("sdGPP","sdReco") %in% names(tmp) ))
}
test_that("partGLInterpolateFluxes runs with rectangular LRCFitter",{
lrcFitter = RectangularLRCFitter(); resEx <- resLRCEx1
testGLInterpolateFluxes( lrcFitter, resEx)
})
test_that("partGLInterpolateFluxes runs with rectangular LRCFitter",{
lrcFitter = NonrectangularLRCFitter(); resEx <- resLRCEx1Nonrectangular
testGLInterpolateFluxes( lrcFitter, resEx)
})
#resLRCEx1
test_that("partitionNEEGL",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
resEx <- resLRCEx1
#DoY.V.n <- as.POSIXlt(dsNEE1$sDateTime)$yday + 1L
#Hour.V.n <- as.POSIXlt(dsNEE1$sDateTime)$hour + as.POSIXlt(dsNEE1$sDateTime)$min/60
#dsNEE1$PotRad_NEW <- fCalcPotRadiation(DoY.V.n, Hour.V.n, LatDeg = 45.0, LongDeg = 1, TimeZoneHour = 0 )
tmp <- partitionNEEGL( dsNEE1,RadVar = 'Rg_f')
#tmp <- partitionNEEGL( dsNEE1,RadVar = 'Rg_f', isVerbose = FALSE) #no messages
#tmp <- partitionNEEGL( dsNEE1,RadVar = 'Rg_f', controlGLPart = partGLControl(nBootUncertainty = 0L, isAssociateParmsToMeanOfValids = FALSE))
expect_equal( nrow(dsNEE1), nrow(tmp) )
#tmp[ is.finite(tmp$FP_beta), ] # note FP_dRecPar is not zero, because iCentralRec != iMeanRec
#iPar <- which(is.finite(tmp$FP_E0))
#plot( tmp$FP_beta[iPar] ~ dsNEE1$sDateTime[iPar],type = "l", ylim = range(c(tmp$FP_beta[iPar],tmp$FP_GPP2000[iPar]))); lines( tmp$FP_GPP2000[iPar] ~ dsNEE1$sDateTime[iPar], col = "red")
#
# now test with different suffix: u50
dsNEE2 <- dsNEE
names(dsNEE2)[ match(c("NEE_f", "NEE_fqc", "NEE_fsd"),names(dsNEE2))] <-
c("NEE_u50_f", "NEE_u50_fqc", "NEE_u50_fsd")
tmp <- partitionNEEGL(
dsNEE2, RadVar = 'Rg_f', suffix = "u50"
, controlGLPart = partGLControl(
nBootUncertainty = 0L, isAssociateParmsToMeanOfValids = FALSE) )
expect_equal( nrow(dsNEE1), nrow(tmp) )
expect_true( all(is.finite(tmp$GPP_DT_u50)))
expect_true( all(tmp$GPP_DT_u50 >= 0))
expect_true( all(tmp$GPP_DT_u50 < 250))
expect_true( all(tmp$Reco_DT_u50 < 10))
expect_true( all(tmp$Reco_DT_u50 > 0))
expect_true( all(tmp$Reco_DT_u50_SD > 0))
expect_true( all(tmp$GPP_DT_u50_SD >= 0))
expect_true( all(abs(diff(tmp$Reco_DT_u50)) < 0.6)) #smooth
# reporting good values at central records
# tmp[resEx$iCentralRec,]
iRowsOpt <- which(is.finite(tmp$FP_E0))
expect_true( length(iRowsOpt) == nrow(resEx) )
#expect_true( all((tmp$FP_alpha[resEx$iCentralRec] - resEx$alpha)[
#resEx$parms_out_range == 0L] < 1e-2) )
.tmp.plot <- function(){
tmp$time <- dsNEE1$sDateTime
plot( Reco_DT_u50 ~ time, tmp)
#plot( diff(Reco_DT_u50) ~ time[-1], tmp)
plot( GPP_DT_u50 ~ time, tmp)
plot( FP_RRef_Night ~ time, tmp)
plot( FP_RRef ~ FP_RRef_Night, tmp)
}
})
test_that("partitionNEEGL with Lasslop options",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
#ds <- data.frame( NEE = dsNEE1$NEE_f, sdNEE = dsNEE1$NEE_fsd, Rg = dsNEE1$Rg_f, VPD = dsNEE1$VPD_f, Temp = dsNEE1$Temp, isDay = dsNEE1$isDay, isNight = dsNEE$isNight )
resEx <- resLRCEx1
dsNEE2 <- dsNEE
partGLControl()
names(dsNEE2)[ match(c("NEE_f", "NEE_fqc", "NEE_fsd"),names(dsNEE2))] <-
c("NEE_u50_f", "NEE_u50_fqc", "NEE_u50_fsd")
dsNEE2$NEE_u50_fsd[24] <- NA
tmp <- partitionNEEGL( dsNEE2, RadVar = 'Rg_f', suffix = "u50"
, controlGLPart = partGLControlLasslopCompatible() )
tmp$sDateTime <- dsNEE2$sDateTime
tmp$iRec <- 1:nrow(tmp)
#subset(tmp, is.finite(FP_errorcode))
expect_equal( nrow(dsNEE1), nrow(tmp) )
expect_true( all(is.finite(tmp$GPP_DT_u50)))
expect_true( all(tmp$GPP_DT_u50 >= 0))
expect_true( all(tmp$GPP_DT_u50 < 250))
expect_true( all(tmp$Reco_DT_u50 < 10))
expect_true( all(tmp$Reco_DT_u50 > 0))
expect_true( all(tmp$Reco_DT_u50_SD > 0))
expect_true( all(tmp$GPP_DT_u50_SD >= 0))
expect_true( all(abs(diff(tmp$Reco_DT_u50)) < 0.6)) #smooth
# reporting good values at central records
# tmp[resEx$iCentralRec,]
#expect_true( sum( is.finite(tmp$FP_alpha) ) == nrow(resEx) )
# options yield large differences, hence do not compare to resEx1
#expect_true( all(abs(tmp$FP_alpha[resEx$iCentralRec] - resEx$alpha)[
#resEx$parms_out_range == 0L & is.finite(tmp$FP_alpha[resEx$iCentralRec])] < 0.05) )
#expect_true( all((is.na(tmp$FP_alpha[resLRCEx1$iFirstRec] - resLRCEx1$a)[
#resLRCEx1$parms_out_range != 0L])) )
expect_true( length(is.finite(tmp$FP_RRef_Night)) > 0 )
.tmp.plot <- function(){
tmp$time <- dsNEE1$sDateTime
plot( Reco_DT_u50 ~ time, tmp)
#plot( diff(Reco_DT_u50) ~ time[-1], tmp)
plot( GPP_DT_u50 ~ time, tmp)
plot( FP_RRef_Night ~ time, tmp)
plot( FP_RRef ~ FP_RRef_Night, tmp)
}
})
isTimeInTestPeriod <- function(sDateTime){
(sDateTime >= as.POSIXct("1998-06-03 00:00:00",tz = tzEx)) &
(sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx) |
sDateTime >= as.POSIXct("1998-06-05 22:00:00",tz = tzEx))
}
.tmp.f <- function(){
isTimeInTestPeriodNoTemp <- function(sDateTime){
# strange: 4 more valid NEE yiels in temperature estimate failing TODO inspect
(sDateTime >= as.POSIXct("1998-06-03 00:00:00",tz = tzEx)) &
(sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx) |
sDateTime >= as.POSIXct("1998-06-06 00:00:00",tz = tzEx))
}
dsNEE2 <- dsNEE
dsNEE2$NEE_fqc[ isTimeInTestPeriodNoTemp(dsNEE$sDateTime) ] <- 2L
#table(dsNEE1$NEE_fqc)
#plot( NEE_f ~ sDateTime, dsNEE1 )
ds2 <- partGLExtractStandardData(dsNEE2)
dsTempSens2 <- REddyProc:::partGLFitNightTimeTRespSens( ds
, nRecInDay = 48L
, controlGLPart = partGLControl()
)
}
test_that("partitionNEEGL sparse data",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
#flag all data except one day bad in order to associate the same
#data with several windows
dsNEE1$NEE_fqc[ isTimeInTestPeriod(dsNEE$sDateTime) ] <- 2L
#table(dsNEE1$NEE_fqc)
#plot( NEE_f ~ sDateTime, dsNEE1 )
ds <- partGLExtractStandardData(dsNEE1)
#
dsTempSens <- REddyProc:::partGLFitNightTimeTRespSens(
ds
, nRecInDay = 48L
, controlGLPart = partGLControl()
)
resLRC <- REddyProc:::partGLFitLRCWindows(
ds, dsTempSens = dsTempSens, lrcFitter = RectangularLRCFitter() )
expect_true( resLRC$iMeanRec[2] == resLRC$iMeanRec[3])
#
tmp <- partitionNEEGL( dsNEE1)
expect_equal( nrow(dsNEE1), nrow(tmp) )
#
expect_true( all(is.finite(tmp$GPP_DT)))
expect_true( all(tmp$GPP_DT >= 0))
expect_true( all(tmp$GPP_DT < 250))
expect_true( all(tmp$Reco_DT < 10))
expect_true( all(tmp$Reco_DT > 0))
expect_true( all(tmp$Reco_DT_SD > 0))
expect_true( all(tmp$GPP_DT_SD >= 0))
#TODO expect_true( all(abs(diff(tmp$Reco_DT)) < 0.6)) #smooth
# reporting good values at first row
expect_true( sum( is.finite(tmp$FP_alpha) ) == sum(is.finite(resLRC$alpha)) )
expect_true( all((is.na(tmp$FP_alpha[resLRC$iCentralRec] - resLRC$alpha)[
resLRC$parms_out_range != 0L & is.finite(tmp$FP_alpha[resLRC$iCentralRec])])) )
.tmp.plot <- function(){
tmp$time <- dsNEE1$sDateTime
plot( Reco_DT ~ time, tmp)
#plot( diff(Reco_DT_u50) ~ time[-1], tmp)
plot( GPP_DT ~ time, tmp)
}
})
test_that("partitionNEEGL isNeglectVPDEffect",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
#flag all VPD data except one day as bad to associate the same data with several windows
dsNEE1$VPD_fqc[ (dsNEE$sDateTime >= as.POSIXct("1998-06-03 00:00:00",tz = tzEx)) &
(dsNEE$sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx) |
dsNEE$sDateTime >= as.POSIXct("1998-06-06 00:00:00",tz = tzEx)) ] <- 2L
#plot(VPD_fqc ~ sDateTime, dsNEE1)
ctrl <- partGLControl(isFilterMeteoQualityFlag = TRUE, isNeglectVPDEffect = TRUE)
ds <- partGLExtractStandardData(dsNEE1, controlGLPart = ctrl)
#plot(VPD ~ sDateTime, ds)
#
dsTempSens <- REddyProc:::partGLFitNightTimeTRespSens( ds
, nRecInDay = 48L
, controlGLPart = ctrl
)
resLRC <- REddyProc:::partGLFitLRCWindows(
ds, dsTempSens = dsTempSens, lrcFitter = RectangularLRCFitter() )
expect_true( resLRC$iMeanRec[2] == resLRC$iMeanRec[3])
#resLRC$nValidRec
#
resLRC2 <- REddyProc:::partGLFitLRCWindows(
ds, dsTempSens = dsTempSens, lrcFitter = RectangularLRCFitter()
,controlGLPart = ctrl )
expect_true( resLRC2$iMeanRec[2] != resLRC2$iMeanRec[3])
# used all records despite missing VPD
expect_true( all( resLRC2$nValidRec >= 50L) )
#resLRC2
#
resLRC2C <- REddyProc:::partGLFitLRCWindows(
ds, dsTempSens = dsTempSens, lrcFitter = RectangularLRCFitterCVersion()
,controlGLPart = ctrl )
# used all records despite missing VPD
expect_true( all( resLRC2C$nValidRec >= 50L) )
# equal unless bootstrapped uncertainty and optimization result
iSd <- grep("_sd$|resOpt$",names(resLRC2))
expect_equal( as.data.frame(resLRC2[,-iSd]), as.data.frame(resLRC2C[,-iSd]) )
#
tmp <- partitionNEEGL( dsNEE1, controlGLPart = ctrl )
#tmp[ resLRC2$iCentralRec,]
expect_equal( tmp[ resLRC2$iCentralRec,"FP_alpha"], resLRC2$alpha )
expect_equal( nrow(dsNEE1), nrow(tmp) )
#
# finite prediction despite missing VPD
expect_true( all(is.finite(tmp$GPP_DT)))
expect_true( all(tmp$GPP_DT >= 0))
expect_true( all(tmp$GPP_DT < 250))
expect_true( all(tmp$Reco_DT < 10))
expect_true( all(tmp$Reco_DT > 0))
expect_true( all(tmp$Reco_DT_SD > 0))
expect_true( all(tmp$GPP_DT_SD >= 0))
})
test_that("partGLPartitionFluxes missing night time data",{
skip_if_not_installed("mlegp")
dsNEE1 <- dsNEE
#set all data except one day to NA to associate the same data
#with several windows
dsNEE1$hourOfDay <- as.POSIXlt(dsNEE1$sDateTime)$hour
dsNEE1$NEE_f[ (dsNEE1$sDateTime >= as.POSIXct("1998-06-01 00:00:00",tz = tzEx)) &
(dsNEE1$sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx)) &
((dsNEE1$hourOfDay >= 18) | (dsNEE1$hourOfDay <= 6))] <- NA
ds <- dsNEE1
ds$NEE <- ds$NEE_f
ds$sdNEE <- ds$NEE_fsd
ds$Temp <- ds$Tair_f
ds$VPD <- ds$VPD_f
ds$Rg <- ds$Rg_f
#
dsTempSens <- REddyProc:::partGLFitNightTimeTRespSens( ds
, nRecInDay = 48L
, controlGLPart = partGLControl()
, winSizeNight = 4L, winExtendSizes = c()
)
expect_true( all( is.finite(dsTempSens$RRef)) )
resLRC <- REddyProc:::partGLFitLRCWindows(
ds, lrcFitter = RectangularLRCFitter(), dsTempSens = dsTempSens )
expect_true( all( is.finite(resLRC$RRef_night)) )
})
test_that("partGLPartitionFluxes filter Meteo flag not enough without VPD",{
skip_if_not_installed("mlegp")
dsNEE1 <- dsNEE
# test setting VPD_fqc to other than zero, for omitting those rows
dsNEE1$VPD_fqc[ isTimeInTestPeriod(dsNEE1$sDateTime) ] <- 1L
#plot( VPD_fqc ~ sDateTime, dsNEE1 )
#plot( NEE_f ~ sDateTime, dsNEE1 )
ds <- dsNEE1
#ds$NEE <- ds$NEE_f
#ds$sdNEE <- ds$NEE_fsd
#ds$Temp <- ds$Tair_f
#ds$VPD <- ds$VPD_f
ds$Rg <- ds$Rg_f
#
.tmp.f <- function(){
tmp0 <- partitionNEEGL(
dsNEE1, controlGLPart = partGLControl(isFilterMeteoQualityFlag = TRUE) )
tmpPar0 <- tmp0[is.finite(tmp0$FP_beta),]
}
tmp <- partitionNEEGL(
ds, controlGLPart = partGLControl(isFilterMeteoQualityFlag = TRUE) )
expect_equal( nrow(ds), nrow(tmp) )
tmpPar <- tmp[is.finite(tmp$FP_RRef_Night),]
# estimate despite missing VPD
expect_equal( sum(is.finite(tmpPar[,"FP_beta"])), 4L )
expect_equal( tmpPar[4L,"FP_k"], 0 ) # k = 0 for missing VPD
#
# now set temperature to bad quality flag,
dsNEE1$Tair_fqc[ isTimeInTestPeriod(dsNEE1$sDateTime) ] <- 1L
#dsNEE1$VPD_fqc[ (dsNEE1$sDateTime > as.POSIXct("1998-06-03 00:00:00",tz = tzEx)) &
#(dsNEE1$sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx) |
#dsNEE1$sDateTime >= as.POSIXct("1998-06-06 00:00:00",tz = tzEx)) ] <- 1L
ds <- dsNEE1
ds$Rg <- ds$Rg_f
tmp <- partitionNEEGL(
ds, controlGLPart = partGLControl(isFilterMeteoQualityFlag = TRUE) )
expect_equal( nrow(ds), nrow(tmp) )
tmpPar <- tmp[is.finite(tmp$FP_RRef_Night),]
# one row less with parameter estimates
expect_equal( sum(is.finite(tmpPar[,"FP_beta"])), 3L )
})
test_that("partGLPartitionFluxes missing prediction VPD",{
skip_if_not_installed("mlegp")
dsNEE1 <- dsNEE
# test setting VPD_fqc to other than zero, for omitting those rows
dsNEE1$VPD_fqc[ (dsNEE1$sDateTime > as.POSIXct("1998-06-03 00:00:00",tz = tzEx)) &
(dsNEE1$sDateTime < as.POSIXct("1998-06-05 00:00:00",tz = tzEx) |
dsNEE1$sDateTime >= as.POSIXct("1998-06-06 00:00:00",tz = tzEx)) ] <- 1L
#plot( VPD_fqc ~ sDateTime, dsNEE1 )
#plot( NEE_f ~ sDateTime, dsNEE1 )
iMissingVPD <- c(10,400)
dsNEE1$VPD_f[iMissingVPD] <- NA
ds <- dsNEE1
#ds$NEE <- ds$NEE_f
#ds$sdNEE <- ds$NEE_fsd
#ds$Temp <- ds$Tair_f
#ds$VPD <- ds$VPD_f
ds$Rg <- ds$Rg_f
#
expect_warning(tmp <- partitionNEEGL(
ds, controlGLPart = partGLControl(
isFilterMeteoQualityFlag = TRUE
, isRefitMissingVPDWithNeglectVPDEffect = FALSE) ))
expect_equal( nrow(ds), nrow(tmp) )
tmpPar <- tmp[is.finite(tmp$FP_RRef_Night),]
# estimate despite missing VPD
expect_equal( sum(is.finite(tmpPar[,"FP_beta"])), 4L )
expect_equal( tmpPar[4L,"FP_k"], 0 ) # k = 0 for missing VPD
expect_true( is.na(tmp$GPP_DT[iMissingVPD[1]]) ) # default could not predict
# second: both parameter sets with k = 0 -> VPD not evaluated for prediction
expect_true( is.finite(tmp$GPP_DT[iMissingVPD[2]]) )
#
# repeat with (default) refitting with neglecting VPD
tmp2 <- partitionNEEGL(
ds, controlGLPart = partGLControl(
isFilterMeteoQualityFlag = TRUE
, isRefitMissingVPDWithNeglectVPDEffect = TRUE) )
expect_equal( tmp2[-iMissingVPD,"GPP_DT"], tmp[-iMissingVPD,"GPP_DT"] )
expect_equal( tmp2[-iMissingVPD,"Reco_DT"], tmp[-iMissingVPD,"Reco_DT"] )
expect_equal( tmp2[iMissingVPD[2],c("GPP_DT","Reco_DT")],
tmp2[iMissingVPD[2],c("GPP_DT","Reco_DT")])
# now also estimate for first missing VPD
expect_true( is.finite(tmp2$GPP_DT[iMissingVPD[1]]) )
})
test_that("partitionNEEGL fixed tempSens",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
#
ds <- dsNEE1
ds$NEE <- ds$NEE_f
ds$sdNEE <- ds$NEE_fsd
ds$Temp <- ds$Tair_f
ds$VPD <- ds$VPD_f
ds$Rg <- ds$Rg_f
#
dsTempSens0 <- REddyProc:::partGLFitNightTimeTRespSens( ds )
startRecs <- REddyProc:::getStartRecsOfWindows( nrow(ds) )
fixedTempSens <- data.frame(E0 = 80, sdE0 = 10 )
tmp <- partitionNEEGL(
dsNEE1, RadVar = "Rg_f", controlGLPart = partGLControl(
fixedTempSens = fixedTempSens) )
expect_equal( nrow(dsNEE1), nrow(tmp) )
#
expect_true( all(is.finite(tmp$GPP_DT)))
expect_true( all(tmp$GPP_DT >= 0))
expect_true( all(tmp$GPP_DT < 250))
expect_true( all(tmp$Reco_DT < 10))
expect_true( all(tmp$Reco_DT > 0))
expect_true( all(tmp$Reco_DT_SD > 0))
expect_true( all(tmp$GPP_DT_SD >= 0))
.tmp.plot <- function(){
tmp$time <- dsNEE1$sDateTime
plot( Reco_DT ~ time, tmp)
#plot( diff(Reco_DT_u50) ~ time[-1], tmp)
plot( GPP_DT ~ time, tmp)
}
})
test_that("partitionNEEGL: missing sdNEE",{
skip_if_not_installed("mlegp")
dsNEE1 <- dsNEE
#summary(dsNEE1)
dsNEE1$NEE_fsd <- NA_real_
#
tmp <- partitionNEEGL(
dsNEE1, RadVar = "Rg_f", controlGLPart = partGLControl() )
tmpWithPar <- tmp[ is.finite(tmp$FP_alpha),]
expect_equal( nrow(tmpWithPar), 4) # fitted with missing fsd
#
rm(tmp)
expect_error(
tmp <- partitionNEEGL(
dsNEE1, RadVar = "Rg_f"
, controlGLPart = partGLControl(replaceMissingSdNEEParms = c(NA,NA)) )
)
})
# see .profilePartGL in develop/profile/profile.R
test_that("partitionNEEGL long gap",{
skip_if_not_installed("mlegp")
skip("cannot repoduce SmoothTempSens failing")
dsNEE1 <- Example_DETha98_Filled
# introduce a long four month gap
bo <- dsNEE1$DoY >= 80 & dsNEE1$DoY <= 360
dsNEE1$NEE_f[bo] <- NA
dsNEE1$NEE_fqc[bo] <- 3
tmp <- partitionNEEGL( dsNEE1)
tmp2 <- cbind(select(dsNEE1, sDateTime), tmp)
plot( FP_qc ~ sDateTime, tmp2)
plot( GPP_DT ~ sDateTime, tmp2)
})
test_that("partitionNEETK",{
skip_if_not_installed("mlegp")
skip_on_cran()
dsNEE1 <- dsNEE
resEx <- resLRCEx1
#DoY.V.n <- as.POSIXlt(dsNEE1$sDateTime)$yday + 1L
#Hour.V.n <- as.POSIXlt(dsNEE1$sDateTime)$hour + as.POSIXlt(dsNEE1$sDateTime)$min/60
#dsNEE1$PotRad_NEW <- fCalcPotRadiation(DoY.V.n, Hour.V.n, LatDeg = 45.0, LongDeg = 1, TimeZoneHour = 0 )
tmp <- partitionNEEGL(
dsNEE1,RadVar = 'Rg_f'
, controlGLPart = partGLControl(useNightimeBasalRespiration = TRUE))
#tmp <- partitionNEEGL( dsNEE1,RadVar = 'Rg_f', controlGLPart = partGLControl(nBootUncertainty = 0L, isAssociateParmsToMeanOfValids = FALSE))
expect_equal( nrow(dsNEE1), nrow(tmp) )
#tmp[ is.finite(tmp$FP_beta), ] # note FP_dRecPar is not zero, because iCentralRec != iMeanRec
#iPar <- which(is.finite(tmp$FP_E0))
#plot( tmp$FP_beta[iPar] ~ dsNEE1$sDateTime[iPar],type = "l", ylim = range(c(tmp$FP_beta[iPar],tmp$FP_GPP2000[iPar]))); lines( tmp$FP_GPP2000[iPar] ~ dsNEE1$sDateTime[iPar], col = "red")
#
# now test with different suffix: u50
dsNEE2 <- dsNEE
names(dsNEE2)[ match(c("NEE_f", "NEE_fqc", "NEE_fsd"),names(dsNEE2))] <-
c("NEE_u50_f", "NEE_u50_fqc", "NEE_u50_fsd")
tmp <- partitionNEEGL(
dsNEE2, RadVar = 'Rg_f', suffix = "u50"
, controlGLPart = partGLControl(
nBootUncertainty = 0L, isAssociateParmsToMeanOfValids = FALSE
,useNightimeBasalRespiration = TRUE )
)
expect_equal( nrow(dsNEE1), nrow(tmp) )
expect_true( all(is.finite(tmp$GPP_DT_u50)))
expect_true( all(tmp$GPP_DT_u50 >= 0))
expect_true( all(tmp$GPP_DT_u50 < 250))
expect_true( all(tmp$Reco_DT_u50 < 10))
expect_true( all(tmp$Reco_DT_u50 > 0))
expect_true( all(tmp$Reco_DT_u50_SD >= 0))
expect_true( all(tmp$GPP_DT_u50_SD >= 0))
# jumps between daytime and nighttime
#expect_true( all(abs(diff(tmp$Reco_DT_u50)) < 0.6)) #smooth
# reporting good values at central records
# tmp[resEx$iCentralRec,]
iRowsOpt <- which(is.finite(tmp$FP_E0))
expect_true( length(iRowsOpt) == nrow(resEx) )
#expect_true( all((tmp$FP_alpha[resEx$iCentralRec] - resEx$alpha)[
#resEx$parms_out_range == 0L] < 1e-2) )
.tmp.plot <- function(){
tmp$time <- dsNEE1$sDateTime
plot( Reco_DT_u50 ~ time, tmp)
#plot( diff(Reco_DT_u50) ~ time[-1], tmp)
plot( GPP_DT_u50 ~ time, tmp)
plot( FP_RRef_Night ~ time, tmp)
plot( FP_RRef ~ FP_RRef_Night, tmp)
plot( Reco_DT_u50_SD ~ time, tmp)
}
})
test_that("report missing VPD_f column in error",{
skip("only interactively test issue #34")
DETha98 <- fConvertTimeToPosix(Example_DETha98, 'YDH', Year = 'Year',
Day = 'DoY', Hour = 'Hour')[-(2:4)]
EProc <- sEddyProc$new('DE-Tha', DETha98,
c('NEE', 'Rg', 'Tair', 'VPD', 'Ustar'))
EProc$sMDSGapFillAfterUstar('NEE', uStarTh = 0.3, FillAll = TRUE)
# missing 'VPD'...
for (i in c('Tair', 'Rg')) EProc$sMDSGapFill(i, FillAll = TRUE)
EProc$sSetLocationInfo(LatDeg = 51.0, LongDeg = 13.6, TimeZoneHour = 1)
# ...produces Error here
# Error should report missing column VPD_f
expect_error(
EProc$sGLFluxPartition(suffix = "uStar")
,"VPD_f"
)
})
test_that("no nighttime data",{
skip_if_not_installed("mlegp")
skip_on_cran()
cols <- c('NEE','Rg','Tair','VPD', 'Ustar'
,"Tair_f","VPD_f","Rg_f","NEE_f"
,"NEE_fqc", "Tair_fqc", "Rg_fqc", "NEE_fsd")
#ds <- Example_DETha98_Filled[,c('DateTime',cols)]
ds <- subset(Example_DETha98_Filled
, sDateTime >= as.POSIXct("1998-05-01 00:15:00",tz = tzEx) &
sDateTime <= as.POSIXct("1998-09-01 21:45:00",tz = tzEx))
ds <- ds[,c('DateTime',cols)]
ds$NEE_f[ds$Rg_f <= 10] <- NA # simulate having no night-time records
EProc <- suppressWarnings(sEddyProc$new('DE-ThaSummer', ds, cols))
EProc$sSetLocationInfo(LatDeg = 51.0, LongDeg = 13.6, TimeZoneHour = 1)
EProc$sGLFluxPartition(controlGLPart = partGLControl(fixedTempSens = data.frame(
E0 = 150, sdE0 = 50, RRef = 20)))
})
.test_sEddyProc_sGLFluxPartitionUStarScens <- function(){}
test_that("sEddyProc_sGLFluxPartitionUStarScens wrong suffix",{
skip_if_not_installed("mlegp")
dsTest <- Example_DETha98_Filled %>% mutate(
NEE_uStar_f = NEE, NEE_uStar_fqc = NEE_fqc, NEE_uStar_fsd = NEE_fsd,
NEE_U50_f = NEE, NEE_U50_fqc = NEE_fqc, NEE_U50_fsd = NEE_fsd
) %>%
select(!"sDateTime")
EProc <- sEddyProc$new(
'DE-Tha', dsTest, setdiff(names(dsTest), c("NEE_fnum","NEE_fwin")))
EProc$sSetUStarSeasons(1L)
EProc$sSetUstarScenarios(data.frame(season = 1L, uStar = 0.28, U50 = 0.38))
EProc$sSetLocationInfo(LatDeg = 51.0, LongDeg = 13.6, TimeZoneHour = 1)
expect_error(
EProc$sGLFluxPartitionUStarScens(uStarScenKeep = "nonexisting")
,"nonexisting")
})
test_that("sEddyProc_sGLFluxPartitionUStarScens",{
skip_if_not_installed("mlegp")
dsTest <- Example_DETha98_Filled %>% mutate(
NEE_uStar_f = NEE, NEE_uStar_fqc = NEE_fqc, NEE_uStar_fsd = NEE_fsd,
NEE_U50_f = NEE, NEE_U50_fqc = NEE_fqc, NEE_U50_fsd = NEE_fsd
) %>%
select(!"sDateTime")
EProc <- sEddyProc$new(
'DE-Tha', dsTest, setdiff(names(dsTest), c("NEE_fnum","NEE_fwin")))
EProc$sSetUStarSeasons(1L)
EProc$sSetUstarScenarios(data.frame(season = 1L, uStar = 0.28, U50 = 0.38))
EProc$sSetLocationInfo(LatDeg = 51.0, LongDeg = 13.6, TimeZoneHour = 1)
#EProc$sMDSGapFill('Tair', FillAll = FALSE, minNWarnRunLength = NA)
#EProc$sMDSGapFill('VPD', FillAll = FALSE, minNWarnRunLength = NA)
#EProc$trace(sGLFluxPartitionUStarScens, recover); # EProc$untrace(sGLFluxPartitionUStarScens)
#EProc$sGLFluxPartitionUStarScens(uStarScenKeep = "U50")
EProc$sGetUstarScenarios()
EProc$sGetUstarSuffixes()
EProc$sGLFluxPartitionUStarScens(warnOnOtherErrors=TRUE, uStarScenKeep="U50")
expect_true(all(c("GPP_DT_uStar","GPP_DT_U50", "GPP_DT_U50_SD") %in%
names(EProc$sTEMP)))
})
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.