scripts/ChinaHW_cluster/Figure_S01_Annual_HWD.R

# 本代码用于检查clusterId年分布变化
source("scripts/main_pkgs.R")
load_all("/mnt/i/Research/phenology/latticeMap.R")
load_all()
InitCluster(6, kill = FALSE)

{
    sp_line <- path.mnt("N:/github/Research/Rcmip5/inst/extdata/shp/bou1_4l_south_sml.shp") %>%
        read_sf() %>%
        as_Spatial() %>%
        list("sp.lines", ., lwd = 0.4)
    date <- seq(as.Date("1961-01-01"), as.Date("2016-12-31"), by = "day")
    grid <- grid_d025
}

files = dir("OUTPUT/clusterId", "*.RDS", full.names = TRUE)#[1]
foreach(file = files, i = icount()) %dopar% {
    runningId(i)
    outfile <- gsub(".RDS$", ".pdf", file)
    if (file.exists(outfile)) return()
    clusterId = readRDS(file)
    # outfile <- glue("Figures/Figure_S01_Tavg_perc{prob*100}_Annual_HWD.pdf")
    
    mat <- apply_3d(!is.na(clusterId), by = year(date), FUN = rowSums2)
    df <- array_3dTo2d(mat) %>%
        as.data.table() %>%
        set_colnames(paste0("x", 1961:2016))
    data <- df %>%
        cbind(I = 1:nrow(.), .) %>%
        melt("I")
    {
        brks <- c(-Inf, 0, 1:5, 10, 20, 30, Inf) # days
        cols <- get_color(rcolors$amwg256, length(brks) - 1)
        cols[1] <- "white"
        p <- sp_plot(grid, data,
            formula = value ~ lon + lat | variable,
            brks = brks, colors = cols, key.num2factor = TRUE,
            colorkey = list(space = "bottom", labels = list(cex = 1.4)),
            aspect = 0.7,
            layout = c(7, 8),
            # sp.layout = sp_line,
            par.settings2 = list(axis.line = list(col = "black"))
        ) +
            layer_title(labels = 1961:2016, hjust = -0.1, vjust = 1.2) +
            theme_lattice(plot.margin = c(1, 1, 2, 1))
        write_fig(p, outfile, 15, 10, show = FALSE)
    }
    gc()
}
CUG-hydro/heatwave documentation built on Dec. 17, 2021, 1:53 p.m.