Simul_DIRSOL_avec_tout_RMQS_LULUCF.R

# Ce script permet de faire tourner les 4 scénarios aux deux résolutions attendues, en prenant en compte
#  tous les points du RMQS et en utilisant les données d'occupation du sol lulucf

titreCarte = "Simulation du réseau de mesures de la SML"

# initialisations --------
source("dirsol/0_InitialisationsLULCF.R")


# chargement des données RMQS -----

load("/home/lpotel/smlsamplingeu_louis/data/data_ETM_C_H1.RData")

#Données ponctuelles

smn = data_ETM_C_H1%>%
  dplyr::select(id_site,x_reel,y_reel) %>%
  st_as_sf(coords = 2:3,
           crs = 2154) %>%
  st_transform(crs(stack))

sample3 <-  terra::extract(stack,
                           smn,
                           xy = TRUE
)

sample3 <- as.data.frame(sample3)
sample3 <- sample3[complete.cases(sample3),]

sample3 <- sample3 %>%
  rename(FR_metzger = climate_metzeger,
         FR_regions = code_supra,
         FR_IGCS = ger,
         FR_WRB = WRBLEV1)
 
#remove the non soil
sample3 <- filter(sample3, !FR_WRB %in% excludeWRB )


# saveRDS( sample3, "~/smlsamplingeu_louis/data/octobre_avec_rmqs_sample3_lulucf.rds")
# sample3 <- readRDS("~/smlsamplingeu_louis/data/octobre_avec_rmqs_sample3_lulucf.rds")
# 



# en forçant la prise en compte de tous les points RMQS 
# Include points from existing surveys (forcing to keep existing points)
# Suppose we have a regular grid we want to preserve in total


##  (S1_lulucf_JRC) ------
scenario = "S1_lulucf_JRC"

ListVarScenario <- c("FR_landuse", 
                     "FR_regions", 
                     "FR_WRB",
                     "FR_metzger")

## 100m-------


domainsS1 <- 1e6 * stackb[["FR_landuse"]] +
  1e4 * stackb[["code_supra"]]  +
  1e2 *stackb[["WRBLEV1"]] +
  stackb[["climate_metzeger"]]


sampleRand <- sample1[rnd.lhs$index_samples,]

source("dirsol/ComputeSampleOptim.R")


# création du masque pour avoir les pixels pris en compte
mask_values <- terra::ifel( domainsS1 %in% as.numeric(tab1$domains[tab1$n>10]),
                            1, 0
                            )

titreCartetot = paste(titreCarte,"\n  ",scenario , "100m")
r = plotResuRMQS( mask_values,
              sample7,
              sampleRand ,
              sample9 ,
              titreCartetot ,
              stack500)


tmap_save(r , file = paste0("dirsol/figures/",scenario,"100.jpeg") )


## 500m ----

domainsS1500 <- 1e6 * stack500[["FR_landuse"]] + 
  1e4 * stack500[["code_supra"]]  +
  1e2 *stack500[["WRBLEV1"]] + 
  stack500[["climate_metzeger"]]

sampleRand <- sample1500[rnd.lhs500$index_samples,]
 

source("dirsol/ComputeSampleOptim.R")


# création du masque pour avoir les pixels pris en compte
mask_values <- terra::ifel( domainsS1500 %in% as.numeric(tab1$domains[tab1$n>10]),
                            1, 0
)

titreCartetot = paste(titreCarte,"\n  ",scenario , "500m")
r = plotResuRMQS( mask_values,
                  sample7,
                  sampleRand ,
                  sample9 ,
                  titreCartetot ,
                  stack500)


tmap_save(r , file = paste0("dirsol/figures/",scenario,"500.jpeg") )



# (S2_lulucf_JRC) ---------
scenario = "S2_lulucf_JRC"
ListVarScenario <- c("FR_landuse", 
                     "FR_regions", 
                     "FR_WRB")

## 100m------

domainsS2 <- 1e4 * stackb[["FR_landuse"]] +
  1e2 * stackb[["code_supra"]]  +
  stackb[["WRBLEV1"]]

sampleRand <- sample1[rnd.lhs$index_samples,]



source("dirsol/ComputeSampleOptim.R")


# création du masque pour avoir les pixels pris en compte
mask_values <- terra::ifel( domainsS2 %in% as.numeric(tab1$domains[tab1$n>10]),
                            1, 0
)

titreCartetot = paste(titreCarte,"\n  ",scenario , "100m")
r = plotResuRMQS( mask_values,
                  sample7,
                  sampleRand ,
                  sample9 ,
                  titreCartetot ,
                  stack500)


