#MODIFIED ON: 06/09/2016            

    #### create a function to generate accuracy layers by tiles

  generate_ac_assessment_layers_by_tile <- function(lf,layers_option,df,df_tile_processed,metric_name,

                                                    var_pred,list_models,use_autokrige,pred_mod_name,

                                                    y_var_name,interpolation_method,

                                                    days_to_process,num_cores,NA_flag_val,file_format,

                                                    out_dir,out_suffix){ 

    #PARAMETERS:

    #lf

    #layers_option

    #df: can be tb,tb_s, data_v or data_s

    #df_tile_processed

    #metric_name <- "rmse" #RMSE, MAE etc.

    #var_pred: e.g. res_mod1, used for kriging of residuals in res_testing or res_training option

    #list_models: NULL then generate the formula for kriging

    #use_autokrige

    #pred_mod_name <- "mod1"

    #y_var_name 

    #interpolation_method #c("gam_CAI") #PARAM3

    #days_to_process <- day_to_mosaic

    #num_cores

    #NA_flag_val <- list_param$NA_flag_val

    #file_format <- list_param$file_format

    #out_dir

    #out_suffix

    #

    #OUTPUT

    #

    #add options to clean up file after use!!

    ###############################

    ### START #####

    out_dir_str <- out_dir

    out_suffix_str <- out_suffix

    #lf <- lf_mosaic

    #Improved by adding multicores option

    num_cores_tmp <- num_cores

    if(layers_option=="ac_training" | layers_option=="ac_testing"){

      list_param_accuracy_metric_raster <- list(lf,df,metric_name,pred_mod_name,y_var_name,interpolation_method,

                                              days_to_process,num_cores_tmp,NA_flag_val,file_format,out_dir_str,out_suffix_str) 

      names(list_param_accuracy_metric_raster) <- c("lf","tb","metric_name","pred_mod_name","y_var_name","interpolation_method",

                                                  "days_to_process","num_cores","NA_flag_val","file_format","out_dir_str","out_suffix_str") 

      list_raster_created_obj <- lapply(1:length(days_to_process),FUN=create_accuracy_metric_raster,

                                      list_param=list_param_accuracy_metric_raster)

      #debug(create_accuracy_metric_raster)

      #list_raster_created_obj <- lapply(1:1,FUN=create_accuracy_metric_raster,

      #                                  list_param=list_param_accuracy_metric_raster)

      #raster_created_obj <- create_accuracy_metric_raster(1, list_param_accuracy_metric_raster)

      #Extract list of files for rmse and date 1 (19920101), there should be 28 raster images

      lf_accuracy_raster <- lapply(1:length(list_raster_created_obj),FUN=function(i){unlist(list_raster_created_obj[[i]]$list_raster_name)}) 

      lf_ac_assessment_layers <- lf_accuracy_raster

    }

    if(layers_option=="res_training" | layers_option=="res_testing"){

      ## Create accuracy surface by kriging

      #num_cores_tmp <-num_cores

      #lf_day_tiles  <- lf_mosaic #list of raster files by dates

      lf_day_tiles  <- lf #list of raster files by dates

      #data_df <- data_day_v # data.frame table/spdf containing stations with residuals and variable

      #df_tile_processed  #tiles processed during assessment usually by region

      #var_pred  #variable being modeled

      #if not list of models is provided generate one

      if(is.null(list_models)){

        list_models <- paste(var_pred,"~","1",sep=" ") #can krige any variable here

      }

      #use_autokrige #if TRUE use automap/gstat package

      #y_var_name  #"dailyTmax" #PARAM2

      #interpolation_method #c("gam_CAI") #PARAM3, need to select reg!!

      #date_processed #can be a monthly layer

      #num_cores #number of cores used

      #NA_flag_val 

      #file_format 

      #out_dir_str <- out_dir #change to specific separate dir??

      #out_suffix_str <- out_suffix 

      #days_to_process <- day_to_mosaic

      #out_suffix_str <- paste("data_day_v_",out_suffix,sep="") 

      out_suffix_str <- paste(layers_option,var_pred,out_suffix,sep="")

      #browser()

      df_tile_processed$path_NEX <- as.character(df_tile_processed$path_NEX) 

      df_tile_processed$reg <- basename(dirname(df_tile_processed$path_NEX))

      ##By regions, selected earlier

      #for(k in 1:length(region_names)){

      df_tile_processed_reg <- subset(df_tile_processed,reg==region_selected)#use reg4

      #i<-1 #loop by days/date to process!!

      #test on the first day 

      list_param_create_accuracy_residuals_raster <- list(lf,df,df_tile_processed_reg,

                                                        var_pred,list_models,use_autokrige,y_var_name,interpolation_method,

                                                        days_to_process,num_cores_tmp,NA_flag_val,file_format,out_dir_str,

                                                        out_suffix_str) 

      names(list_param_create_accuracy_residuals_raster) <- c("lf_day_tiles","data_df","df_tile_processed_reg",

                                                            "var_pred","list_models","use_autokrige","y_var_name","interpolation_method",

                                                            "days_to_process","num_cores_tmp","NA_flag_val","file_format","out_dir_str",

                                                            "out_suffix_str") 

      #browser()  

      list_create_accuracy_residuals_raster_obj <- lapply(1:length(day_to_mosaic),FUN=create_accuracy_residuals_raster,

                                                        list_param=list_param_create_accuracy_residuals_raster)

      #undebug(create_accuracy_residuals_raster)

      #list_create_accuracy_residuals_raster_obj <- lapply(1:1,FUN=create_accuracy_residuals_raster,

      #                                list_param=list_param_create_accuracy_residuals_raster)

      #create_accuracy_residuals_raster_obj <- create_accuracy_residuals_raster(1, list_param_create_accuracy_residuals_raster_obj)

      #note that three tiles did not produce a residuals surface!!! find out more later, join the output

      #to df_raste_tile to keep track of which one did not work...

      #lf_accuracy_residuals_raster <- as.character(unlist(lapply(1:length(list_create_accuracy_residuals_raster_obj),FUN=function(i,x){unlist(extract_from_list_obj(x[[i]]$list_pred_res_obj,"raster_name"))},x=list_create_accuracy_residuals_raster_obj))) 

      lf_accuracy_residuals_raster <- lapply(1:length(list_create_accuracy_residuals_raster_obj),FUN=function(i,x){as.character(unlist(extract_from_list_obj(x[[i]]$list_pred_res_obj,"raster_name")))},x=list_create_accuracy_residuals_raster_obj)

      lf_ac_assessment_layers <- lf_accuracy_residuals_raster

    }

    return(lf_ac_assessment_layers)

  }

  browser()

  #try running the function now:

  lf <- lf_mosaic[i]

  df <- tb #for ac_testing

  days_to_process <- day_to_mosaic[i]

  #browser()

  debug(generate_ac_assessment_layers_by_tile)

  generate_ac_assessment_layers_by_tile(lf,layers_option,df,df_tile_processed,metric_name,

                                                    var_pred,list_models,use_autokrige,pred_mod_name,

                                                    y_var_name,interpolation_method,

                                                    days_to_process,num_cores,NA_flag_val,file_format,

                                                    out_dir,out_suffix)   #### create a function to generate accuracy layers by tiles