Revision ac6ed65d
Added by Benoit Parmentier about 9 years ago
| climate/research/oregon/interpolation/global_run_scalingup_mosaicing.R | ||
|---|---|---|
| 5 | 5 |
#Analyses, figures, tables and data are also produced in the script. |
| 6 | 6 |
#AUTHOR: Benoit Parmentier |
| 7 | 7 |
#CREATED ON: 04/14/2015 |
| 8 |
#MODIFIED ON: 10/05/2015
|
|
| 8 |
#MODIFIED ON: 10/09/2015
|
|
| 9 | 9 |
#Version: 5 |
| 10 | 10 |
#PROJECT: Environmental Layers project |
| 11 | 11 |
#COMMENTS: analyses run for reg4 1992 for test of mosaicing using 1500x4500km and other tiles |
| ... | ... | |
| 112 | 112 |
|
| 113 | 113 |
### Start new function here |
| 114 | 114 |
|
| 115 |
## Add numcores |
|
| 116 |
## use mclapply |
|
| 115 | 117 |
create_accuracy_metric_raster <- function(i, list_param){
|
| 116 | 118 |
#This function generates weights from a point location on a raster layer. |
| 117 | 119 |
#Note that the weights are normatlized on 0-1 scale using max and min values. |
| ... | ... | |
| 180 | 182 |
df_centroids <- merge(df_centroids,tb_date,by="tile_coord") |
| 181 | 183 |
|
| 182 | 184 |
#r1 <- raster(lf[i]) |
| 183 |
r1 <- raster(df_centroids$files[i]) |
|
| 184 |
r1[] <- df_centroids[[metric_name]][i] #improve this |
|
| 185 |
#set1f <- function(x){rep(NA, x)}
|
|
| 186 |
#r_init <- init(r_in, fun=set1f) |
|
| 187 |
|
|
| 188 |
raster_name_tmp <- paste("r_",metric_name,"_",interpolation_method,"_",date_processed,"_",out_suffix_str,sep="")
|
|
| 189 |
raster_name <- file.path(out_dir_str,raster_name_tmp) |
|
| 190 |
writeRaster(r1, NAflag=NA_flag_val,filename=raster_name,overwrite=TRUE) |
|
| 185 |
|
|
| 186 |
#function(j,df_centroids,metric_name,date_processed,interpolation,) |
|
| 187 |
#loop over files, make this a function later, works for now |
|
| 188 |
#use mclapply |
|
| 189 |
list_raster_name <- vector("list",length=length(lf))
|
|
| 190 |
|
|
| 191 |
for(j in 1:length(lf)){
|
|
| 192 |
|
|
| 193 |
inFilename <- df_centroids$files[j] |
|
| 194 |
r1 <- raster(inFilename) |
|
| 195 |
r1[] <- df_centroids[[metric_name]][j] #improve this |
|
| 196 |
#set1f <- function(x){rep(NA, x)}
|
|
| 197 |
#r_init <- init(r_in, fun=set1f) |
|
| 198 |
lf_tmp <- gsub(file_format,"",lf) |
|
| 199 |
|
|
| 200 |
extension_str <- extension(inFilename) |
|
| 201 |
raster_name_tmp <- gsub(extension_str,"",basename(inFilename)) |
|
| 202 |
outFilename <- file.path(out_dir,paste(raster_name_tmp,"_",metric_name,"_",out_suffix,file_format,sep="")) #for use in function later... |
|
| 203 |
|
|
| 204 |
#raster_name <- file.path(out_dir_str,paste(raster_name_tmp,"_",out_suffix,file_format,sep=""))#output file |
|
| 205 |
|
|
| 206 |
#raster_name_tmp <- paste("r_",metric_name,"_",interpolation_method,"_",date_processed,"_",out_suffix_str,sep="")
|
|
| 207 |
#raster_name <- file.path(out_dir_str,raster_name_tmp) |
|
| 208 |
writeRaster(r1, NAflag=NA_flag_val,filename=outFilename,overwrite=TRUE) |
|
| 209 |
list_raster_name[[j]] <- outFilename |
|
| 210 |
} |
|
| 211 |
|
|
| 191 | 212 |
|
| 192 |
raster_created_obj <- list(raster_name,df_centroids) |
|
| 193 |
names(raster_created_obj) <- c("raster_name","df_centroids")
|
|
| 213 |
raster_created_obj <- list(list_raster_name,df_centroids)
|
|
| 214 |
names(raster_created_obj) <- c("list_raster_name","df_centroids")
|
|
| 194 | 215 |
return(raster_created_obj) |
| 195 | 216 |
|
| 196 | 217 |
} |
| 197 | 218 |
|
| 198 | 219 |
#### end of function |
| 199 | 220 |
|
| 200 |
|
|
| 201 |
|
|
| 202 | 221 |
lf_mosaic1 <-list.files(path=file.path(in_dir_tiles), |
| 203 | 222 |
pattern=paste(".*.",day_to_mosaic[1],".*.tif$",sep=""),full.names=T) #choosing date 2...20100901
|
| 204 | 223 |
|
| ... | ... | |
| 223 | 242 |
names(list_param_accuracy_metric_raster) <- c("lf","tb","metric_name","pred_mod_name","y_var_name","interpolation_method",
|
| 224 | 243 |
"days_to_process","NA_flag_val","file_format","out_dir_str","out_suffix_str") |
| 225 | 244 |
debug(create_accuracy_metric_raster) |
| 226 |
test <- create_accuracy_metric_raster(1, list_param_accuracy_metric_raster) |
|
| 245 |
raster_created_obj <- create_accuracy_metric_raster(1, list_param_accuracy_metric_raster) |
|
| 246 |
|
|
| 247 |
lf_accuracy_raster <- unlist(raster_created_obj$list_raster_name) |
|
| 248 |
|
|
| 249 |
lf_r_rmse <- mclapply(1:length(days_to_process),FUN=create_accuracy_metric_raster,list_param=list_param_accuracy_metric_raster,mc.preschedule=FALSE,mc.cores = 3) |
|
| 227 | 250 |
|
| 228 | 251 |
#lf_mosaic <- lf_mosaic[1:20] |
| 229 | 252 |
r1 <- raster(lf_mosaic1[1]) |
| ... | ... | |
| 243 | 266 |
out_suffix_str <- paste(day_to_mosaic[i],out_suffix,sep="_") |
| 244 | 267 |
#undebug(mosaicFiles) |
| 245 | 268 |
#can also loop through methods!!! |
| 246 |
mosaic_edge_obj <- mosaicFiles(lf_mosaic1,mosaic_method="use_edge_weights", |
|
| 269 |
mosaic_edge_obj_prediction <- mosaicFiles(lf_mosaic1,mosaic_method="use_edge_weights",
|
|
| 247 | 270 |
num_cores=num_cores, |
| 248 | 271 |
python_bin=NULL, |
| 249 | 272 |
df_points=NULL,NA_flag=NA_flag_val, |
| 250 | 273 |
file_format=file_format,out_suffix=out_suffix_str, |
| 251 | 274 |
out_dir=out_dir) |
| 252 | 275 |
|
| 276 |
mosaic_method <- "use_edge_weights" #this is distance from edge |
|
| 277 |
out_suffix_str <- paste(day_to_mosaic[i],out_suffix,sep="_") |
|
| 278 |
#undebug(mosaicFiles) |
|
| 279 |
#can also loop through methods!!! |
|
| 280 |
mosaic_edge_obj_accuracy <- mosaicFiles(lf_accuracy_raster,mosaic_method="use_edge_weights", |
|
| 281 |
num_cores=num_cores, |
|
| 282 |
python_bin=NULL, |
|
| 283 |
df_points=NULL,NA_flag=NA_flag_val, |
|
| 284 |
file_format=file_format,out_suffix=out_suffix_str, |
|
| 285 |
out_dir=out_dir) |
|
| 253 | 286 |
|
| 254 | 287 |
#mosaic_unweighted_obj <- mosaicFiles(lf_mosaic1,mosaic_method="unweighted", |
| 255 | 288 |
# num_cores=num_cores, |
| ... | ... | |
| 259 | 292 |
# out_dir=out_dir) |
| 260 | 293 |
|
| 261 | 294 |
#list_mosaic_obj[[i]] <- list(unweighted=mosaic_unweighted_obj,edge=mosaic_edge_obj) |
| 262 |
list_mosaic_obj[[i]] <- list(unweighted=mosaic_unweighted_obj,edge=mosaic_edge_obj) |
|
| 295 |
#list_mosaic_obj[[i]] <- list(prediction=mosaic_edge_obj_accuracy,accuracy=mosaic_edge_obj_accuracy) |
|
| 296 |
|
|
| 297 |
list_mosaic_obj[[i]] <- list(prediction=mosaic_edge_obj_prediction,accuracy=mosaic_edge_obj_accuracy) |
|
| 263 | 298 |
} |
| 264 | 299 |
|
| 265 | 300 |
##################### |
| ... | ... | |
| 290 | 325 |
list_param= list_param_plot_mosaic, |
| 291 | 326 |
mc.preschedule=FALSE,mc.cores = num_cores) |
| 292 | 327 |
|
| 293 |
#################### |
|
| 294 |
#### Now difference figures... |
|
| 295 |
|
|
| 296 |
lf_obj1 <- list.files(path=out_dir,pattern="*unweighted.*.RData") |
|
| 297 |
lf_obj2 <- list.files(path=out_dir,pattern="*edge_.*.RData") |
|
| 298 |
|
|
| 299 |
lf1 <- unlist(lapply(lf_obj1,function(x){load_obj(x)[["mean_mosaic"]]}))
|
|
| 300 |
lf2 <- unlist(lapply(lf_obj2,function(x){load_obj(x)[["mean_mosaic"]]}))
|
|
| 301 |
|
|
| 302 |
out_suffix_str <- paste(paste(mosaicing_method,collapse="_"),day_to_mosaic,out_suffix,sep="_") |
|
| 303 |
|
|
| 304 |
list_param_plot_diff <- list(lf1=lf1,lf2=lf2,out_suffix=out_suffix_str) |
|
| 305 |
|
|
| 306 |
#debug(plot_diff_raster) |
|
| 307 |
#plot_diff_raster(1,list_param=list_param_plot_diff) |
|
| 308 |
|
|
| 309 | 328 |
num_cores <- 2 |
| 310 | 329 |
l_diff_png_files <- mclapply(1:length(lf1),FUN=plot_diff_raster,list_param= list_param_plot_diff, |
| 311 | 330 |
mc.preschedule=FALSE,mc.cores = num_cores) |
| 312 | 331 |
|
| 313 |
|
|
| 314 | 332 |
############### |
| 315 |
##### Get all the tiles togheter merged |
|
| 316 |
|
|
| 317 |
#ls -ltr ./reg*/*/*mean*.tif | wc |
|
| 318 |
in_dir1 <- "/data/project/layers/commons/NEX_data/mosaicing_data_test" |
|
| 319 |
lf_unweighted_20100831 <- list.files(path=in_dir1,pattern="r_m_mean_20100831.*.tif",recursive=T,full.names=T) |
|
| 320 |
lf_edge_weighted_20100831 <- list.files(path=in_dir1,pattern="r_.*.edge.*._mean_20100831.*.tif",recursive=T,full.names=T) |
|
| 321 |
lf_unweighted_20100901 <- list.files(path=in_dir1,pattern="r_m_mean_20100901.*.tif",recursive=T,full.names=T) |
|
| 322 |
lf_edge_weighted_20100901 <- list.files(path=in_dir1,pattern="r_.*.edge.*.mean_20100901.*.tif",recursive=T,full.names=T) |
|
| 323 |
|
|
| 324 |
output_fnames <- c("mean_unweighted_world_20100831_global_analyses_07012015.tif",
|
|
| 325 |
"mean_edge_weighted_world_20100831_global_analyses_07012015.tif", |
|
| 326 |
"mean_unweighted_world_20100901_global_analyses_07012015.tif", |
|
| 327 |
"mean_edge_weighted_world_20100901_global_analyses_07012015.tif" |
|
| 328 |
) |
|
| 329 |
|
|
| 330 |
list_lf <- list(lf_unweighted_20100831,lf_edge_weighted_20100831,lf_unweighted_20100901,lf_edge_weighted_20100901) |
|
| 331 |
out_dir_str <- "/data/project/layers/commons/NEX_data/mosaicing_data_test/mosaic_world_07012015" |
|
| 332 |
list_output_fnames <- vector("list",length=length(output_fnames))
|
|
| 333 |
for(i in 1:length(output_fnames)){
|
|
| 334 |
rast_ref <- file.path(out_dir_str,output_fnames[i]) #this is a the ref ouput file |
|
| 335 |
lf_to_mosaic <- list_lf[[i]] |
|
| 336 |
cmd_str <- paste("python","/usr/bin/gdal_merge.py","-o ",rast_ref,paste(lf_to_mosaic,collapse=" "))
|
|
| 337 |
system(cmd_str) |
|
| 338 |
list_output_fnames[[i]] <- rast_ref |
|
| 339 |
} |
|
| 340 |
|
|
| 341 |
#list_lf_m <- mixedsort(list.files(path=out_dir_str,pattern="mean.*.world.*.global_analyses_07012015.tif",full.names=T)) |
|
| 342 |
list_lf_m <- unlist(list_output_fnames) |
|
| 343 |
reg_name <- "world" |
|
| 344 |
out_suffix_str <- "mosaic_07092015" |
|
| 345 |
l_dates <- c("unweighted_20100831","edge_weighted_20100831","unweighted_20100901","edge_weighted_20100901")
|
|
| 346 |
list_param_plot_daily_mosaics <- list(list_lf_m,reg_name,out_dir_str,out_suffix_str,l_dates) |
|
| 347 |
names(list_param_plot_daily_mosaics) <- c("lf_m","reg_name","out_dir_str","out_suffix","l_dates")
|
|
| 348 |
|
|
| 349 |
#undebug(plot_daily_mosaics) |
|
| 350 |
#test<- plot_daily_mosaics(1,list_param_plot_daily_mosaics) |
|
| 351 |
num_cores <- 4 |
|
| 352 |
lf_plot <- mclapply(1:length(l_dates),FUN=plot_daily_mosaics,list_param=list_param_plot_daily_mosaics, |
|
| 353 |
mc.preschedule=FALSE,mc.cores = num_cores) |
|
| 354 |
|
|
| 355 |
r_test <- raster("r_m_use_edge_weighted_mean_19920101_reg4_mosaic_run10_1500x4500_global_analyses_10052015.tif")
|
|
| 356 |
plot(r_test) |
|
| 357 | 333 |
|
| 358 | 334 |
##################### END OF SCRIPT ###################### |
| 359 | 335 |
|
Also available in: Unified diff
testing accuracy metric function and mosaic