tmap_save(r , file = paste0("dirsol/figures/",scenario,"100.jpeg") )





##500m----

domainsS2500 <- 1e4 * stack500[["FR_landuse"]] + 
  1e2 * stack500[["code_supra"]]  +
  stack500[["WRBLEV1"]] 

sampleRand <- sample1500[rnd.lhs500$index_samples,]
 

source("dirsol/ComputeSampleOptim.R")


# création du masque pour avoir les pixels pris en compte
mask_values <- terra::ifel( domainsS2500 %in% as.numeric(tab1$domains[tab1$n>10]),
                            1, 0
)

titreCartetot = paste(titreCarte,"\n  ",scenario , "500m")
r = plotResuRMQS( mask_values,
                  sample7,
                  sampleRand ,
                  sample9 ,
                  titreCartetot ,
                  stack500)


tmap_save(r , file = paste0("dirsol/figures/",scenario,"500.jpeg") )


# (S1_lulucf_IGCS) -----

scenario = "S1_lulucf_IGCS"
ListVarScenario <- c("FR_landuse", "FR_regions", "FR_IGCS", "FR_metzger")

### 100m-------

domainsS3 <- 1e6 * stackb[["FR_landuse"]] + 1e4 * stackb[["code_supra"]]  +
  1e2 *stackb[["ger"]] + stackb[["climate_metzeger"]]

sampleRand <- sample1[rnd.lhs$index_samples,]
 
source("dirsol/ComputeSampleOptim.R")


# création du masque pour avoir les pixels pris en compte
mask_values <- terra::ifel( domainsS3 %in% as.numeric(tab1$domains[tab1$n>10]),
                            1, 0
)

titreCartetot = paste(titreCarte,"\n  ",scenario , "100m")
r = plotResuRMQS( mask_values,
                  sample7,
                  sampleRand ,
                  sample9 ,
                  titreCartetot ,
                  stack500)


tmap_save(r , file = paste0("dirsol/figures/",scenario,"100.jpeg") )

## 500m-----

domainsS3500 <- 1e6 * stack500[["FR_landuse"]] + 1e4 * stack500[["code_supra"]]  +
  1e2 *stack500[["ger"]] + stack500[["climate_metzeger"]]

sampleRand <- sample1500[rnd.lhs500$index_samples,]

source("dirsol/ComputeSampleOptim.R")


# création du masque pour avoir les pixels pris en compte
mask_values <- terra::ifel( domainsS3500 %in% as.numeric(tab1$domains[tab1$n>10]),
                            1, 0
)

titreCartetot = paste(titreCarte,"\n  ",scenario , "500m")
r = plotResuRMQS( mask_values,
                  sample7,
                  sampleRand ,
                  sample9 ,
                  titreCartetot ,
                  stack500)


tmap_save(r , file = paste0("dirsol/figures/",scenario,"500.jpeg") )


# (S2_lulucf_IGCS) ----------
scenario = "S2_lulucf_IGCS"
ListVarScenario <- c("FR_landuse", "FR_regions", "FR_IGCS")

## 100m-------

domainsS4 <- 1e4 * stackb[["FR_landuse"]] +
  1e2 * stackb[["code_supra"]]  +
  stackb[["ger"]]





sampleRand <- sample1[rnd.lhs$index_samples,]

source("dirsol/ComputeSampleOptim.R")


# création du masque pour avoir les pixels pris en compte
mask_values <- terra::ifel( domainsS3500 %in% as.numeric(tab1$domains[tab1$n>10]),
                            1, 0
)

titreCartetot = paste(titreCarte,"\n  ",scenario , "100m")
r = plotResuRMQS( mask_values,
                  sample7,
                  sampleRand ,
                  sample9 ,
                  titreCartetot ,
                  stack500)


tmap_save(r , file = paste0("dirsol/figures/",scenario,"100.jpeg") )

## 500m------------

domainsS4500 <- 1e4 * stack500[["FR_landuse"]] + 1e2 * stack500[["code_supra"]]  +
  stack500[["ger"]] 

sampleRand <- sample1500[rnd.lhs500$index_samples,]
 
source("dirsol/ComputeSampleOptim.R")


# création du masque pour avoir les pixels pris en compte
mask_values <- terra::ifel( domainsS3500 %in% as.numeric(tab1$domains[tab1$n>10]),
                            1, 0
)

titreCartetot = paste(titreCarte,"\n  ",scenario , "500m")
r = plotResuRMQS( mask_values,
                  sample7,
                  sampleRand ,
                  sample9 ,
                  titreCartetot ,
                  stack500)


tmap_save(r , file = paste0("dirsol/figures/",scenario,"500.jpeg") )