Mercurial > repos > ecology > srs_diversity_maps
comparison functions.r @ 0:9adccd3da70c draft default tip
planemo upload for repository https://github.com/Marie59/Sentinel_2A/srs_tools commit b32737c1642aa02cc672534e42c5cb4abe0cd3e7
author | ecology |
---|---|
date | Mon, 09 Jan 2023 13:37:37 +0000 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:9adccd3da70c |
---|---|
1 #Rscript | |
2 | |
3 ######################### | |
4 ## Functions ## | |
5 ######################### | |
6 | |
7 #####Packages : raster | |
8 # sp | |
9 # ggplot2 | |
10 | |
11 ####Set paramaters for all tools using BiodivMapR | |
12 | |
13 # path for the Mask raster corresponding to image to process | |
14 # expected to be in ENVI HDR format, 1 band, integer 8bits | |
15 # expected values in the raster: 0 = masked, 1 = selected | |
16 # set to FALSE if no mask available | |
17 input_mask_file <- FALSE | |
18 | |
19 # relative or absolute path for the Directory where results will be stored | |
20 # For each image processed, a subdirectory will be created after its name | |
21 output_dir <- "RESULTS" | |
22 | |
23 # SPATIAL RESOLUTION | |
24 # resolution of spatial units for alpha and beta diversity maps (in pixels), relative to original image | |
25 # if Res.Map = 10 for images with 10 m spatial resolution, then spatial units will be 10 pixels x 10m = 100m x 100m surfaces | |
26 # rule of thumb: spatial units between 0.25 and 4 ha usually match with ground data | |
27 # too small window_size results in low number of pixels per spatial unit, hence limited range of variation of diversity in the image | |
28 window_size <- 10 | |
29 | |
30 # PCA FILTERING: Set to TRUE if you want second filtering based on PCA outliers to be processed. Slower | |
31 filterpca <- TRUE | |
32 | |
33 ################################################################################ | |
34 ## DEFINE PARAMETERS FOR METHOD ## | |
35 ################################################################################ | |
36 nbcpu <- 4 | |
37 maxram <- 0.5 | |
38 nbclusters <- 50 | |
39 | |
40 ################################################################################ | |
41 ## PROCESS IMAGE ## | |
42 ################################################################################ | |
43 # 1- Filter data in order to discard non vegetated / shaded / cloudy pixels | |
44 ndvi_thresh <- 0.5 | |
45 blue_thresh <- 500 | |
46 nir_thresh <- 1500 | |
47 continuum_removal <- TRUE | |
48 | |
49 | |
50 | |
51 #### Convert raster to dataframe | |
52 | |
53 # Convert raster to SpatialPointsDataFrame | |
54 convert_raster <- function(data_raster) { | |
55 r_pts <- raster::rasterToPoints(data_raster, spatial = TRUE) | |
56 | |
57 # reproject sp object | |
58 geo_prj <- "+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0" | |
59 r_pts <- sp::spTransform(r_pts, sp::CRS(geo_prj)) | |
60 | |
61 | |
62 # Assign coordinates to @data slot, display first 6 rows of data.frame | |
63 r_pts@data <- data.frame(r_pts@data, longitude = sp::coordinates(r_pts)[, 1], | |
64 latitude = sp::coordinates(r_pts)[, 2]) | |
65 | |
66 return(r_pts@data) | |
67 } | |
68 | |
69 | |
70 #### Potting indices | |
71 | |
72 plot_indices <- function(data, titre) { | |
73 graph_indices <- ggplot2::ggplot(data) + | |
74 ggplot2::geom_point(ggplot2::aes_string(x = data[, 2], y = data[, 3], color = data[, titre]), shape = "square", size = 2) + ggplot2::scale_colour_gradient(low = "blue", high = "orange", na.value = "grey50") + | |
75 ggplot2::xlab("Longitude") + ggplot2::ylab("Latitude") + | |
76 ggplot2::theme(axis.text.x = ggplot2::element_text(angle = 90, vjust = 0.5, hjust = 1), plot.title = ggplot2::element_text(color = "black", size = 12, face = "bold")) + ggplot2::ggtitle(titre) | |
77 | |
78 ggplot2::ggsave(paste0(titre, ".png"), graph_indices, width = 12, height = 10, units = "cm") | |
79 } |