Error of `phenofit-MOD13A2_HeNan.R`

liuh1995 commented 1 year ago

In a Windows platform

phenofit-MOD13A2_HeNan.txt

Error in { : task 1 failed - "could not find function "runningId"" Timing stopped at: 67.56 11.98 85.7

next, i changed the code

Is NULL normal?

I want to konw how to solve it. How long does the code run once. Thank you!

kongdd commented 1 year ago

如注释中所写，并行请用linux或wsl，windows并行比较复杂。如果一定是windows系统，请自行搜索R语言如何并行，以及并行如何共享变量和函数。 %dopar%改写成%do%可关掉并行模式。
关于计算时间，phenofit效率没有这么高。

PS:

Windows terminal+WSL搭建Linux开发环境, https://zhuanlan.zhihu.com/p/209032586

liuh1995 commented 1 year ago

以后同样的问题写在一个issue里，谢谢老师答复，NULL是正常的吧

kongdd commented 1 year ago

bug已修复，详见https://github.com/eco-hydro/phenofit-scripts/blob/master/multiGS_MODIS_Henan/phenofit-MOD13A2_HeNan.md，另外添加了并行的示例

liuh1995 commented 1 year ago

谢谢老师，已解决！

cymkG commented 1 year ago

source("C:/UserData/lium/Downloads/phenofit-scripts-master/main_pkgs.R", encoding = "UTF-8")
#> Warning: package 'phenofit' was built under R version 4.2.1
#> Registered S3 method overwritten by 'Ipaper':
#>   method           from      
#>   print.data.table data.table
#> Warning: package 'magrittr' was built under R version 4.2.1
#> 
#> Attaching package: 'Ipaper'
#> The following object is masked from 'package:phenofit':
#> 
#>     melt_list
#> 
#> Attaching package: 'sf2'
#> The following object is masked from 'package:graphics':
#> 
#>     plot
#> The following object is masked from 'package:base':
#> 
#>     plot
#> Warning: package 'ggplot2' was built under R version 4.2.2
#> Warning: package 'ggnewscale' was built under R version 4.2.1
#> Loading required package: timechange
#> Warning: package 'timechange' was built under R version 4.2.2
#> 
#> Attaching package: 'lubridate'
#> The following objects are masked from 'package:base':
#> 
#>     date, intersect, setdiff, union
#> Warning: package 'zeallot' was built under R version 4.2.1
#> Warning: package 'stringr' was built under R version 4.2.2
#> Warning: package 'dplyr' was built under R version 4.2.2
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union
#> Warning: package 'job' was built under R version 4.2.1
#> Warning: package 'sp' was built under R version 4.2.1
#> Warning: package 'data.table' was built under R version 4.2.2
#> 
#> Attaching package: 'data.table'
#> The following objects are masked from 'package:dplyr':
#> 
#>     between, first, last
#> The following objects are masked from 'package:lubridate':
#> 
#>     hour, isoweek, mday, minute, month, quarter, second, wday, week,
#>     yday, year
#> The following objects are masked from 'package:Ipaper':
#> 
#>     first, last, transpose
#> Warning: package 'rcolors' was built under R version 4.2.1
#> Registered S3 method overwritten by 'lattice.layers':
#>   method    from        
#>   +.trellis latticeExtra
#> 
#> Attaching package: 'lattice.layers'
#> The following object is masked from 'package:rcolors':
#> 
#>     get_color
#> The following object is masked from 'package:Ipaper':
#> 
#>     dev_off
#> The following objects are masked from 'package:phenofit':
#> 
#>     get_options, set_options
#> Warning: package 'terra' was built under R version 4.2.1
#> terra 1.6.17
#> 
#> Attaching package: 'terra'
#> The following objects are masked from 'package:lattice.layers':
#> 
#>     rotate, wrap
#> The following object is masked from 'package:data.table':
#> 
#>     shift
#> The following object is masked from 'package:grid':
#> 
#>     depth
#> The following object is masked from 'package:Ipaper':
#> 
#>     clamp
#> The following objects are masked from 'package:magrittr':
#> 
#>     extract, inset

# infile = "C:/UserData/lium/Downloads/phenofit-scripts-master/data/MOD13A2_Henan_2015_2020.rda"
infile = "C:/UserData/lium/Downloads/phenofit-scripts-master/MEE2022_phenofit/MOD13A2_Henan_2015_2020.rda"

if (!file.exists(infile)) {
    indir = path.mnt("C:/UserData/lium/Downloads/phenofit-scripts-master/MEE2022_phenofit/phenofit_MEE2022_INPUTs_NorthChina")
    file_vi  = glue("{indir}/MOD13A2_Henan_2016_2020_EVI.tif")
    file_qc  = glue("{indir}/MOD13A2_Henan_2016_2020_SummaryQA.tif")
    file_doy = glue("{indir}/MOD13A2_Henan_2016_2020_DayOfYear.tif")

    # dates = terra::rast(infile) %>% names() %>%
    #     substr(1, 10) %>% as.Date("%Y_%m_%d")
    dates = modis_date("2015-01-01", "2020-12-31", 16)
    l = read_rast(file_vi) %>% rast2mat()
    l_qc = read_rast(file_qc) %>% rast2mat()
    l_doy = read_rast(file_doy) %>% rast2mat()

    I_grid = l$I_grid
    d_coord = l$d_coord[I_grid, ]
    data <- listk(
        VI = l$mat[I_grid, ] / 1e4,
        QC = l_qc$mat[I_grid, ],
        DOY = l_doy$mat[I_grid, ],
        dates = l$dates, qcFUN = qc_summary)
    save(data, d_coord, file = infile)
} else {
    load(infile)
    n = nrow(data$QC)
}

## Parameter setting
# If provided, LCs should have the same length as `data$VI`
LCs <- NULL
nptperyear <- 23
opt_old <- phenofit::get_options()

# Different parameters for different land covers (LCs) refered by TIMESAT
# Name of `list_options` corresponds to `LCs` code
list_options <- list(
    "1" = list(),
    "2" = list(),
    default = list(
        wFUN = wTSM, wmin = 0.2,
        verbose = FALSE,
        season = list(
            # rFUN = "smooth_wWHIT", lambda = 0.5,
            rFUN = "smooth_wHANTS", nf = 6,
            maxExtendMonth = 12, r_min = 0.0),
        # fine fitting parameters
        fitting = list(nextend = 1, minExtendMonth = 0, maxExtendMonth = 0.5,
            methods = c("AG", "Zhang", "Beck", "Elmore", "Gu")[3], #,"klos",, 'Gu'
            iters = 1,
            minPercValid = 0)
    )
)
do.call(phenofit::set_options, list_options$default)

# # 1. check the performance at random sampled points
# # only run in debug mode
# if (0) {
#     source("scripts/main_pkgs.R", encoding = "UTF-8")
#     set.seed(1)
#     inds = sample(1:nrow(data$VI), 10) %>% sort()
#     # id_bads = c(19564, 25769, 68760, 79072, 79073, 83378, 83379, 91274, 99172, 99173,
#     #             113144, 119239, 125169, 125810, 132338, 136962, 136963, 138684)
#     # inds = c(1049, 1270, 1271, 1273, 1391, 1511, 1512, 1513, 1514, 1642, 1644,
#     #                 1654, 2240, 2614, 2645, 2649, 2650, 2651, 2652, 2825)
#     # i <- inds[1]
#     ofile = "test_v0.3.5_Henan_MODIS.pdf"
#     dev_open(ofile, 10, 3, use.cairo_pdf = FALSE)
#     # library(proffer)
#     do.call(phenofit::set_options, list_options$default)
# 
#     # proffer::pprof({
#         res = foreach(i = inds, k = icount(10)) %do% {
#             runningId(k, 1)
#             d = get_input(i, data)
#             tryCatch({
#                 r <- phenofit_point(d, plot = T, period = c(2015, 2020))#$pheno
#                 if (k < length(inds)) grid.newpage()
#             }, error = function(e) {
#                 message(sprintf('%s', e$message))
#             })
#         }
#     # })
#     dev_off()
#     Ipaper::pdf_view(ofile)
#     # set_options(opt_old)
# }

# If need parallel mode, should be `linux` or `wsl` system
# library(Ipaper)
# InitCluster(10, kill = FALSE)
# inds = 1:n
# t = system.time({
#     res = foreach(i = inds %>% set_names(., .), icount()) %dopar% {
#         runningId(i, 1000)
#         if (!is.null(LCs)) {
#             LC = LCs[i]
#             do.call(set_options, list_options[[LC]])
#         }
#         tryCatch({
#             d = get_input(i, data)
#             r <- phenofit_point(d, plot = TRUE, period = c(2015, 2020))$pheno
#             phenofit_point(d)$pheno
#         }, error = function(e) {
#             message(sprintf('%s', e$message))
#         })
#     }
# })

main_phenofit <- function(n_run = 10, step = 1000, outfile = NULL, .parallel = FALSE) {
  `%dof%` <- ifelse(.parallel, foreach::`%dopar%`, foreach::`%do%`)
  inds <- 1:n

  t <- system.time({
    res <- foreach(i = inds %>% set_names(., .), icount(n_run)) %dof% {
      runningId(i, step)

      if (!is.null(LCs)) {
        LC <- LCs[i]
        do.call(set_options, list_options[[LC]])
      }
      tryCatch({
        d <- get_input(i, data)
        r <- phenofit_point(d, plot = FALSE, period = c(2015, 2020))$pheno
      }, error = function(e) {
        message(sprintf('%s', e$message))
      })
    }
  })

  if (is.null(outfile)) {
    res
  } else save(res, t, file = outfile)
}

# save(res, t, file = "C:/UserData/lium/Downloads/phenofit-scripts-master/data/pheno_Henan_MODIS_V6.1.rda")

res = main_phenofit(10, 1)
#> [] running 1 ...
#> <environment: 0x0000019b93a085c8> 
#> NULL
#> non-numeric argument to binary operator
#> [] running 2 ...
#> <environment: 0x0000019b93a085c8> 
#> NULL
#> non-numeric argument to binary operator
#> [] running 3 ...
#> <environment: 0x0000019b93a085c8> 
#> NULL
#> non-numeric argument to binary operator
#> [] running 4 ...
#> <environment: 0x0000019b93a085c8> 
#> NULL
#> non-numeric argument to binary operator
#> [] running 5 ...
#> <environment: 0x0000019b93a085c8> 
#> NULL
#> non-numeric argument to binary operator
#> [] running 6 ...
#> <environment: 0x0000019b93a085c8> 
#> NULL
#> non-numeric argument to binary operator
#> [] running 7 ...
#> <environment: 0x0000019b93a085c8> 
#> NULL
#> non-numeric argument to binary operator
#> [] running 8 ...
#> <environment: 0x0000019b93a085c8> 
#> NULL
#> non-numeric argument to binary operator
#> [] running 9 ...
#> <environment: 0x0000019b93a085c8> 
#> NULL
#> non-numeric argument to binary operator
#> [] running 10 ...
#> <environment: 0x0000019b93a085c8> 
#> NULL
#> non-numeric argument to binary operator

## 2. Visualization ------------------------------------------------------------
# growing season dividing
# Ipaper::write_fig({
#     par(cex = 1.1)
#     plot_season(INPUT, brks, ylab = "EVI", margin = 0.2, show.shade = FALSE)
# }, "Figure4_seasons.pdf", 9, 3.8)

^{Created on 2022-11-08 with reprex v2.0.2} 运行了非并行版本，前10个点还是显示null infile文件是老师给的百度网盘中的文件，没有调用程序开始生成的，但是和老师上述示例输出不一样

kongdd commented 1 year ago

https://github.com/eco-hydro/phenofit-scripts/blob/master/main_pkgs.R#L74

cymkG commented 1 year ago

also, in the above code.
nptperyear <- 23
nptperyear refers to the number of VI observations per year and per pixel. this variable was not quoted in the remaining code. Is it a redundant line? thank you very much for the answer

cymkG commented 1 year ago

刚好看到回复，谢谢老师, nptperyear has been quoted in phenofit::get_options(). problem has been solved

liuh1995 commented 1 year ago

孔老师，根据您提供的代码和原始数据（MEE2022_phenofit），生成新的的MOD13A2_Henan_2015_2020.rda，dates 有问题，导致后面无法有结果。您原始提供的MOD13A2_Henan_2015_2020.rda 中变量dates 本身是没有问题的。请问下是代码问题还是原始tif读取进来前需要做额外处理

将原始代码l$dates修改成将 dates= dates1 ，即运行正常了 dates1 = modis_date("2015-01-01", "2020-12-31", 16) l = read_rast(file_vi) %>% rast2mat() l_qc = read_rast(file_qc) %>% rast2mat() l_doy = read_rast(file_doy) %>% rast2mat()

I_grid = l$I_grid
d_coord = l$d_coord[I_grid, ]
data <- listk(
    VI = l$mat[I_grid, ] / 1e4,
    QC = l_qc$mat[I_grid, ],
    DOY = l_doy$mat[I_grid, ],
    dates = l$dates, qcFUN = qc_summary)       
save(data, d_coord, file = infile)

kongdd commented 1 year ago

@liuh1995 确实发现read_rast函数有一些问题，已经修复，https://github.com/eco-hydro/phenofit-scripts/commit/de67875c20088b2d155c946428210f3527d60d2e。

但是不会产生你说描述的错误，我无法复现。按照样例数据格式准备输入数据即可，没必要原封不动照搬代码。

eco-hydro / phenofit-scripts

Error of `phenofit-MOD13A2_HeNan.R` #6