Mercurial > repos > ecology > srs_spectral_indices
comparison pca_raster.r @ 0:a8dabbf47e15 draft
planemo upload for repository https://github.com/Marie59/Sentinel_2A/srs_tools commit b32737c1642aa02cc672534e42c5cb4abe0cd3e7
| author | ecology |
|---|---|
| date | Mon, 09 Jan 2023 13:39:08 +0000 |
| parents | |
| children |
comparison
equal
deleted
inserted
replaced
| -1:000000000000 | 0:a8dabbf47e15 |
|---|---|
| 1 #Rscript | |
| 2 | |
| 3 ########################################### | |
| 4 ## Getting PCA raster ## | |
| 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 ## remotes::install_github("jbferet/biodivMapR") | |
| 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 typepca <- as.character(args[4]) | |
| 33 source(args[5]) | |
| 34 } | |
| 35 | |
| 36 ################################################################################ | |
| 37 ## DEFINE PARAMETERS FOR DATASET TO BE PROCESSED ## | |
| 38 ################################################################################ | |
| 39 # expected to be in ENVI HDR | |
| 40 | |
| 41 if (data_raster == "") { | |
| 42 #Create a directory where to unzip your folder of data | |
| 43 dir.create("data_dir") | |
| 44 unzip(data, exdir = "data_dir") | |
| 45 # Path to raster | |
| 46 data_raster <- list.files("data_dir/results/Reflectance", pattern = "_Refl") | |
| 47 input_image_file <- file.path("data_dir/results/Reflectance", data_raster[1]) | |
| 48 input_header_file <- file.path("data_dir/results/Reflectance", data_raster[2]) | |
| 49 | |
| 50 } else { | |
| 51 input_image_file <- file.path(getwd(), data_raster, fsep = "/") | |
| 52 input_header_file <- file.path(getwd(), rasterheader, fsep = "/") | |
| 53 } | |
| 54 | |
| 55 ################################################################################ | |
| 56 ## PROCESS IMAGE ## | |
| 57 ################################################################################ | |
| 58 # 1- Filter data in order to discard non vegetated / shaded / cloudy pixels | |
| 59 print("PERFORM PCA ON RASTER") | |
| 60 pca_output <- biodivMapR::perform_PCA(Input_Image_File = input_image_file, Input_Mask_File = input_mask_file, | |
| 61 Output_Dir = output_dir, TypePCA = typepca, FilterPCA = filterpca, nbCPU = nbcpu, MaxRAM = maxram) | |
| 62 | |
| 63 | |
| 64 pca_path <- file.path(output_dir, basename(data_raster), typepca, "PCA", "OutputPCA_8_PCs") | |
| 65 pca_raster <- raster::raster(pca_path) | |
| 66 get_pca <- convert_raster(pca_raster) | |
| 67 | |
| 68 colnames(get_pca) <- c("PCA", "longitude", "latitude") | |
| 69 plot_indices(get_pca, titre = "PCA") | |
| 70 | |
| 71 write.table(get_pca, file = "PCA.tabular", sep = "\t", dec = ".", na = " ", row.names = FALSE, col.names = TRUE, quote = FALSE) | |
| 72 #### Get the raster layer files | |
| 73 pca_files <- file.path("RESULTS", basename(data_raster), typepca, "PCA") | |
| 74 to_dir_short <- output_dir | |
| 75 file.copy(pca_files, to_dir_short) #copy files from long to short paths |