Revision 7ff19746
Added by Benoit Parmentier over 9 years ago
| climate/research/oregon/interpolation/global_run_scalingup_assessment_part2.R | ||
|---|---|---|
| 5 | 5 |
#Analyses, figures, tables and data are also produced in the script. |
| 6 | 6 |
#AUTHOR: Benoit Parmentier |
| 7 | 7 |
#CREATED ON: 03/23/2014 |
| 8 |
#MODIFIED ON: 04/15/2015
|
|
| 8 |
#MODIFIED ON: 04/27/2015
|
|
| 9 | 9 |
#Version: 4 |
| 10 | 10 |
#PROJECT: Environmental Layers project |
| 11 |
#COMMENTS: analyses for run 10 global analyses,all regions 1500x4500km and other tiles
|
|
| 11 |
#COMMENTS: analyses for run 10 global analyses,all regions 1500x4500km with additional tiles to increase overlap
|
|
| 12 | 12 |
#TODO: |
| 13 | 13 |
#1) Split functions and master script |
| 14 | 14 |
#2) Make this is a script/function callable from the shell/bash |
| ... | ... | |
| 381 | 381 |
interpolation_method <- c("gam_CAI") #PARAM2
|
| 382 | 382 |
#out_suffix<-"run10_global_analyses_01282015" |
| 383 | 383 |
#out_suffix <- "output_run10_1000x3000_global_analyses_02102015" |
| 384 |
out_suffix <- "run10_1500x4500_global_analyses_pred_2003_04102015" #PARAM3
|
|
| 385 |
out_dir <- "/data/project/layers/commons/NEX_data/output_run10_1500x4500_global_analyses_pred_2003_04102015" #PARAM4
|
|
| 384 |
out_suffix <- "run10_1500x4500_global_analyses_04172015" #PARAM3
|
|
| 385 |
out_dir <- "/data/project/layers/commons/NEX_data/output_run10_1500x4500_global_analyses_04172015" #PARAM4
|
|
| 386 | 386 |
create_out_dir_param <- FALSE #PARAM 5 |
| 387 | 387 |
|
| 388 | 388 |
mosaic_plot <- FALSE #PARAM6 |
| 389 | 389 |
|
| 390 | 390 |
#if daily mosaics NULL then mosaicas all days of the year |
| 391 |
day_to_mosaic <- c("20030101","20030102","20030103","20030104","20030105",
|
|
| 392 |
"20030301","20030302","20030303","20030304","20030305", |
|
| 393 |
"20030501","20030502","20030503","20030504","20030505", |
|
| 394 |
"20030701","20030702","20030703","20030704","20030705", |
|
| 395 |
"20030901","20030902","20030903","20030904","20030905", |
|
| 396 |
"20031101","20031102","20031103","20031104","20031105") #PARAM7 |
|
| 391 |
|
|
| 392 |
day_to_mosaic <- c("20100101","20100102","20100103","20100104","20100105",
|
|
| 393 |
"20100301","20100302","20100303","20100304","20100305", |
|
| 394 |
"20100501","20100502","20100503","20100504","20100505", |
|
| 395 |
"20100701","20100702","20100703","20100704","20100705", |
|
| 396 |
"20100901","20100902","20100903","20100904","20100905", |
|
| 397 |
"20101101","20101102","20101103","20101104","20101105") #PARAM7 |
|
| 398 |
|
|
| 397 | 399 |
|
| 398 | 400 |
#CRS_locs_WGS84 <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +towgs84=0,0,0") #Station coords WGS84
|
| 399 | 401 |
CRS_WGS84 <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +towgs84=0,0,0") #Station coords WGS84 #CONSTANT1
|
| ... | ... | |
| 410 | 412 |
region_name <- "world" #PARAM 13 |
| 411 | 413 |
plot_region <- TRUE |
| 412 | 414 |
num_cores <- 10 #PARAM 14 |
| 415 |
reg_modified <- TRUE |
|
| 413 | 416 |
|
| 414 | 417 |
########################## START SCRIPT ############################## |
| 415 | 418 |
|
| ... | ... | |
| 435 | 438 |
tb_s <- read.table(file=file.path(out_dir,paste("tb_diagnostic_s_NA","_",out_suffix,".txt",sep="")),sep=",")
|
| 436 | 439 |
|
| 437 | 440 |
tb_month_s <- read.table(file=file.path(out_dir,paste("tb_month_diagnostic_s_NA","_",out_suffix,".txt",sep="")),sep=",")
|
| 438 |
#pred_data_month_info <- read.table(file=file.path(out_dir,paste("pred_data_month_info_",out_suffix,".txt",sep="")),sep=",")
|
|
| 441 |
pred_data_month_info <- read.table(file=file.path(out_dir,paste("pred_data_month_info_",out_suffix,".txt",sep="")),sep=",")
|
|
| 439 | 442 |
#pred_data_day_info <- read.table(file=file.path(out_dir,paste("pred_data_day_info_",out_suffix,".txt",sep="")),sep=",")
|
| 440 | 443 |
df_tile_processed <- read.table(file=file.path(out_dir,paste("df_tile_processed_",out_suffix,".txt",sep="")),sep=",")
|
| 441 | 444 |
|
| ... | ... | |
| 466 | 469 |
tb_all <- tb |
| 467 | 470 |
|
| 468 | 471 |
summary_metrics_v_all <- summary_metrics_v |
| 472 |
#deal with additional tiles... |
|
| 473 |
if(reg_modified==T){
|
|
| 474 |
|
|
| 475 |
summary_metrics_v_tmp <- summary_metrics_v |
|
| 476 |
summary_metrics_v_tmp$reg[summary_metrics_v_tmp$reg=="reg_1b"] <- "reg1" |
|
| 477 |
summary_metrics_v_tmp$reg[summary_metrics_v_tmp$reg=="reg_1c"] <- "reg1" |
|
| 478 |
summary_metrics_v_tmp$reg[summary_metrics_v_tmp$reg=="reg_3b"] <- "reg3" |
|
| 479 |
summary_metrics_v_tmp$reg_all <- summary_metrics_v$reg |
|
| 480 |
### |
|
| 481 |
summary_metrics_v<- summary_metrics_v_tmp |
|
| 482 |
|
|
| 483 |
### |
|
| 484 |
tb_tmp <- tb |
|
| 485 |
tb_tmp$reg[tb_tmp$reg=="reg_1b"] <- "reg1" |
|
| 486 |
tb_tmp$reg[tb_tmp$reg=="reg_1c"] <- "reg1" |
|
| 487 |
tb_tmp$reg[tb_tmp$reg=="reg_3b"] <- "reg3" |
|
| 488 |
### |
|
| 489 |
tb <- tb_tmp |
|
| 490 |
} |
|
| 491 |
|
|
| 492 |
table(summary_metrics_v_all$reg) |
|
| 493 |
table(summary_metrics_v$reg) |
|
| 494 |
table(tb_all$reg) |
|
| 495 |
table(tb$reg) |
|
| 469 | 496 |
|
| 470 | 497 |
############ PART 2: PRODUCE FIGURES ################ |
| 471 | 498 |
|
| ... | ... | |
| 480 | 507 |
list_shp_reg_files <- file.path("/data/project/layers/commons/NEX_data/",out_dir,
|
| 481 | 508 |
"shapefiles",basename(list_shp_reg_files)) |
| 482 | 509 |
|
| 510 |
#table(summary_metrics_v$reg) |
|
| 511 |
#table(summary_metrics_v$reg) |
|
| 483 | 512 |
|
| 484 | 513 |
### Potential function starts here: |
| 485 | 514 |
#function(in_dir,out_dir,list_shp_reg_files,title_str,region_name,num_cores,out_suffix,out_suffix) |
| ... | ... | |
| 871 | 900 |
|
| 872 | 901 |
test |
| 873 | 902 |
|
| 874 |
unique(test$tile_id) #72 tiles |
|
| 903 |
as.character(unique(test$tile_id)) #141 tiles |
|
| 904 |
|
|
| 875 | 905 |
dim(subset(test,test$predicted==365 & test$pred_mod=="mod1")) |
| 876 | 906 |
histogram(subset(test, test$pred_mod=="mod1")$predicted) |
| 877 | 907 |
unique(subset(test, test$pred_mod=="mod1")$predicted) |
| ... | ... | |
| 898 | 928 |
##### Figure 8: Breaking down accuracy by regions!! ##### |
| 899 | 929 |
|
| 900 | 930 |
#summary_metrics_v <- merge(summary_metrics_v,df_tile_processed,by="tile_id") |
| 901 |
table(summary_metrics_v$reg) |
|
| 902 | 931 |
|
| 903 | 932 |
## Figure 8a |
| 904 | 933 |
res_pix <- 480 |
| ... | ... | |
| 1052 | 1081 |
#undebug(plot_screen_raster_val) |
| 1053 | 1082 |
|
| 1054 | 1083 |
#world_m_list1<- plot_screen_raster_val(1,list_param_plot_screen_raster) |
| 1055 |
#world_m_list <- mclapply(1:10, list_param=list_param_plot_screen_raster, plot_screen_raster_val,mc.preschedule=FALSE,mc.cores = 5) #This is the end bracket from mclapply(...) statement
|
|
| 1056 |
world_m_list <- mclapply(1:length(l_dates), list_param=list_param_plot_screen_raster, plot_screen_raster_val,mc.preschedule=FALSE,mc.cores = 6) #This is the end bracket from mclapply(...) statement
|
|
| 1084 |
#world_m_list <- mclapply(11:30, list_param=list_param_plot_screen_raster, plot_screen_raster_val,mc.preschedule=FALSE,mc.cores = num_cores) #This is the end bracket from mclapply(...) statement
|
|
| 1085 |
world_m_list <- mclapply(1:length(l_dates), list_param=list_param_plot_screen_raster, plot_screen_raster_val,mc.preschedule=FALSE,mc.cores = num_cores) #This is the end bracket from mclapply(...) statement
|
|
| 1057 | 1086 |
|
| 1058 | 1087 |
#lf_world_mask_reg <- vector("list",length=length(lf_world_pred))
|
| 1059 | 1088 |
#for(i in 1:length(lf_world_pred)){
|
Also available in: Unified diff
global assessment part 2, changes to accomodate additional tiles in figure production