| 20 |
20 |
# 13) out_prefix: output suffix added to output names--it is the same in the interpolation script
|
| 21 |
21 |
#
|
| 22 |
22 |
#
|
| 23 |
|
# 14)
|
|
23 |
# 14) subampling:
|
| 24 |
24 |
# 15)..
|
|
25 |
#
|
|
26 |
#
|
|
27 |
|
| 25 |
28 |
#The output is a list of four shapefile names produced by the function:
|
|
29 |
|
| 26 |
30 |
#1) loc_stations: locations of stations as shapefile in EPSG 4326
|
| 27 |
31 |
#2) loc_stations_ghcn: ghcn daily data for the year range of interpolation (locally projected)
|
| 28 |
32 |
#3) daily_query_ghcn_data: ghcn daily data from daily query before application of quality flag
|
| ... | ... | |
| 85 |
89 |
qc_flags_stations <- list_param_prep$qc_flags_stations #flags allowed for the query from the GHCND??
|
| 86 |
90 |
out_prefix<-list_param_prep$out_prefix #"_365d_GAM_fus5_all_lstd_03012013" #User defined output prefix
|
| 87 |
91 |
|
| 88 |
|
## New parameters added for sub samplineg in areas with important density of meteo stations
|
|
92 |
## New parameters added for sub sampling in areas with important density of meteo stations
|
|
93 |
|
| 89 |
94 |
sub_sampling <- list_param$sub_sampling #if TRUE then monthly stations data are resampled
|
| 90 |
95 |
sub_sample_rnd <- list_param$sub_sample_rnd #if TRUE use random sampling in addition to spatial sub-sampling
|
|
96 |
|
| 91 |
97 |
target_range_nb <- list_param$target_range_nb # number of stations desired as min and max, convergence to min for now
|
|
98 |
#needs to be added in master script!!
|
|
99 |
target_range_daily_nb <- list_param$target_range_daily_nb #desired number range of daily stations
|
|
100 |
|
| 92 |
101 |
dist_range <- list_param$dist_range #distance range for pruning, usually (0,5) in km or 0,0.009*5 for degreee
|
| 93 |
102 |
step_dist <- list_param$step_dist #stepping distance used in pruning spatially, use 1km or 0.009 for degree data
|
| 94 |
103 |
|
| 95 |
104 |
## working directory is the same for input and output for this function
|
| 96 |
105 |
#setwd(in_path)
|
| 97 |
106 |
setwd(out_path)
|
|
107 |
|
| 98 |
108 |
##### STEP 1: Select station in the study area
|
| 99 |
109 |
|
| 100 |
110 |
filename<-sub(".shp","",fixed=TRUE,infile_reg_outline) #Removing the extension from file.
|
| ... | ... | |
| 212 |
222 |
data_RST_SDF$value<-data_RST_SDF$value/10 #TMax is the average max temp for monthy data
|
| 213 |
223 |
}
|
| 214 |
224 |
|
|
225 |
|
|
226 |
## Adding subsampling for daily stations...
|
|
227 |
|
|
228 |
#This must be set up in master script
|
|
229 |
#target_max_nb <- 100,000 #this is not actually used yet in the current implementation,can be set to very high value...
|
|
230 |
#target_min_nb <- 600 #this is the target number of stations we would like for daily and 1000x3000 tiles
|
|
231 |
#to be set by Alberto...
|
|
232 |
##max_dist <- 1000 # the maximum distance used for pruning ie removes stations that are closer than 1000m, this in degree...?
|
|
233 |
#max_dist <- 0.009*5 #5km in degree
|
|
234 |
#min_dist <- 0 #minimum distance to start with
|
|
235 |
#step_dist <- 0.009 #iteration step to remove the stations
|
|
236 |
|
|
237 |
#test5 <- sub_sampling_by_dist_nb_stat(target_range_nb=target_range_nb,dist_range=dist_range,step_dist=step_dist,data_in=data_month,sampling=T,combined=F)
|
|
238 |
if(sub_sampling_day==TRUE){
|
|
239 |
|
|
240 |
sub_sampling_obj <- sub_sampling_by_dist_nb_stat(target_range_nb=target_range_day_nb,dist_range=dist_range,step_dist=step_dist,data_in=data_RST_SDF,sampling=T,combined=F)
|
|
241 |
data_RST_SDF <- sub_sampling_obj$data #get sub-sampled data...for monhtly stations
|
|
242 |
#save the information for later use (validation at monthly step!!)
|
|
243 |
save(sub_sampling_obj,file= file.path(out_path,paste("sub_sampling_obj_","dayly_",interpolation_method,"_", out_prefix,".RData",sep="")))
|
|
244 |
}
|
|
245 |
|
|
246 |
#Make sure this is still a shapefile...!! This might need to be uncommented...
|
|
247 |
|
|
248 |
#coordinates(data_RST_SDF)<-cbind(data_RST_SDF$x,data_RST_SDF$y) #Transforming data_RST_SDF into a spatial point dataframe
|
|
249 |
#CRS_reg<-proj4string(data_RST_SDF)
|
|
250 |
#proj4string(data_RST_SDF)<-CRS_reg #Need to assign coordinates...
|
|
251 |
|
| 215 |
252 |
#write out a new shapefile (including .prj component)
|
| 216 |
253 |
outfile4<-file.path(out_path,paste("daily_covariates_ghcn_data_",var,"_",range_years[1],"_",range_years[2],out_prefix,".shp",sep="")) #Name of the file
|
| 217 |
254 |
writeOGR(data_RST_SDF,dsn= dirname(outfile4),layer= sub(".shp","",basename(outfile4)), driver="ESRI Shapefile",overwrite_layer=TRUE)
|
| ... | ... | |
| 318 |
355 |
|
| 319 |
356 |
#test5 <- sub_sampling_by_dist_nb_stat(target_range_nb=target_range_nb,dist_range=dist_range,step_dist=step_dist,data_in=data_month,sampling=T,combined=F)
|
| 320 |
357 |
if(sub_sampling==TRUE){
|
| 321 |
|
sub_sampling_obj <- sub_sampling_by_dist_nb_stat(target_range_nb=target_range_nb,dist_range=dist_range,step_dist=step_dist,data_in=data_month,sampling=T,combined=F)
|
|
358 |
sub_sampling_obj <- sub_sampling_by_dist_nb_stat(target_range_nb=target_range_nb,dist_range=dist_range,step_dist=step_dist,data_in=dst,sampling=T,combined=F)
|
| 322 |
359 |
dst <- sub_sampling_obj$data #get sub-sampled data...for monhtly stations
|
| 323 |
360 |
#save the information for later use (validation at monthly step!!)
|
| 324 |
361 |
save(sub_sampling_obj,file= file.path(out_path,paste("sub_sampling_obj_",interpolation_method,"_", out_prefix,".RData",sep="")))
|
| 325 |
362 |
}
|
| 326 |
363 |
|
|
364 |
#Make sure this is still a shapefile...!! This might need to be uncommented...
|
|
365 |
|
|
366 |
#coordinates(dst)<-cbind(dst$x,dst$y) #Transforming data_RST_SDF into a spatial point dataframe
|
|
367 |
#CRS_reg<-proj4string(data_reg)
|
|
368 |
#proj4string(dst)<-CRS_reg #Need to assign coordinates...
|
|
369 |
|
| 327 |
370 |
####
|
| 328 |
371 |
#write out a new shapefile (including .prj component)
|
| 329 |
372 |
dst$OID<-1:nrow(dst) #need a unique ID?
|
data preparation adding subsampling for daily stations