Revision b4a28e3d
Added by Benoit Parmentier about 10 years ago
| climate/research/oregon/interpolation/global_run_scalingup_assessment_part1.R | ||
|---|---|---|
| 5 | 5 |
#Part 1 create summary tables and inputs for figure in part 2 and part 3. |
| 6 | 6 |
#AUTHOR: Benoit Parmentier |
| 7 | 7 |
#CREATED ON: 03/23/2014 |
| 8 |
#MODIFIED ON: 06/27/2014
|
|
| 8 |
#MODIFIED ON: 08/14/2014
|
|
| 9 | 9 |
#Version: 3 |
| 10 | 10 |
#PROJECT: Environmental Layers project |
| 11 | 11 |
#TO DO: |
| ... | ... | |
| 15 | 15 |
# - generate data_s and data_v inputs as giant table |
| 16 | 16 |
# - generate accuracy for mosaic (part 2 and part3) |
| 17 | 17 |
# - clean up |
| 18 |
|
|
| 19 |
#First source file: |
|
| 20 |
#source /nobackupp4/aguzman4/climateLayers/sharedModules/etc/environ.sh |
|
| 18 | 21 |
################################################################################################# |
| 19 | 22 |
|
| 20 | 23 |
### Loading R library and packages |
| ... | ... | |
| 46 | 49 |
|
| 47 | 50 |
#### FUNCTION USED IN SCRIPT |
| 48 | 51 |
|
| 49 |
function_analyses_paper1 <-"contribution_of_covariates_paper_interpolation_functions_05212014.R" #first interp paper
|
|
| 50 |
function_analyses_paper2 <-"multi_timescales_paper_interpolation_functions_05052014.R"
|
|
| 52 |
function_analyses_paper1 <-"contribution_of_covariates_paper_interpolation_functions_07182014.R" #first interp paper
|
|
| 53 |
function_analyses_paper2 <-"multi_timescales_paper_interpolation_functions_08132014.R"
|
|
| 51 | 54 |
|
| 52 | 55 |
load_obj <- function(f) |
| 53 | 56 |
{
|
| ... | ... | |
| 310 | 313 |
#### Parameters and constants |
| 311 | 314 |
|
| 312 | 315 |
#in_dir1 <- "/data/project/layers/commons/NEX_data/test_run1_03232014/output" #On Atlas |
| 313 |
in_dir1 <- "/nobackupp4/aguzman4/climateLayers/output10Deg/reg1/" #On NEX |
|
| 316 |
#in_dir1 <- "/nobackupp4/aguzman4/climateLayers/output10Deg/reg1/" #On NEX |
|
| 317 |
in_dir1 <- "/nobackupp4/aguzman4/climateLayers/output20Deg/" |
|
| 318 |
#/nobackupp4/aguzman4/climateLayers/output20Deg/reg5/20.0_0.0/ |
|
| 319 |
#/nobackupp4/aguzman4/climateLayers/output20Deg/reg3/-20.0_-70.0/ |
|
| 320 |
#/nobackupp4/aguzman4/climateLayers/output20Deg/reg5/20.0_30.0/ |
|
| 321 |
#/nobackupp4/aguzman4/climateLayers/output20Deg/reg4/40.0_0.0/ |
|
| 322 |
#/nobackupp4/aguzman4/climateLayers/output20Deg/reg5/20.0_-10.0/ |
|
| 323 |
#/nobackupp4/aguzman4/climateLayers/output20Deg/reg4/50.0_0.0/ |
|
| 324 |
#/nobackupp4/aguzman4/climateLayers/output20Deg/reg6/60.0_40.0/ |
|
| 325 |
#/nobackupp4/aguzman4/climateLayers/output20Deg/reg6/30.0_40.0/ |
|
| 326 |
|
|
| 327 |
|
|
| 314 | 328 |
#/nobackupp4/aguzman4/climateLayers/output10Deg/reg1/finished.txt |
| 315 |
in_dir_list <- list.dirs(path=in_dir1,recursive=FALSE) #get the list of directories with resutls by 10x10 degree tiles |
|
| 329 |
#in_dir_list <- list.dirs(path=in_dir1,recursive=FALSE) #get the list of directories with resutls by 10x10 degree tiles |
|
| 330 |
in_dir_list <- c("/nobackupp4/aguzman4/climateLayers/output20Deg/reg5/20.0_0.0/",
|
|
| 331 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg3/-20.0_-70.0/", |
|
| 332 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg5/20.0_30.0/", |
|
| 333 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg4/40.0_0.0/", |
|
| 334 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg5/20.0_-10.0/", |
|
| 335 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg4/50.0_0.0/", |
|
| 336 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg6/60.0_40.0/", |
|
| 337 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg6/30.0_40.0/") |
|
| 338 |
|
|
| 316 | 339 |
#use subset for now: |
| 317 | 340 |
|
| 318 | 341 |
#in_dir_list <- c( |
| ... | ... | |
| 322 | 345 |
#in_dir_list <- file.path(in_dir1,read.table(file.path(in_dir1,"processed.txt"))$V1) |
| 323 | 346 |
#in_dir_list <- as.list(in_dir_list[-1]) |
| 324 | 347 |
#in_dir_list <- in_dir_list[grep("bak",basename(basename(in_dir_list)),invert=TRUE)] #the first one is the in_dir1
|
| 325 |
#in_dir_shp <- in_dir_list[grep("shapefiles",basename(in_dir_list),invert=FALSE)] #select directory with shapefiles...
|
|
| 326 |
in_dir_shp <- "/nobackupp4/aguzman4/climateLayers/output10Deg/reg1/subset/shapefiles/" |
|
| 327 |
in_dir_shp_list <- list.files(in_dir_shp,".shp") |
|
| 348 |
in_dir_shp <- in_dir_list[grep("shapefiles",basename(in_dir_list),invert=FALSE)] #select directory with shapefiles...
|
|
| 349 |
in_dir_shp <- c( |
|
| 350 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg3/subset/shapefiles/", |
|
| 351 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg5/subset/shapefiles/", |
|
| 352 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg4/subset/shapefiles/", |
|
| 353 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg5/subset/shapefiles/", |
|
| 354 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg4/subset/shapefiles/", |
|
| 355 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg6/subset/shapefiles/", |
|
| 356 |
"/nobackupp4/aguzman4/climateLayers/output20Deg/reg6/subset/shapefiles/") |
|
| 357 |
#in_dir_shp <- "/nobackupp4/aguzman4/climateLayers/output10Deg/reg1/subset/shapefiles/" |
|
| 358 |
#in_dir_shp <- "/nobackupp4/aguzman4/climateLayers/output20Deg/reg2/subset/shapefiles" |
|
| 359 |
in_dir_shp_list <- list.files(in_dir_shp,".shp",full.names=T) |
|
| 328 | 360 |
|
| 329 | 361 |
#in_dir_list <- in_dir_list[grep("shapefiles",basename(in_dir_list),invert=TRUE)]
|
| 330 | 362 |
#the first one is the in_dir1 |
| 331 | 363 |
# the last directory contains shapefiles |
| 332 | 364 |
y_var_name <- "dailyTmax" |
| 333 | 365 |
interpolation_method <- c("gam_CAI")
|
| 334 |
out_prefix<-"run3_global_analyses_06192014"
|
|
| 366 |
out_prefix<-"run4_global_analyses_08142014"
|
|
| 335 | 367 |
|
| 336 | 368 |
#out_dir<-"/data/project/layers/commons/NEX_data/" #On NCEAS Atlas |
| 337 | 369 |
out_dir <- "/nobackup/bparmen1/" #on NEX |
| ... | ... | |
| 371 | 403 |
|
| 372 | 404 |
lf_covar_obj <- lapply(in_dir_list,FUN=function(x){list.files(path=x,pattern="covar_obj.*.RData",full.names=T)})
|
| 373 | 405 |
lf_covar_tif <- lapply(in_dir_list,FUN=function(x){list.files(path=x,pattern="covar.*.tif",full.names=T)})
|
| 406 |
lf_diagnostic_obj <- lapply(in_dir_list,FUN=function(x){list.files(path=x,pattern="diagnostic.*.RData",full.names=T)})
|
|
| 374 | 407 |
|
| 375 | 408 |
########################## START SCRIPT ############################## |
| 376 | 409 |
|
| ... | ... | |
| 483 | 516 |
###################################################### |
| 484 | 517 |
####### PART 4: Get shapefile tiling with centroids ### |
| 485 | 518 |
|
| 486 |
#in_dir_shp <- "/nobackupp4/aguzman4/climateLayers/output4/subset/shapefiles/" |
|
| 487 |
|
|
| 488 | 519 |
#get shape files for the region being assessed: |
| 489 |
list_shp_global_tiles_files <- list.files(path=in_dir_shp,pattern="*.shp") |
|
| 490 |
l_shp<-lapply(1:length(in_dir_shp_list),FUN=function(i){paste(strsplit(in_dir_shp_list[i],"_")[[1]][2:3],collapse="_")})
|
|
| 491 |
list_shp_global_tiles_files <- l_shp |
|
| 492 |
pattern_str <- basename(in_dir_list) |
|
| 493 |
list_shp_global_tiles_files |
|
| 494 |
list_shp_reg_files <- lapply(pattern_str,function(x){list_shp_global_tiles_files[grep(x,invert=FALSE,list_shp_global_tiles_files)]}) #select directory with shapefiles...
|
|
| 495 |
df_tile_processed$shp_files <- unlist(list_shp_reg_files) |
|
| 520 |
|
|
| 521 |
list_shp_world <- list.files(path=in_dir_shp,pattern=".*.shp",full.names=T) |
|
| 522 |
l_shp <- unlist(lapply(1:length(list_shp_world),FUN=function(i){paste(strsplit(list_shp_world[i],"_")[[1]][2:3],collapse="_")}))
|
|
| 523 |
matching_index <- match(basename(in_dir_list),l_shp) |
|
| 524 |
list_shp_reg_files <- list_shp_world[matching_index] |
|
| 525 |
df_tile_processed$shp_files <-list_shp_world[matching_index] |
|
| 496 | 526 |
|
| 497 | 527 |
tx<-strsplit(as.character(df_tile_processed$tile_coord),"_") |
| 498 | 528 |
lat<- as.numeric(lapply(1:length(tx),function(i,x){x[[i]][1]},x=tx))
|
| ... | ... | |
| 500 | 530 |
df_tile_processed$lat <- lat |
| 501 | 531 |
df_tile_processed$lon <- long |
| 502 | 532 |
|
| 503 |
list_shp_world <- list.files(in_dir_shp,".shp") |
|
| 504 |
l_shp <- unlist(lapply(1:length(list_shp_world),FUN=function(i){paste(strsplit(list_shp_world[i],"_")[[1]][2:3],collapse="_")}))
|
|
| 505 |
list_shp_reg_files <- as.character(df_tile_processed$tile_coord) |
|
| 506 |
#matching_index <- match(l_shp,list_shp_reg_files) |
|
| 507 |
matching_index <- match(list_shp_reg_files,l_shp) |
|
| 508 |
df_tile_processed$shp_files <-list_shp_world[matching_index] |
|
| 509 |
|
|
| 510 |
df_tile_processed <- na.omit(df_tile_processed) #remove other list of folders irrelevant |
|
| 511 |
list_shp_reg_files <- df_tile_processed$shp_files |
|
| 512 |
list_shp_reg_files<- file.path(in_dir_shp,list_shp_reg_files) |
|
| 513 |
|
|
| 514 | 533 |
#put that list in the df_processed and also the centroids!! |
| 515 | 534 |
write.table(df_tile_processed, |
| 516 | 535 |
file=file.path(out_dir,paste("df_tile_processed_",out_prefix,".txt",sep="")),sep=",")
|
| ... | ... | |
| 774 | 793 |
|
| 775 | 794 |
### This assumes the tree structure has been replicated on Atlas: |
| 776 | 795 |
#for i in 1:length(df_tiled_processed$tile_coord) |
| 777 |
output_atlas_dir <- "/data/project/layers/commons/NEX_data/output_run3_global_analyses_06192014/output10Deg/reg1" |
|
| 796 |
#output_atlas_dir <- "/data/project/layers/commons/NEX_data/output_run3_global_analyses_06192014/output10Deg/reg1" |
|
| 797 |
output_atlas_dir <- "/data/project/layers/commons/NEX_data/output_run4_global_analyses_08142014/output20Deg" |
|
| 798 |
|
|
| 799 |
#Make directories on ATLAS |
|
| 778 | 800 |
#for (i in 1:length(df_tiled_processed$tile_coord)){
|
| 779 | 801 |
# create_dir_fun(file.path(output_atlas_dir,as.character(df_tiled_processed$tile_coord[i])),out_suffix=NULL) |
| 780 | 802 |
#} |
| 781 | 803 |
|
| 804 |
#Make directories on ATLAS for shapefiles |
|
| 805 |
#for (i in 1:length(df_tiled_processed$tile_coord)){
|
|
| 806 |
# create_dir_fun(file.path(output_atlas_dir,as.character(df_tiled_processed$tile_coord[i]),"/shapefiles"),out_suffix=NULL) |
|
| 807 |
#} |
|
| 808 |
|
|
| 809 |
|
|
| 810 |
#Copy summary textfiles and mosaic back to atlas |
|
| 811 |
|
|
| 812 |
Atlas_dir <- file.path("/data/project/layers/commons/NEX_data/",basename(out_dir))#,"output/subset/shapefiles")
|
|
| 813 |
Atlas_hostname <- "parmentier@atlas.nceas.ucsb.edu" |
|
| 814 |
lf_cp_f <- list.files(out_dir,full.names=T)#copy all files can filter later |
|
| 815 |
filenames_NEX <- paste(lf_cp_f,collapse=" ") #copy raster prediction object |
|
| 816 |
cmd_str <- paste("scp -p",filenames_NEX,paste(Atlas_hostname,Atlas_dir,sep=":"), sep=" ")
|
|
| 817 |
system(cmd_str) |
|
| 818 |
|
|
| 819 |
#system("scp -p ./*.txt parmentier@atlas.nceas.ucsb.edu:/data/project/layers/commons/NEX_data/output_run2_global_analyses_05122014")
|
|
| 820 |
#system("scp -p ./*.txt ./*.tif parmentier@atlas.nceas.ucsb.edu:/data/project/layers/commons/NEX_data/output_run2_global_analyses_05122014")
|
|
| 821 |
|
|
| 822 |
#### COPY SHAPEFILES, TIF MOSAIC, COMBINED TEXT FILES etc... |
|
| 823 |
|
|
| 824 |
#copy shapefiles defining regions |
|
| 825 |
Atlas_dir <- file.path("/data/project/layers/commons/NEX_data/",basename(out_dir),"output/subset/shapefiles")
|
|
| 826 |
Atlas_hostname <- "parmentier@atlas.nceas.ucsb.edu" |
|
| 827 |
lf_cp_shp <- df_tile_processed$shp_files #get all the files... |
|
| 828 |
|
|
| 829 |
#lf_cp_shp <- list.files(in_dir_shp, ".shp",full.names=T) |
|
| 830 |
list_tile_scp <- 1:8 |
|
| 831 |
|
|
| 832 |
for (j in 1:length(list_tile_scp)){
|
|
| 833 |
tile_nb <- list_tile_scp[j] |
|
| 834 |
|
|
| 835 |
in_dir_tile <-dirname(df_tile_processed$shp_files[tile_nb]) |
|
| 836 |
#/data/project/layers/commons/NEX_data/output_run2_05122014/output |
|
| 837 |
#output_atlas_dir |
|
| 838 |
#Atlas_dir <- file.path(file.path("/data/project/layers/commons/NEX_data/",basename(out_dir),"output"),in_dir_tile)
|
|
| 839 |
Atlas_dir <- file.path(output_atlas_dir,as.character(df_tile_processed$tile_coord[j]),"/shapefiles") |
|
| 840 |
|
|
| 841 |
Atlas_hostname <- "parmentier@atlas.nceas.ucsb.edu" |
|
| 842 |
|
|
| 843 |
filenames_NEX <- paste(lf_cp_shp,collapse=" ") #copy raster prediction object |
|
| 844 |
cmd_str <- paste("scp -p",filenames_NEX,paste(Atlas_hostname,Atlas_dir,sep=":"), sep=" ")
|
|
| 845 |
system(cmd_str) |
|
| 846 |
} |
|
| 847 |
|
|
| 782 | 848 |
#### FIRST COPY DATA FOR SPECIFIC TILES ##### |
| 783 | 849 |
#Copy specific tiles info back...This assumes that the tree structre |
| 784 | 850 |
#has been created on ATLAS: |
| 785 | 851 |
#../$out_dir/ouput/tile_coord |
| 786 | 852 |
|
| 787 |
list_tile_scp <- c(1,2) |
|
| 853 |
#list_tile_scp <- c(1,2) |
|
| 854 |
list_tile_scp <- 1:8 |
|
| 788 | 855 |
|
| 789 | 856 |
for (j in 1:length(list_tile_scp)){
|
| 790 | 857 |
tile_nb <- list_tile_scp[j] |
| ... | ... | |
| 827 | 894 |
lf_cp_month[[i]] <- unlist(lapply(1:nb_mod,FUN=function(x){lf_clim_tif[[x]][[index]][[tile_nb]]}))
|
| 828 | 895 |
} |
| 829 | 896 |
##Add RData object for specified tile... |
| 830 |
lf_cp_RData_tif <- c(lf_covar_obj[tile_nb],lf_covar_tif[tile_nb],list_raster_obj_files[[tile_nb]]) |
|
| 897 |
lf_cp_RData_tif <- c(lf_covar_obj[tile_nb],lf_covar_tif[tile_nb],list_raster_obj_files[[tile_nb]],lf_diagnostic_obj[[tile_nb]])
|
|
| 831 | 898 |
#unlist(lf_cp_RData_tif) |
| 832 | 899 |
lf_cp <- unlist(c(lf_cp_day,lf_cp_month,lf_cp_RData_tif)) |
| 833 | 900 |
#lf_cp <- c(unlist(c(lf_cp_day,lf_cp_month)),list_raster_obj_files[tile_nb]) |
| ... | ... | |
| 837 | 904 |
system(cmd_str) |
| 838 | 905 |
} |
| 839 | 906 |
|
| 840 |
#### COPY SHAPEFILES, TIF MOSAIC, COMBINED TEXT FILES etc... |
|
| 841 |
|
|
| 842 |
#copy shapefiles defining regions |
|
| 843 |
Atlas_dir <- file.path("/data/project/layers/commons/NEX_data/",basename(out_dir),"output/subset/shapefiles")
|
|
| 844 |
Atlas_hostname <- "parmentier@atlas.nceas.ucsb.edu" |
|
| 845 |
lf_cp_shp <- list.files(in_dir_shp,full.names=T) #get all the files... |
|
| 846 |
#lf_cp_shp <- list.files(in_dir_shp, ".shp",full.names=T) |
|
| 847 |
|
|
| 848 |
filenames_NEX <- paste(lf_cp_shp,collapse=" ") #copy raster prediction object |
|
| 849 |
cmd_str <- paste("scp -p",filenames_NEX,paste(Atlas_hostname,Atlas_dir,sep=":"), sep=" ")
|
|
| 850 |
system(cmd_str) |
|
| 851 |
|
|
| 852 |
#Copy summary textfiles and mosaic back to atlas |
|
| 853 |
|
|
| 854 |
Atlas_dir <- file.path("/data/project/layers/commons/NEX_data/",basename(out_dir))#,"output/subset/shapefiles")
|
|
| 855 |
Atlas_hostname <- "parmentier@atlas.nceas.ucsb.edu" |
|
| 856 |
lf_cp_f <- list.files(out_dir,full.names=T)#copy all files can filter later |
|
| 857 |
filenames_NEX <- paste(lf_cp_f,collapse=" ") #copy raster prediction object |
|
| 858 |
cmd_str <- paste("scp -p",filenames_NEX,paste(Atlas_hostname,Atlas_dir,sep=":"), sep=" ")
|
|
| 859 |
system(cmd_str) |
|
| 860 |
|
|
| 861 |
#system("scp -p ./*.txt parmentier@atlas.nceas.ucsb.edu:/data/project/layers/commons/NEX_data/output_run2_global_analyses_05122014")
|
|
| 862 |
#system("scp -p ./*.txt ./*.tif parmentier@atlas.nceas.ucsb.edu:/data/project/layers/commons/NEX_data/output_run2_global_analyses_05122014")
|
|
| 863 |
|
|
| 864 |
system("scp -p /nobackupp4/aguzman4/climateLayers/output4/subset/shapefiles/* parmentier@atlas.nceas.ucsb.edu:/data/project/layers/commons/NEX_data/shapefiles")
|
|
| 865 | 907 |
|
| 866 | 908 |
##################### END OF SCRIPT ###################### |
Also available in: Unified diff
run 4 assessment part1, gam fitting comparison k for 6 specific tiles