R Under development (unstable) (2026-06-25 r90191 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2026 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> ## Testing apsim_met suite of functions
> ## Some observations:
> ##
> ## 1. check should never fail. It should always return warnings
> ## 2. napad should work for many different tricky situations
> ## 3. plotting and different combinations
> 
> require(apsimx)
Loading required package: apsimx
APSIM(X) not found.
                        If APSIM(X) is installed in an alternative location,
                        set paths manually using 'apsimx_options' or 'apsim_options'.
                        You can still try as the package will look into the registry (under Windows)
> require(mgcv)
Loading required package: mgcv
Loading required package: nlme
This is mgcv 1.9-4. For overview type '?mgcv'.
> require(ggplot2)
Loading required package: ggplot2
> packageVersion("apsimx")
[1] '2.8.270'
> packageVersion("GSODR")
[1] '5.0.1'
> packageVersion("data.table")
[1] '1.18.4'
> 
> run.apsim.met <- get(".run.local.tests", envir = apsimx.options)
> 
> if(.Platform$OS.type == "unix"){
+   internet <- !is.null(nsl("google.com"))
+ }else{
+   internet <- FALSE  
+ }
> run.apsim.met <- FALSE  ## This does not work because GSODR does not work
> 
> if(run.apsim.met && internet){
+   
+   ## Testing napad and impute
+   ## Buenos Aires 34d 36 12, 58 22 54 W
+   bsas.lat <- unit_conv("34d 36' 12\" W", from = "degrees", to = "decimal")
+   bsas.lon <- unit_conv("58d 22' 54\" S", from = "degrees", to = "decimal")
+   
+   pwr <- get_power_apsim_met(lonlat = c(bsas.lon, bsas.lat), dates = c("2010-01-01", "2015-12-31"))
+   dim(pwr)
+   check_apsim_met(pwr)
+   gsd <- get_gsod_apsim_met(lonlat = c(bsas.lon, bsas.lat), 
+                             dates = c("2010-01-01", "2015-12-31"))
+   dim(gsd)
+   check_apsim_met(gsd)
+   gsd$radn <- pwr$radn
+   check_apsim_met(gsd)
+   gsd.impt <- impute_apsim_met(gsd)
+   check_apsim_met(gsd.impt)
+   
+   ## Test for Brazil (Sao Paulo)
+   ## 23 33′S 46 38′W
+   sp.lat <- unit_conv("23d 33' 0\" W", from = "degrees", to = "decimal")
+   sp.lon <- unit_conv("46d 38' 0\" S", from = "degrees", to = "decimal")
+   
+   pwr <- get_power_apsim_met(lonlat = c(sp.lon, sp.lat), dates = c("2010-01-01", "2015-12-31"))
+   check_apsim_met(pwr)
+   
+   plot(pwr, met.var = "rain", years = 2010:2015, cumulative = TRUE, climatology = TRUE)
+   
+   gsd <- get_gsod_apsim_met(lonlat = c(sp.lon, sp.lat), dates = c("2010-01-01", "2015-12-31"), fillin.radn = TRUE)
+   check_apsim_met(gsd)
+   gsd.np <- napad_apsim_met(gsd)
+   gsd.imp <- impute_apsim_met(gsd.np)
+   check_apsim_met(gsd.imp) ## Rain variable is missing
+   plot(gsd.imp, met.var = "radn", years = 2010:2015, cumulative = TRUE, climatology = TRUE)
+   
+   ## Test for Daymet data
+   dmt <- get_daymet2_apsim_met(lonlat = c(-93, 42), years = c(2000, 2020))
+   check_apsim_met(dmt)
+   dmt.np <- napad_apsim_met(dmt)  
+   dmt.imp <- impute_apsim_met(dmt.np)
+   check_apsim_met(dmt.imp)
+   
+   summary(dmt.imp)
+   plot(dmt.imp)
+   plot(dmt.imp, met.var = "rain", years = 2010:2015, cumulative = TRUE, climatology = TRUE)
+ }
> 
> ## Testing the feature for adding a column to a met file
> if(run.apsim.met){
+   
+   ## Testing adding a column
+   extd.dir <- system.file("extdata", package = "apsimx")
+   ames <- read_apsim_met("Ames.met", src.dir = extd.dir)
+   
+   ## ames$vp <- 5 This does not work and it is not supposed to 
+   
+   vp <- data.frame(vp = abs(rnorm(nrow(ames), 2)))
+   attr(vp, "units") <- "(hPa)"
+   
+   ames$vp <- vp ## This works now (2026-6-18)
+   ames[["vp"]] <- vp ## This works now (2026-6-18)
+ 
+   val <- abs(rnorm(nrow(ames), 2))
+   nm <- "vp"
+   ames <- add_column_apsim_met(ames, value = val, name = "vp", units = "(hPa)")
+   
+   ## Extra tests
+   ames <- ames[ames$year < 2018, ]
+   ames <- apsimx:::pp_apsim_met(ames)
+   ames <- tt_apsim_met(ames, dates = c("01-05", "30-10"), method = "Classic_TT")
+   
+   plot(ames, met.var = "photoperiod")
+   plot(ames, met.var = "Classic_TT")
+     
+ }
> 
> if(run.apsim.met && internet){
+  
+   iem <- get_iem_apsim_met(lonlat = c(-93.77, 42.02), 
+                            dates = c("2013-01-01","2013-12-31"))
+   
+   pwr <- get_power_apsim_met(lonlat = c(-93.77, 42.02), 
+                            dates = c("2013-01-01","2013-12-31"))
+   
+   dmt <- get_daymet_apsim_met(lonlat = c(-93.77, 42.02), 
+                               years = 2013)
+   
+   cmp <- compare_apsim_met(iem, pwr[, 1:6], dmt[, 1:6], 
+                            labels = c("IEM", "POWER", "DAYMET"))
+   
+   cmp <- compare_apsim_met(iem, pwr[, 1:6], dmt[, 1:6], 
+                            met.var = "rain",
+                            labels = c("IEM", "POWER", "DAYMET"))
+ 
+   plot(cmp)
+   plot(cmp, plot.type = "ts")
+   plot(cmp, plot.type = "ts", met.var = "rain", cumulative = TRUE)
+ 
+ }
> 
> #### Different graph tests ----
> if(run.apsim.met && internet){
+   
+   ## Testing the graph functionality
+   pwr <- get_power_apsim_met(lonlat = c(-93.77, 42.02), 
+                              dates = c("1990-01-01","2020-12-31"))
+   
+   ### default plotting options
+   plot(pwr)
+   plot(pwr, climatology = TRUE)
+   plot(pwr, climatology = TRUE, years = 2015:2018)
+   plot(pwr, climatology = TRUE, cumulative = TRUE)
+   plot(pwr, climatology = TRUE, years = 2015:2018, cumulative = TRUE)
+   
+   plot(pwr, summary = TRUE)
+   plot(pwr, summary = TRUE, climatology = TRUE)
+   plot(pwr, summary = TRUE, met.var = "rain", climatology = TRUE)
+   
+   ## Plotting using the density option ----
+   plot(pwr, plot.type = "density")
+   plot(pwr, plot.type = "density", met.var = "radn")
+   plot(pwr, plot.type = "density", met.var = "rain")
+   
+   ## Plotting using the density option
+   plot(pwr, plot.type = "density", climatology = TRUE)
+   plot(pwr, plot.type = "density", met.var = "radn", climatology = TRUE)
+   plot(pwr, plot.type = "density", met.var = "rain", climatology = TRUE)
+   
+   plot(pwr, plot.type = "density", climatology = TRUE, years = c(2012, 2019))
+   plot(pwr, plot.type = "density", met.var = "radn", climatology = TRUE, years = c(2012, 2019))
+   plot(pwr, plot.type = "density", met.var = "rain", climatology = TRUE, years = c(2012, 2019))
+ 
+   plot(pwr, plot.type = "density", summary = TRUE)
+   plot(pwr, plot.type = "density", met.var = "radn_sum", summary = TRUE)
+   plot(pwr, plot.type = "density", met.var = "rain", summary = TRUE)
+   
+   ## Testing the options with individual years
+   plot(pwr, plot.type = "density", summary = TRUE, years = c(2012:2015))
+   plot(pwr, plot.type = "density", met.var = "radn", summary = TRUE, years = c(2012, 2019))
+   plot(pwr, plot.type = "density", met.var = "rain", summary = TRUE, years = c(2015, 2017))
+   
+   ## Testing with climatology and years
+   plot(pwr, plot.type = "density", met.var = "rain", summary = TRUE, climatology = TRUE, years = c(2012, 2019))
+   plot(pwr, plot.type = "density", met.var = "radn_sum", summary = TRUE, climatology = TRUE, years = c(2012, 2019))
+   plot(pwr, plot.type = "density", met.var = "avg_maxt", summary = TRUE, climatology = TRUE, years = c(2012, 2019))
+ 
+ }
> 
> 
> #### Testing anomaly features ----
> if(run.apsim.met && internet){
+   
+   ## Testing the graph functionality
+   pwr <- get_power_apsim_met(lonlat = c(-93.77, 42.02), 
+                              dates = c("1990-01-01","2020-12-31"))
+   
+   apwr1 <- summary(pwr, anomaly = TRUE)
+   apwr2 <- summary(pwr, anomaly = "rain_sum")
+   apwr3 <- summary(pwr, anomaly = "rain")
+   apwr4 <- summary(pwr, anomaly = "maxt")
+   apwr5 <- summary(pwr, anomaly = c("maxt", "rain"))
+   
+   plot(pwr, plot.type = "anomaly")
+   plot(pwr, plot.type = "anomaly", summary = TRUE, met.var = "rain")
+   plot(pwr, plot.type = "anomaly", summary = TRUE, years = 2012:2015)
+   
+   plot(pwr, plot.type = "anomaly", summary = TRUE, met.var = c("avg_maxt", "rain_sum"))
+   plot(pwr, plot.type = "anomaly", summary = TRUE, years = 2012:2015, met.var = c("avg_maxt", "rain_sum"))
+   plot(pwr, plot.type = "anomaly", summary = TRUE, years = 2012:2015, met.var = c("rain_sum", "avg_maxt"))
+   
+   pwr <- tt_apsim_met(pwr, dates = c("01-05", "30-10"), method = "Classic_TT")
+   
+   summary(pwr)
+   summary(pwr, compute.frost = TRUE)
+   
+   plot(pwr, met.var = "Classic_TT")
+   plot(pwr, met.var = "Classic_TT", plot.type = "density")
+   plot(pwr, met.var = "Classic_TT", plot.type = "anomaly", years = 2012:2019)
+   
+   plot(pwr, plot.type = "anomaly", summary = TRUE, met.var = c("rain", "high_classic_tt"))
+   plot(pwr, plot.type = "anomaly", summary = TRUE, met.var = c("rain", "high_classic_tt"),
+        years = 2012:2016)
+ }
> 
> if(run.apsim.met && internet){
+   
+   #### Will use this code to identify the potential issue with 
+   #### modifying a iem.met object
+   iem.met <- get_iem_apsim_met(lonlat = c(-93.77, 42.02), 
+                                 dates = c("2000-01-01","2020-12-31"))
+   
+   iem.met <- tt_apsim_met(iem.met, method = "Classic_TT", dates = c("1-5", "31-10"))
+   
+   names(iem.met)
+   ### This works as expected now and it removes the units properly
+   iem.met$Classic_TT <- NULL 
+     
+   iem.met <- tt_apsim_met(iem.met, method = "Classic_TT", dates = c("1-5", "31-10"))
+     
+   iem.met <- remove_column_apsim_met(iem.met, "Classic_TT")
+ 
+   ### This should return an error  
+   try(iem.met$hola <- 100, silent = TRUE)
+   
+ }
> 
> ### Next test should be for creating the 'by' argument
> 
> #### Comparing met files with different year ranges ----
> 
> if(run.apsim.met){
+   
+   pwr1 <- get_power_apsim_met(lonlat = c(-93, 42), dates = c("2000-01-01", "2010-12-31"))
+   pwr2 <- get_power_apsim_met(lonlat = c(-93, 42), dates = c("2011-01-01", "2024-12-31"))
+   
+   ### This has a better error message
+   merge.test1 <- try(compare_apsim_met(pwr1, pwr2), silent = TRUE)
+   
+   ### What about one year in common?
+   
+   pwr3 <- get_power_apsim_met(lonlat = c(-93, 42), dates = c("2000-01-01", "2011-12-31"))
+   
+   ### This should work but it should give a message
+   merge.test2 <- compare_apsim_met(pwr2, pwr3)
+   
+   dmt <- get_daymet_apsim_met(lonlat = c(-93, 42), years = c(2000, 2023))
+   
+   merge.test3 <- try(compare_apsim_met(pwr3, dmt, labels = c("POWER", "DAYMET")), silent = TRUE)
+   ### The number of columns is the same
+   ### This gives an error because some column names are different
+   
+   merge.test4 <- compare_apsim_met(pwr3[, 1:6], dmt[, 1:6], labels = c("POWER", "DAYMET"))
+   plot(merge.test4, met.var = "rain", plot.type = "ts", cumulative = TRUE)
+   plot(merge.test4, plot.type = "ts")
+   
+   #### How can I fix the previous gam fit?
+   mm4 <- merge.test4$met.mrg
+   
+   ggplot(data = mm4, aes(x = dates, y = radn.1)) + 
+     geom_point() + 
+     geom_smooth(formula = y ~ s(x, bs = "cs", k = 70, fx = TRUE))
+   
+   
+   merge.test5 <- compare_apsim_met(pwr3[, 1:6], dmt[, 1:6], labels = c("POWER", "DAYMET"),
+                                    check = FALSE)
+   
+   plot(merge.test5, met.var = "rain", plot.type = "ts", cumulative = TRUE)
+   
+   ### Scramble the order of columns
+   dmt.attr <- attributes(dmt)
+   dmt2 <- subset(dmt, select = names(dmt)[8:1])
+   attr(dmt2, "filename") <- dmt.attr$filename
+   attr(dmt2, "site") <- dmt.attr$site
+   attr(dmt2, "latitude") <- dmt.attr$latitude
+   attr(dmt2, "longitude") <- dmt.attr$longitude
+   attr(dmt2, "tav") <- dmt.attr$tav
+   attr(dmt2, "colnames") <- rev(dmt.attr$colnames)
+   attr(dmt2, "units") <- dmt.attr$units
+   attr(dmt2, "comments") <- dmt.attr$comments
+   attr(dmt2, "class") <- dmt.attr$class
+   attr(dmt2, "amp") <- dmt.attr$amp
+   head(dmt2)
+ 
+   merge.test6 <- try(compare_apsim_met(pwr3, dmt2, labels = c("POWER", "DAYMET")), silent = TRUE)
+   ### this results in an error as it should
+ 
+   merge.test7 <- compare_apsim_met(pwr3[, 1:6], dmt2[, 3:8], labels = c("POWER", "DAYMET"))  
+   plot(merge.test7)
+   plot(merge.test7, met.var = "rain", plot.type = "ts", cumulative = TRUE)
+   
+   summary(pwr3)
+   summary(dmt2)
+   
+   plot(pwr3, summary = TRUE, climatology = TRUE)
+   plot(dmt2, summary = TRUE, climatology = TRUE)
+ 
+   plot(pwr3, met.var = "rain", summary = TRUE, climatology = TRUE)
+   plot(dmt2, met.var = "rain", summary = TRUE, climatology = TRUE)
+   
+   ### Testing the extractor operation, this works now
+   ### DO I need to create other methods?
+   ### I have created `[<-.met` and `[[<-.met`
+   ### which basically send an error message and ask instead to use 'add_column_apsim_met'
+   dmt5 <- dmt[, 1:4]
+   class(dmt5)
+   dmt6 <- dmt[1:365, 1:4]
+   class(dmt6)
+   
+   summary(dmt6)
+   
+ }
> 
> if(FALSE){
+   
+   ### How the heck does the 'Extract' work in 'met' files?
+   ## Testing adding a column
+   extd.dir <- system.file("extdata", package = "apsimx")
+   ames <- read_apsim_met("Ames.met", src.dir = extd.dir)
+   
+   ## ames$vp <- 5 This does not work and it is not supposed to 
+   
+   vp <- data.frame(vp = abs(rnorm(nrow(ames), 2)))
+   attr(vp, "units") <- "(hPa)"
+   
+   ames$vp <- vp ## This works now (2026-6-18)
+   ames[["vp"]] <- vp ## This works now (2026-6-18)
+   
+   
+ }
> 
> proc.time()
   user  system elapsed 
   2.17    0.29    2.45