Revision 87cbaa98
Added by Benoit Parmentier over 11 years ago
| climate/research/oregon/interpolation/master_script_temp.R | ||
|---|---|---|
| 10 | 10 |
#STAGE 5: Output analyses-visualization of results for specific dates... |
| 11 | 11 |
# |
| 12 | 12 |
#AUTHOR: Benoit Parmentier |
| 13 |
#DATE: 04/02/2013
|
|
| 13 |
#DATE: 04/05/2013
|
|
| 14 | 14 |
|
| 15 | 15 |
#PROJECT: NCEAS INPLANT: Environment and Organisms --TASK#363, TASK$568-- |
| 16 | 16 |
|
| ... | ... | |
| 40 | 40 |
library(plotrix) |
| 41 | 41 |
|
| 42 | 42 |
### Parameters and arguments |
| 43 |
|
|
| 44 |
script_path<-"/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/" |
|
| 43 |
script_path<-"/home/parmentier/Data/IPLANT_project/env_layers_scripts/" |
|
| 44 |
#script_path<-"/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/"
|
|
| 45 | 45 |
#list_script_files<- |
| 46 | 46 |
#stages_to_run<-c(1,2,3,4,5) #May decide on antoher strategy later on... |
| 47 | 47 |
stages_to_run<-c(0,0,3,4,5) #May decide on antoher strategy later on... |
| ... | ... | |
| 54 | 54 |
|
| 55 | 55 |
#/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/climatology_03192013.py |
| 56 | 56 |
source(file.path(script_path,"covariates_production_temperatures_03212013.R")) |
| 57 |
source(file.path(script_path,"Database_stations_covariates_processing_function_03272013.R"))
|
|
| 58 |
source(file.path(script_path,"GAM_fusion_analysis_raster_prediction_multisampling_03272013.R"))
|
|
| 57 |
source(file.path(script_path,"Database_stations_covariates_processing_function_04042013.R"))
|
|
| 58 |
source(file.path(script_path,"GAM_fusion_analysis_raster_prediction_multisampling_04052013.R"))
|
|
| 59 | 59 |
source(file.path(script_path,"results_interpolation_date_output_analyses_03272013.R")) |
| 60 | 60 |
#source(file.path(script_path,"results_covariates_database_stations_output_analyses_04012013.R")) |
| 61 | 61 |
|
| ... | ... | |
| 75 | 75 |
############ STAGE 2: Covariate production ################ |
| 76 | 76 |
|
| 77 | 77 |
##Paths to inputs and output |
| 78 |
var<-"TMIN" |
|
| 79 |
in_path <- "/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/input_data/" |
|
| 80 |
out_path<- "/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/output_data/" |
|
| 78 |
var<-"TMAX" |
|
| 79 |
#in_path <- "/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/input_data/" |
|
| 80 |
#out_path<- "/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/output_data/" |
|
| 81 |
in_path <-"/home/parmentier/Data/IPLANT_project/Oregon_interpolation/Oregon_03142013/input" |
|
| 82 |
out_path<-"/home/parmentier/Data/IPLANT_project/Oregon_interpolation/Oregon_03142013/output" |
|
| 83 |
|
|
| 81 | 84 |
|
| 82 | 85 |
lc_path<-"/home/layers/data/land-cover/lc-consensus-global" |
| 83 |
infile_modis_grid<-"modis_sinusoidal_grid_world.shp" |
|
| 86 |
infile_modis_grid<-"modis_sinusoidal_grid_world.shp" #Give path!!! NEED TO CHANGE THIS... |
|
| 87 |
|
|
| 84 | 88 |
infile_elev<-"/home/layers/data/terrain/dem-cgiar-srtm-1km-tif/srtm_1km.tif" #this is the global file: replace later with the input produced by the DEM team |
| 85 |
infile_canheight<-"Simard_Pinto_3DGlobalVeg_JGR.tif" #Canopy height |
|
| 86 |
list_tiles_modis = c('h11v08','h11v07','h12v07','h12v08','h10v07','h10v08') #tile for Venezuel and surrounding area
|
|
| 87 |
infile_reg_outline="" #input region outline defined by polygon |
|
| 88 |
CRS_interp<-"+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs"; |
|
| 89 |
infile_canheight<-"/home/layers/data/land-cover/treeheight-simard2011/Simard_Pinto_3DGlobalVeg_JGR.tif" #Canopy height |
|
| 90 |
#list_tiles_modis <- c('h11v08','h11v07','h12v07','h12v08','h10v07','h10v08') #tile for Venezuel and surrounding area
|
|
| 91 |
list_tiles_modis <- c("h08v04","h09v04") #tiles for Oregon
|
|
| 92 |
|
|
| 93 |
#infile_reg_outline="" #input region outline defined by polygon: none for Venezuel |
|
| 94 |
infile_reg_outline <- "OR83M_state_outline.shp" #input region outline defined by polygon: Oregon |
|
| 95 |
|
|
| 96 |
#CRS_interp<-"+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs"; |
|
| 97 |
CRS_interp <-"+proj=lcc +lat_1=43 +lat_2=45.5 +lat_0=41.75 +lon_0=-120.5 +x_0=400000 +y_0=0 +ellps=GRS80 +units=m +no_defs"; |
|
| 98 |
|
|
| 89 | 99 |
CRS_locs_WGS84<-CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +towgs84=0,0,0") #Station coords WGS84
|
| 90 |
out_region_name<-"_venezuela_region" #generated on the fly |
|
| 91 |
out_suffix<-"_VE_03192013" |
|
| 92 |
ref_rast_name<-"" #local raster name defining resolution, exent, local projection--. set on the fly?? |
|
| 100 |
out_region_name<-"_oregon_region" #generated on the fly |
|
| 101 |
out_suffix<-"_OR_04042013" |
|
| 102 |
#ref_rast_name<-"" #local raster name defining resolution, exent, local projection--. set on the fly?? |
|
| 103 |
ref_rast_name<-"mean_day244_rescaled.rst" #local raster name defining resolution, exent: oregon |
|
| 104 |
|
|
| 93 | 105 |
#The names of covariates can be changed...these names should be output/input from covar script!!! |
| 94 |
rnames<-c("x","y","lon","lat","N","E","N_w","E_w","elev","slope","aspect","CANHEIGHT","DISTOC")
|
|
| 95 |
lc_names<-c("LC1","LC2","LC3","LC4","LC5","LC6","LC7","LC8","LC9","LC10","LC11","LC12")
|
|
| 106 |
rnames<-c("x","y","lon","lat","N","E","N_w","E_w","elev_s","slope","aspect","CANHEIGHT","DISTOC")
|
|
| 107 |
#lc_names<-c("LC1","LC2","LC3","LC4","LC5","LC6","LC7","LC8","LC9","LC10","LC11","LC12")
|
|
| 108 |
lc_names<-c("LC1","LC2","LC3","LC4","LC5","LC6","LC7","LC8","LC9","LC10") #use older version for continuity check to be changed
|
|
| 96 | 109 |
lst_names<-c("mm_01","mm_02","mm_03","mm_04","mm_05","mm_06","mm_07","mm_08","mm_09","mm_10","mm_11","mm_12",
|
| 97 | 110 |
"nobs_01","nobs_02","nobs_03","nobs_04","nobs_05","nobs_06","nobs_07","nobs_08", |
| 98 | 111 |
"nobs_09","nobs_10","nobs_11","nobs_12") |
| ... | ... | |
| 119 | 132 |
#var <- "TMIN" # name of the variables to keep: TMIN, TMAX or PRCP --already set up earlier |
| 120 | 133 |
|
| 121 | 134 |
#infile_covariates<-"covariates__venezuela_region__VE_01292013.tif" #this is an output from covariate script and used in stage 3 and stage 4 |
| 122 |
infile_covariates<-"covariates__venezuela_region_TMIN__VE_03192013.tif" #covariates stack for TMIN |
|
| 135 |
#infile_covariates<-"covariates__venezuela_region_TMIN__VE_03192013.tif" #covariates stack for TMIN |
|
| 136 |
infile_covariates<- "covariates_Oregon_region_TMAX__OR_04052013.tif" #Oregon covar TMAX from earlier codes...for continuity |
|
| 123 | 137 |
|
| 124 | 138 |
CRS_locs_WGS84<-CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +towgs84=0,0,0") #Station coords WGS84: same as earlier
|
| 125 |
infile1<- "outline_venezuela_region__VE_01292013.shp" #This is the shape file of outline of the study area #It is an input/output of the covariate script |
|
| 139 |
#infile1<- "outline_venezuela_region__VE_01292013.shp" #This is the shape file of outline of the study area #It is an input/output of the covariate script |
|
| 140 |
#infile_reg_outline <- "OR83M_state_outline.shp" #input region outline defined by polygon: Oregon |
|
| 141 |
infile1 <-"OR83M_state_outline.shp" #remove this parameter!!! |
|
| 142 |
|
|
| 126 | 143 |
#covar_names see stage 2 |
| 127 | 144 |
|
| 128 | 145 |
#specific to this stage |
| 129 | 146 |
db.name <- "ghcn" # name of the Postgres database |
| 130 |
range_years<-c("2000","2001") #right bound not included in the range!!
|
|
| 131 |
range_years_clim<-c("1981","2011") #right bound not included in the range!!
|
|
| 147 |
range_years<-c("2010","2011") #right bound not included in the range!!
|
|
| 148 |
range_years_clim<-c("2000","2011") #right bound not included in the range!!
|
|
| 132 | 149 |
infile2<-"/home/layers/data/climate/ghcn/v2.92-upd-2012052822/ghcnd-stations.txt" #This is the textfile of station locations from GHCND |
| 133 |
in_path <- "/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/input_data/" |
|
| 134 |
out_prefix<-"_365d_GAM_fus5_all_lstd_04022013" #User defined output prefix |
|
| 135 |
qc_flags_stations<-c("0","S") #flags allowed for screening after the query from the GHCND??
|
|
| 150 |
#in_path <- "/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/input_data/" |
|
| 151 |
out_prefix<-"_365d_GAM_fus5_all_lstd_04042013" #User defined output prefix |
|
| 152 |
#qc_flags_stations<-c("0","S") #flags allowed for screening after the query from the GHCND??
|
|
| 153 |
qc_flags_stations<-c("0") #flags allowed for screening after the query from the GHCND??
|
|
| 136 | 154 |
|
| 137 | 155 |
#list of 12 parameters for input in the function... |
| 138 | 156 |
|
| ... | ... | |
| 171 | 189 |
dates_selected<-"" # if empty string then predict for the full year specified earlier |
| 172 | 190 |
|
| 173 | 191 |
#Models to run...this can be change for each run |
| 174 |
list_models<-c("y_var ~ s(elev_1)",
|
|
| 192 |
list_models<-c("y_var ~ s(elev_s)",
|
|
| 175 | 193 |
"y_var ~ s(LST)", |
| 176 |
"y_var ~ s(elev_1,LST)",
|
|
| 177 |
"y_var ~ s(lat) + s(lon)+ s(elev_1)",
|
|
| 178 |
"y_var ~ s(lat,lon,elev_1)",
|
|
| 179 |
"y_var ~ s(lat,lon) + s(elev_1) + s(N_w,E_w) + s(LST)",
|
|
| 180 |
"y_var ~ s(lat,lon) + s(elev_1) + s(N_w,E_w) + s(LST) + s(LC2)",
|
|
| 181 |
"y_var ~ s(lat,lon) + s(elev_1) + s(N_w,E_w) + s(LST) + s(LC6)",
|
|
| 182 |
"y_var ~ s(lat,lon) + s(elev_1) + s(N_w,E_w) + s(LST) + s(DISTOC)")
|
|
| 194 |
"y_var ~ s(elev_s,LST)",
|
|
| 195 |
"y_var ~ s(lat) + s(lon)+ s(elev_s)",
|
|
| 196 |
"y_var ~ s(lat,lon,elev_s)",
|
|
| 197 |
"y_var ~ s(lat,lon) + s(elev_s) + s(N_w,E_w) + s(LST)",
|
|
| 198 |
"y_var ~ s(lat,lon) + s(elev_s) + s(N_w,E_w) + s(LST) + s(LC1)",
|
|
| 199 |
"y_var ~ s(lat,lon) + s(elev_s) + s(N_w,E_w) + s(LST) + s(LC2)",
|
|
| 200 |
"y_var ~ s(lat,lon) + s(elev_s) + s(N_w,E_w) + s(LST) + s(DISTOC)")
|
|
| 183 | 201 |
|
| 184 | 202 |
#Choose interpolation method... |
| 185 | 203 |
interpolation_method<-c("gam_fusion","gam_CAI") #other otpions to be added later
|
| 186 | 204 |
|
| 187 | 205 |
#Default name of LST avg to be matched |
| 188 | 206 |
lst_avg<-c("mm_01","mm_02","mm_03","mm_04","mm_05","mm_06","mm_07","mm_08","mm_09","mm_10","mm_11","mm_12")
|
| 189 |
in_path<-"/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/input_data" |
|
| 207 |
#in_path<-"/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/input_data"
|
|
| 190 | 208 |
#Create on the fly output folder... |
| 191 |
out_path<-"/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/output_data" |
|
| 192 |
script_path<-"/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/" |
|
| 209 |
#out_path<-"/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/output_data" |
|
| 210 |
#script_path<-"/home/parmentier/Data/IPLANT_project/Venezuela_interpolation/Venezuela_01142013/" |
|
| 211 |
|
|
| 193 | 212 |
|
| 194 | 213 |
#Collect all parameters in a list |
| 195 | 214 |
list_param_raster_prediction<-list(list_param_data_prep, |
Also available in: Unified diff
master script run first test for Oregon interpolation TMAX using covar stack