0
|
1 #Rscript
|
|
2
|
|
3 ###########################################
|
|
4 ## Mapping alpha and beta diversity ##
|
|
5 ###########################################
|
|
6
|
|
7 #####Packages : stars
|
|
8 # utils
|
|
9 # biodivmapr
|
|
10 # raster
|
|
11 # sf
|
|
12 # mapview
|
|
13 # leafpop
|
|
14 # RColorBrewer
|
|
15 # labdsv
|
|
16 # rgdal
|
|
17 # ggplot2
|
|
18 # gridExtra
|
|
19
|
|
20 #####Load arguments
|
|
21
|
|
22 args <- commandArgs(trailingOnly = TRUE)
|
|
23
|
|
24 #####Import the S2 data
|
|
25
|
|
26 if (length(args) < 1) {
|
|
27 stop("This tool needs at least 1 argument")
|
|
28 }else {
|
|
29 data_raster <- args[1]
|
|
30 rasterheader <- args[2]
|
|
31 data <- args[3]
|
|
32 # type of PCA:
|
|
33 # PCA: no rescaling of the data
|
|
34 # SPCA: rescaling of the data
|
|
35 typepca <- as.character(args[4])
|
|
36 alpha <- as.logical(args[5])
|
|
37 beta <- as.logical(args[6])
|
|
38 funct <- as.logical(args[7])
|
|
39 all <- as.logical(args[8])
|
|
40 nbcpu <- as.integer(args[9])
|
|
41 source(args[10])
|
|
42 }
|
|
43
|
|
44 ################################################################################
|
|
45 ## DEFINE PARAMETERS FOR DATASET TO BE PROCESSED ##
|
|
46 ################################################################################
|
|
47 if (data_raster == "") {
|
|
48 #Create a directory where to unzip your folder of data
|
|
49 dir.create("data_dir")
|
|
50 unzip(data, exdir = "data_dir")
|
|
51 # Path to raster
|
|
52 data_raster <- list.files("data_dir/results/Reflectance", pattern = "_Refl")
|
|
53 input_image_file <- file.path("data_dir/results/Reflectance", data_raster[1])
|
|
54 input_header_file <- file.path("data_dir/results/Reflectance", data_raster[2])
|
|
55
|
|
56 } else {
|
|
57 input_image_file <- file.path(getwd(), data_raster, fsep = "/")
|
|
58 input_header_file <- file.path(getwd(), rasterheader, fsep = "/")
|
|
59 }
|
|
60
|
|
61 ################################################################################
|
|
62 ## PROCESS IMAGE ##
|
|
63 ################################################################################
|
|
64 # 1- Filter data in order to discard non vegetated / shaded / cloudy pixels
|
|
65
|
|
66 print("PERFORM PCA ON RASTER")
|
|
67 pca_output <- biodivMapR::perform_PCA(Input_Image_File = input_image_file, Input_Mask_File = input_mask_file,
|
|
68 Output_Dir = output_dir, TypePCA = typepca, FilterPCA = filterpca, nbCPU = nbcpu, MaxRAM = maxram)
|
|
69
|
|
70 pca_files <- pca_output$PCA_Files
|
|
71 pix_per_partition <- pca_output$Pix_Per_Partition
|
|
72 nb_partitions <- pca_output$nb_partitions
|
|
73 # path for the updated mask
|
|
74 input_mask_file <- pca_output$MaskPath
|
|
75
|
|
76
|
|
77 selected_pcs <- seq(1, dim(raster::stack(input_image_file))[3])
|
|
78
|
|
79 selected_pcs <- all(selected_pcs)
|
|
80 ################################################################################
|
|
81 ## MAP ALPHA AND BETA DIVERSITY ##
|
|
82 ################################################################################
|
|
83 print("MAP SPECTRAL SPECIES")
|
|
84
|
|
85 kmeans_info <- biodivMapR::map_spectral_species(Input_Image_File = input_image_file, Output_Dir = output_dir, PCA_Files = pca_files, Input_Mask_File = input_mask_file, SelectedPCs = selected_pcs, Pix_Per_Partition = pix_per_partition, nb_partitions = nb_partitions, TypePCA = typepca, nbCPU = nbcpu, MaxRAM = maxram, nbclusters = nbclusters)
|
|
86
|
|
87 image_name <- tools::file_path_sans_ext(basename(input_image_file))
|
|
88 if (alpha == TRUE || beta == TRUE || all == TRUE) {
|
|
89 ## alpha
|
|
90 print("MAP ALPHA DIVERSITY")
|
|
91 index_alpha <- c("Shannon")
|
|
92 alpha_div <- biodivMapR::map_alpha_div(Input_Image_File = input_image_file, Output_Dir = output_dir, TypePCA = typepca, window_size = window_size, nbCPU = nbcpu, MaxRAM = maxram, Index_Alpha = index_alpha, nbclusters = nbclusters, FullRes = TRUE, LowRes = FALSE, MapSTD = FALSE)
|
|
93
|
|
94 alpha_zip <- file.path(output_dir, image_name, typepca, "ALPHA", "Shannon_10_Fullres.zip")
|
|
95 alpha_path <- file.path(output_dir, image_name, typepca, "ALPHA")
|
|
96 unzip(alpha_zip, exdir = alpha_path)
|
|
97 alpha_path <- file.path(output_dir, image_name, typepca, "ALPHA", "Shannon_10_Fullres")
|
|
98 alpha_raster <- raster::raster(alpha_path)
|
|
99 get_alpha <- convert_raster(alpha_raster)
|
|
100
|
|
101 if (alpha == TRUE || all == TRUE) {
|
|
102 colnames(get_alpha) <- c("Alpha", "longitude", "latitude")
|
|
103 plot_indices(get_alpha, titre = "Alpha")
|
|
104
|
|
105 write.table(get_alpha, file = "alpha.tabular", sep = "\t", dec = ".", na = " ", row.names = FALSE, col.names = TRUE, quote = FALSE)
|
|
106 }
|
|
107 if (beta == TRUE || all == TRUE) {
|
|
108 ## beta
|
|
109 print("MAP BETA DIVERSITY")
|
|
110 beta_div <- biodivMapR::map_beta_div(Input_Image_File = input_image_file, Output_Dir = output_dir, TypePCA = typepca, window_size = window_size, nb_partitions = nb_partitions, nbCPU = nbcpu, MaxRAM = maxram, nbclusters = nbclusters)
|
|
111
|
|
112 beta_path <- file.path(output_dir, image_name, typepca, "BETA", "BetaDiversity_BCdiss_PCO_10")
|
|
113 beta_raster <- raster::raster(beta_path)
|
|
114 get_beta <- convert_raster(beta_raster)
|
|
115
|
|
116 colnames(get_beta) <- c("Beta", "longitude", "latitude")
|
|
117 plot_indices(get_beta, titre = "Beta")
|
|
118
|
|
119 write.table(get_beta, file = "beta.tabular", sep = "\t", dec = ".", na = " ", row.names = FALSE, col.names = TRUE, quote = FALSE)
|
|
120 }
|
|
121 }
|
|
122
|
|
123
|
|
124 ################################################################################
|
|
125 ## COMPUTE ALPHA AND BETA DIVERSITY FROM FIELD PLOTS ##
|
|
126 ################################################################################
|
|
127
|
|
128 if (funct == TRUE || all == TRUE) {
|
|
129 mapper <- biodivMapR::map_functional_div(Original_Image_File = input_image_file, Functional_File = pca_files, Selected_Features = selected_pcs, Output_Dir = output_dir, window_size = window_size, nbCPU = nbcpu, MaxRAM = maxram, TypePCA = typepca, FullRes = TRUE, LowRes = FALSE, MapSTD = FALSE)
|
|
130
|
|
131 funct_zip <- file.path(output_dir, image_name, typepca, "FUNCTIONAL", "FunctionalDiversity_Map_MeanFilter_Fullres.zip")
|
|
132 funct_path <- file.path(output_dir, image_name, typepca, "FUNCTIONAL")
|
|
133 unzip(funct_zip, exdir = funct_path)
|
|
134 funct_path <- file.path(output_dir, image_name, typepca, "FUNCTIONAL", "FunctionalDiversity_Map_MeanFilter_Fullres")
|
|
135 funct_raster <- raster::raster(funct_path)
|
|
136 get_funct <- convert_raster(funct_raster)
|
|
137
|
|
138 colnames(get_funct) <- c("Functionnal", "longitude", "latitude")
|
|
139 plot_indices(get_funct, titre = "Functionnal")
|
|
140
|
|
141 write.table(get_funct, file = "Functionnal.tabular", sep = "\t", dec = ".", na = " ", row.names = FALSE, col.names = TRUE, quote = FALSE)
|
|
142 }
|