Skip to content
GitLab
Select Git revision
  • master default
Find file BlameHistoryPermalink
map.R 37.91 KiB
1 
# \\\
2 
# Copyright 2021-2022 Louis Héraut*1
3 
#
4 
# *1   INRAE, France
5 
#      louis.heraut@inrae.fr
6 
#
7 
# This file is part of ash R toolbox.
8 
#
9 
# ash R toolbox is free software: you can redistribute it and/or modify
10 
# it under the terms of the GNU General Public License as published by
11 
# the Free Software Foundation, either version 3 of the License, or (at
12 
# your option) any later version.
13 
#
14 
# ash R toolbox is distributed in the hope that it will be useful, but
15 
# WITHOUT ANY WARRANTY; without even the implied warranty of
16 
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
17 
# General Public License for more details.
18 
#
19 
# You should have received a copy of the GNU General Public License
20 
# along with ash R toolbox.  If not, see <https://www.gnu.org/licenses/>.
21 
# ///
22 
#
23 
#
24 
# plotting/map.R
25 
#
26 
# Deals with the creation of a map for presenting the trend analysis of hydrological variables
27
28
29 
## 1. MAP PANEL ______________________________________________________
30 
# Generates a map plot of the tendancy of a hydrological variable
31 
map_panel = function (list_df2plot, df_meta, df_shapefile, idPer_trend=1,
32 
                      trend_period, mean_period, colorForce=FALSE,
33 
                      outdirTmp='', codeLight=NULL,
34 
                      margin=NULL, showSea=TRUE,
35 
                      foot_note=FALSE,
36 
                      foot_height=0, resources_path=NULL,
37 
                      logo_dir=NULL,
38 
                      AEAGlogo_file=NULL, INRAElogo_file=NULL,
39 
                      FRlogo_file=NULL, df_page=NULL,
40 
                      verbose=TRUE) {
41
42 
    # Extract shapefiles
43 
    df_france = df_shapefile$france
44 
    df_bassin = df_shapefile$bassin
45 
    df_subbassin = df_shapefile$subbassin
46 
    df_river = df_shapefile$river
47
48 
    # Number of variable/plot
49 
    nbp = length(list_df2plot)
50
51 
    # Get all different stations code
52 
    Code = levels(factor(df_meta$code))
53 
    nCode = length(Code)
54
55 
    if (!is.null(trend_period)) {
56 
        # Convert 'trend_period' to list
57 
        trend_period = as.list(trend_period)
58 
        # Number of trend period
59 
        nPeriod_trend = length(trend_period)
60
61 
        # Extracts the min and the max of the mean trend for all the station
62 
        res = short_trendExtremes(list_df2plot, Code, nPeriod_trend, nbp, nCode, colorForce)
63 
        minTrendValue = res$min
64 
        maxTrendValue = res$max
65 
    }
66
67 
    # If there is a 'mean_period'
68 
    if (!is.null(mean_period)) {
69 
        # Convert 'mean_period' to list
70 
        mean_period = as.list(mean_period)
7172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140
# Number of mean period nPeriod_mean = length(mean_period) res = short_meanExtremes(list_df2plot, Code, nPeriod_mean, nbp, nCode) minBreakValue = res$min maxBreakValue = res$max breakValue_code = res$value } else { nPeriod_mean = 1 } # Number of ticks for the colorbar nbTick = 10 for (j in 1:nPeriod_mean) { # For all variable for (i in 1:nbp) { # If there is a specified station code to highlight (mini map) # and there has already been one loop if ((i > 1 | j > 1) & !is.null(codeLight)) { # Stop the for loop over the variable break } # Extracts the variable of the plot var = list_df2plot[[i]]$var # Extracts the type of variable of the plot type = list_df2plot[[i]]$type # Explanations about the variable glose = list_df2plot[[i]]$glose # Createsa name for the map if (j > 1) { outname = paste('map_d', var, sep='') } else { outname = paste('map_', var, sep='') } n_loop = i + nbp*(j-1) N_loop = nbp*nPeriod_mean # If there is the verbose option if (verbose) { if (j > 1) { mapName = 'difference' } else { mapName = 'tendence' } # Prints the name of the map print(paste('Map of ', mapName, ' for : ', var, " (", round(n_loop/N_loop*100, 0), " %)", sep='')) } # If there is no specified station code to highlight # (mini map) if (is.null(codeLight)) { # Sets the size of the countour sizefr = 0.45 sizebs = 0.4 sizerv = 0.3 } else { sizefr = 0.35 sizebs = 0.3 sizerv = 0.2 } # Stores the coordonate system cf = coord_fixed() # Makes it the default one to remove useless warning
141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210
cf$default = TRUE # Open a new plot with the personalise theme map = ggplot() + theme_void() + # theme( # # Ticks marker # axis.ticks.x=element_line(color='grey75', size=0.3), # axis.ticks.y=element_line(color='grey75', size=0.3), # # Ticks label # axis.text.x=element_text(color='grey40'), # axis.text.y=element_text(color='grey40')) + # Fixed coordinate system (remove useless warning) cf + # Plot the background of France geom_polygon(data=df_france, aes(x=long, y=lat, group=group), color=NA, fill="grey97") # If the river shapefile exists if (!is.null(df_river)) { # Plot the river map = map + geom_path(data=df_river, aes(x=long, y=lat, group=group), color="grey85", size=sizerv) } map = map + # Plot the hydrological basin geom_polygon(data=df_bassin, aes(x=long, y=lat, group=group), color="grey70", fill=NA, size=sizebs) + # Plot the hydrological sub-basin geom_polygon(data=df_subbassin, aes(x=long, y=lat, group=group), color="grey70", fill=NA, size=sizebs) + # Plot the countour of France geom_polygon(data=df_france, aes(x=long, y=lat, group=group), color="grey40", fill=NA, size=sizefr) if (is.null(codeLight)) { xBasin = c(410000, 520000, 630000, 620000, 510000, 450000, 390000) yBasin = c(6280000, 6290000, 6320000, 6385000, 6450000, 6530000, 6361000) nameBasin = c('Adour', 'Garonne', 'Tarn-Aveyron', 'Lot', 'Dordogne', 'Charente', 'Fleuves-\nCôtiers') nBasin = length(xBasin) plot_basin = tibble(x=xBasin, y=yBasin, label=nameBasin) map = map + geom_shadowtext(data=plot_basin, aes(x=x, y=y, label=label), fontface="bold", color="grey85", bg.colour="grey97", hjust=0.5, vjust=0.5, size=5) } # If the sea needs to be shown if (showSea) { # Leaves space around the France xlim = c(295000, 790000) ylim = c(6125000, 6600000)
211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280
# Otherwise } else { # Leaves minimal space around France xlim = c(305000, 790000) ylim = c(6135000, 6600000) } # If there is no specified station code to highlight (mini map) if (is.null(codeLight)) { # Sets a legend scale start xmin = gpct(4, xlim, shift=TRUE) # Sets graduations xint = c(0, 10*1E3, 50*1E3, 100*1E3) # Sets the y postion ymin = gpct(5, ylim, shift=TRUE) # Sets the height of graduations ymax = ymin + gpct(1, ylim) # Size of the value size = 3 # Size of the 'km' unit sizekm = 2.5 # If there is a specified station code } else { # Same but with less graduation and smaller size xmin = gpct(2, xlim, shift=TRUE) xint = c(0, 100*1E3) ymin = gpct(1, ylim, shift=TRUE) ymax = ymin + gpct(3, ylim) size = 2 sizekm = 1.5 } map = map + # Adds the base line of the scale geom_line(aes(x=c(xmin, max(xint)+xmin), y=c(ymin, ymin)), color="grey40", size=0.2) + # Adds the 'km' unit annotate("text", x=max(xint)+xmin+gpct(1, xlim), y=ymin, vjust=0, hjust=0, label="km", color="grey40", size=sizekm) # For all graduations for (x in xint) { map = map + # Draws the tick annotate("segment", x=x+xmin, xend=x+xmin, y=ymin, yend=ymax, color="grey40", size=0.2) + # Adds the value annotate("text", x=x+xmin, y=ymax+gpct(0.5, ylim), vjust=0, hjust=0.5, label=x/1E3, color="grey40", size=size) } map = map + # Allows to crop shapefile without graphical problem coord_sf(xlim=xlim, ylim=ylim, expand=FALSE) # If there is no margins specified if (is.null(margin)) { # Sets all margins to 0 map = map + theme(plot.margin=margin(t=0, r=0, b=0, l=0, unit="mm")) # Otherwise } else { # Sets margins to the given ones
281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350
map = map + theme(plot.margin=margin) } # Blank vector to store data about station lon = c() lat = c() fill = c() shape = c() Value = c() OkVal = c() # For all code for (k in 1:nCode) { # Gets the code code = Code[k] if (j > 1) { value = breakValue_code[j, i, k] minValue = minBreakValue[j, i] maxValue = maxBreakValue[j, i] pVal = 0 } else if (is.null(trend_period)) { value = NA minValue = NULL maxValue = NULL pVal = 0 } else { # Extracts the data corresponding to the # current variable df_data = list_df2plot[[i]]$data # Extracts the trend corresponding to the # current variable df_trend = list_df2plot[[i]]$trend # Gets the risk of the test alpha = list_df2plot[[i]]$alpha # Extracts the data corresponding to the code df_data_code = df_data[df_data$code == code,] # Extracts the trend corresponding to the code df_trend_code = df_trend[df_trend$code == code,] # Extract start and end of trend periods Start = df_trend_code$period_start[idPer_trend] End = df_trend_code$period_end[idPer_trend] # Extracts the corresponding data for the period df_data_code_per = df_data_code[df_data_code$Date >= Start & df_data_code$Date <= End,] # Same for trend df_trend_code_per = df_trend_code[df_trend_code$period_start == Start & df_trend_code$period_end == End,] # Computes the number of trend analysis selected Ntrend = nrow(df_trend_code_per) # If there is more than one trend on the same period if (Ntrend > 1) { # Takes only the first because they are similar df_trend_code_per = df_trend_code_per[1,] } # If it is a flow variable if (type == 'sévérité') { # Computes the mean of the data on the period dataMean = mean(df_data_code_per$Value, na.rm=TRUE)