From 2b170c4c2cd85df822a4179d427d58529de1c8f5 Mon Sep 17 00:00:00 2001 From: "louis.heraut" <louis.heraut@inrae.fr> Date: Fri, 7 Jan 2022 20:44:05 +0100 Subject: [PATCH] Cleaned and commented --- plotting/info.R | 140 ------- plotting/layout.R | 245 ++++++++---- plotting/map.R | 1 + plotting/matrix.R | 8 +- plotting/other.R | 239 ------------ plotting/time.R | 949 ---------------------------------------------- 6 files changed, 177 insertions(+), 1405 deletions(-) delete mode 100644 plotting/info.R delete mode 100644 plotting/other.R delete mode 100644 plotting/time.R diff --git a/plotting/info.R b/plotting/info.R deleted file mode 100644 index 9759946..0000000 --- a/plotting/info.R +++ /dev/null @@ -1,140 +0,0 @@ -# \\\ -# Copyright 2021-2022 Louis Héraut*1 -# -# *1 INRAE, France -# louis.heraut@inrae.fr -# -# This file is part of ash R toolbox. -# -# ash R toolbox is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or (at -# your option) any later version. -# -# ash R toolbox is distributed in the hope that it will be useful, but -# WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -# General Public License for more details. -# -# You should have received a copy of the GNU General Public License -# along with ash R toolbox. If not, see <https://www.gnu.org/licenses/>. -# /// -# -# -# plotting/info.R -# -# - - -info_panel = function(list_df2plot, df_meta, df_shapefile, codeLight, df_data_code=NULL) { - - if (!is.null(df_data_code)) { - hyd = hydrogramme(df_data_code, - margin=margin(t=3, r=0, b=0, l=5, unit="mm")) - } else { - hyd = void - } - - # databin = list_df2plot[[1]]$data[list_df2plot[[1]]$data$code == codeLight,] - # yearLast = format(databin$Date[nrow(databin)], "%Y") - # yearFirst = format(databin$Date[1], "%Y") - # Nyear = as.numeric(yearLast) - as.numeric(yearFirst) + 1 - - map = map_panel(list_df2plot, - df_meta, - df_shapefile=df_shapefile, - codeLight=codeLight, - margin=margin(t=5, r=2, b=0, l=0, unit="mm"), - showSea=FALSE, - verbose=FALSE) - - df_meta_code = df_meta[df_meta$code == codeLight,] - - nom = df_meta_code$nom - nom = gsub("-", "- ", nom) - - duration = as.numeric(format(as.Date(df_meta_code$fin), "%Y")) - - as.numeric(format(as.Date(df_meta_code$debut), "%Y")) - debut = format(as.Date(df_meta_code$debut), "%d/%m/%Y") - fin = format(as.Date(df_meta_code$fin), "%d/%m/%Y") - - text1 = paste( - "<b>", codeLight, '</b> - ', nom, - sep='') - - text2 = paste( - "<b>", - "Gestionnaire : ", df_meta_code$gestionnaire, "<br>", - "Région hydro : ", df_meta_code$region_hydro, - "</b>", - sep='') - - text3 = paste( - "<b>", - "Superficie : ", df_meta_code$surface_km2_BH, " [km<sup>2</sup>] <br>", - "Altitude : ", df_meta_code$altitude_m_BH, " [m]<br>", - "X = ", df_meta_code$L93X_m_BH, " [m ; Lambert 93]<br>", - "Y = ", df_meta_code$L93Y_m_BH, " [m ; Lambert 93]", - "</b>", - sep='') - - text4 = paste( - "<b>", - "Date de début : ", debut, "<br>", - "Date de fin : ", fin, "<br>", - "Nombre d'années : ", duration, " [ans]", "<br>", - "Taux de lacunes : ", signif(df_meta_code$tLac100, 2), " [%]", - "</b>", - sep='') - - gtext1 = richtext_grob(text1, - x=0, y=1, - margin=unit(c(t=5, r=5, b=10, l=5), "mm"), - hjust=0, vjust=1, - gp=gpar(col="#00A3A8", fontsize=14)) - - gtext2 = richtext_grob(text2, - x=0, y=1, - margin=unit(c(t=0, r=0, b=0, l=5), "mm"), - hjust=0, vjust=1, - gp=gpar(col="grey20", fontsize=8)) - - gtext3 = richtext_grob(text3, - x=0, y=1, - margin=unit(c(t=0, r=0, b=0, l=5), "mm"), - hjust=0, vjust=1, - gp=gpar(col="grey20", fontsize=9)) - - gtext4 = richtext_grob(text4, - x=0, y=1, - margin=unit(c(t=0, r=0, b=0, l=0), "mm"), - hjust=0, vjust=1, - gp=gpar(col="grey20", fontsize=9)) - - - void = void + - theme(plot.background=element_rect(fill=NA, color="#EC4899"), - plot.margin=margin(t=0, r=0, b=0, l=0, unit="mm")) - - P = list(gtext1, gtext2, gtext3, gtext4, hyd, map) - # P = list(void, void, void, void, void, void, void) - - LM = matrix(c(1, 1, 1, 6, - # 7, 7, 7, 6, - 2, 2, 5, 6, - # 2, 2, 5, 6, - 3, 4, 5, 6, - 3, 4, 5, 6), - nrow=4, - byrow=TRUE) - - heights = rep(1, times=nrow(LM)) - # heights[2] = 0.1 - heights[2] = 0.8 - - plot = grid.arrange(grobs=P, - layout_matrix=LM, - heights=heights) - - return(plot) -} diff --git a/plotting/layout.R b/plotting/layout.R index 2cdc9a5..9c53437 100644 --- a/plotting/layout.R +++ b/plotting/layout.R @@ -40,11 +40,10 @@ library(rgdal) library(shadowtext) # Sourcing R file -source('plotting/time.R', encoding='latin1') -source('plotting/info.R', encoding='latin1') +source('plotting/datasheet.R', encoding='latin1') source('plotting/map.R', encoding='latin1') source('plotting/matrix.R', encoding='latin1') -source('plotting/other.R', encoding='latin1') +source('plotting/break.R', encoding='latin1') ## 1. PERSONALISATION @@ -84,22 +83,64 @@ theme_ash = ) ### 1.2. Color palette -palette_perso = c('#0f3b57', +palette_perso = c('#0f3b57', # cold '#1d7881', '#80c4a9', - '#e2dac6', #mid + '#e2dac6', # mid '#fadfad', '#d08363', - '#7e392f') + '#7e392f') # hot -## 2. LAYOUT -datasheet_layout = function (df_data, df_meta, layout_matrix, isplot=c('datasheet', 'matrix', 'map'), figdir='', filedir_opt='', filename_opt='', variable='', df_trend=NULL, p_threshold=0.1, unit2day=365.25, type='', trend_period=NULL, mean_period=NULL, axis_xlim=NULL, missRect=FALSE, time_header=NULL, info_header=TRUE, info_ratio=1, time_ratio=2, var_ratio=3, df_shapefile=NULL) { - +## 2. USEFUL GENERICAL PLOT +### 2.1. Void plot +# A plot completly blank +void = ggplot() + geom_blank(aes(1,1)) + + theme( + plot.background = element_blank(), + panel.grid.major = element_blank(), + panel.grid.minor = element_blank(), + panel.border = element_blank(), + panel.background = element_blank(), + axis.title.x = element_blank(), + axis.title.y = element_blank(), + axis.text.x = element_blank(), + axis.text.y = element_blank(), + axis.ticks = element_blank(), + axis.line = element_blank() + ) + +### 2.2. Circle +# Allow to draw circle in ggplot2 with a radius and a center position +gg_circle = function(r, xc, yc, color="black", fill=NA, ...) { + x = xc + r*cos(seq(0, pi, length.out=100)) + ymax = yc + r*sin(seq(0, pi, length.out=100)) + ymin = yc + r*sin(seq(0, -pi, length.out=100)) + annotate("ribbon", x=x, ymin=ymin, ymax=ymax, color=color, + fill=fill, ...) +} + + +## 3. LAYOUT +# Generates a PDF that gather datasheets, map and summarize matrix about the trend analyses realised on selected stations +datasheet_layout = function (df_data, df_meta, layout_matrix, + isplot=c('datasheet', 'matrix', 'map'), + figdir='', filedir_opt='', filename_opt='', + variable='', df_trend=NULL, + p_threshold=0.1, unit2day=365.25, type='', + trend_period=NULL, mean_period=NULL, + axis_xlim=NULL, missRect=FALSE, + time_header=NULL, info_header=TRUE, + info_ratio=1, time_ratio=2, var_ratio=3, + df_shapefile=NULL) { + + # Name of the document outfile = "Panels" + # If there is an option to mention in the filename it adds it if (filename_opt != '') { outfile = paste(outfile, '_', filename_opt, sep='') } + # Add the 'pdf' extensionto the name outfile = paste(outfile, '.pdf', sep='') # If there is not a dedicated figure directory it creats one @@ -108,16 +149,22 @@ datasheet_layout = function (df_data, df_meta, layout_matrix, isplot=c('datashee dir.create(outdir) } + # Names of a temporary directory to store all the independent pages outdirTmp = file.path(outdir, 'tmp') + # Creates it if it does not exist if (!(file.exists(outdirTmp))) { dir.create(outdirTmp) + # If it already exists it deletes the pre-existent directory + # and recreates one } else { unlink(outdirTmp, recursive=TRUE) dir.create(outdirTmp) } + # Number of variable studied nbp = length(df_data) + # Convert data tibble to list of tibble if it is not the case if (all(class(df_data) != 'list')) { df_data = list(df_data) } @@ -125,15 +172,20 @@ datasheet_layout = function (df_data, df_meta, layout_matrix, isplot=c('datashee if (all(class(df_trend) != 'list')) { df_trend = list(df_trend) if (length(df_trend) == 1) { + df_trend = replicate(nbp, df_trend) }} if (all(class(p_threshold) != 'list')) { p_threshold = list(p_threshold) + # If there is only one value if (length(p_threshold) == 1) { + # Replicates the value the number of times that there + # is of studied variables p_threshold = replicate(nbp, p_threshold) }} - + + # Same if (all(class(unit2day) != 'list')) { unit2day = list(unit2day) if (length(unit2day) == 1) { @@ -152,38 +204,34 @@ datasheet_layout = function (df_data, df_meta, layout_matrix, isplot=c('datashee missRect = replicate(nbp, missRect) }} + # Creates a blank list to store all the data of each type of plot list_df2plot = vector(mode='list', length=nbp) - # minTrend = c() - # maxTrend = c() + # For all the type of graph / number of studied variables for (i in 1:nbp) { - + # Creates a list that gather all the info for one type of graph df2plot = list(data=df_data[[i]], trend=df_trend[[i]], p_threshold=p_threshold[[i]], unit2day=unit2day[[i]], type=type[[i]], missRect=missRect[[i]]) - - # okTrend = df_trend[[i]]$trend[df_trend[[i]]$p <= p_threshold[[i]]] - - # minTrend[i] = min(okTrend, na.rm=TRUE) - # maxTrend[i] = max(okTrend, na.rm=TRUE) - + # Stores it list_df2plot[[i]] = df2plot } - + # If datasheets needs to be plot if ('datasheet' %in% isplot) { - time_panel(list_df2plot, df_meta, trend_period, info_header=info_header, time_header=time_header, layout_matrix=layout_matrix, info_ratio=info_ratio, time_ratio=time_ratio, var_ratio=var_ratio, outdirTmp=outdirTmp) + datasheet_panel(list_df2plot, df_meta, trend_period, info_header=info_header, time_header=time_header, layout_matrix=layout_matrix, info_ratio=info_ratio, time_ratio=time_ratio, var_ratio=var_ratio, outdirTmp=outdirTmp) } - + # If summarize matrix needs to be plot if ('matrix' %in% isplot) { matrix_panel(list_df2plot, df_meta, trend_period, mean_period, slice=12, outdirTmp=outdirTmp, A3=TRUE) } - + + # If map needs to be plot if ('map' %in% isplot) { map_panel(list_df2plot, df_meta, @@ -193,170 +241,225 @@ datasheet_layout = function (df_data, df_meta, layout_matrix, isplot=c('datashee margin=margin(t=5, r=0, b=5, l=5, unit="mm")) } - # PDF combine + # Combine independant pages into one PDF pdf_combine(input=file.path(outdirTmp, list.files(outdirTmp)), output=file.path(outdir, outfile)) - # unlink(outdirTmp, recursive=TRUE) } -## 3. COLOR MANAGEMENT -### 3.1. Color on colorbar +## 4. COLOR MANAGEMENT +### 4.1. Color on colorbar +# Returns a color of a palette corresponding to a value included +# between the min and the max of the variable get_color = function (value, min, max, ncolor=256, palette_name='perso', reverse=FALSE) { + + # If the value is a NA return NA color + if (is.na(value)) { + return (NA) + } + # If the palette chosen is the personal ones if (palette_name == 'perso') { colorList = palette_perso + # Else takes the palette corresponding to the name given } else { colorList = brewer.pal(11, palette_name) } + # Gets the number of discrete colors in the palette nSample = length(colorList) - + # Recreates a continuous color palette palette = colorRampPalette(colorList)(ncolor) + # Separates it in the middle to have a cold and a hot palette Sample_hot = 1:(as.integer(nSample/2)+1) Sample_cold = (as.integer(nSample/2)+1):nSample - palette_hot = colorRampPalette(colorList[Sample_hot])(ncolor) palette_cold = colorRampPalette(colorList[Sample_cold])(ncolor) + # Reverses the palette if it needs to be if (reverse) { palette = rev(palette) palette_hot = rev(palette_hot) palette_cold = rev(palette_cold) } - if (is.na(value)) { - return (NA) - } - + # Computes the absolute max maxAbs = max(abs(max), abs(min)) - + + # If the value is negative if (value < 0) { + # Gets the relative position of the value in respect + # to its span idNorm = (value + maxAbs) / maxAbs + # The index corresponding id = round(idNorm*(ncolor - 1) + 1, 0) + # The associated color color = palette_cold[id] + # Same if it is a positive value } else { idNorm = value / maxAbs id = round(idNorm*(ncolor - 1) + 1, 0) color = palette_hot[id] } - return(color) } -### 3.2. Colorbar +### 4.2. Colorbar +# Returns the colorbar but also positions, labels and colors of some +# ticks along it get_palette = function (min, max, ncolor=256, palette_name='perso', reverse=FALSE, nbTick=10) { + # If the palette chosen is the personal ones if (palette_name == 'perso') { colorList = palette_perso + # Else takes the palette corresponding to the name given } else { colorList = brewer.pal(11, palette_name) } + # Gets the number of discrete colors in the palette nSample = length(colorList) - + # Recreates a continuous color palette palette = colorRampPalette(colorList)(ncolor) + # Separates it in the middle to have a cold and a hot palette Sample_hot = 1:(as.integer(nSample/2)+1) Sample_cold = (as.integer(nSample/2)+1):nSample - palette_hot = colorRampPalette(colorList[Sample_hot])(ncolor) palette_cold = colorRampPalette(colorList[Sample_cold])(ncolor) + # Reverses the palette if it needs to be if (reverse) { palette = rev(palette) palette_hot = rev(palette_hot) palette_cold = rev(palette_cold) - } + } + # If the min and the max are below zero if (min < 0 & max < 0) { + # The palette show is only the cold one paletteShow = palette_cold + # If the min and the max are above zero } else if (min > 0 & max > 0) { + # The palette show is only the hot one paletteShow = palette_hot + # Else it is the entire palette that is shown } else { paletteShow = palette } + # The position of ticks is between 0 and 1 posTick = seq(0, 1, length.out=nbTick) - + # Blank vector to store corresponding labels and colors labTick = c() colTick = c() + # For each tick for (i in 1:nbTick) { + # Computes the graduation between the min and max lab = (i-1)/(nbTick-1) * (max - min) + min - labTick = c(labTick, lab) + # Gets the associated color col = get_color(lab, min=min, max=max, ncolor=ncolor, palette_name=palette_name, reverse=reverse) - colTick = c(colTick, col) + # Stores them + labTick = c(labTick, lab) + colTick = c(colTick, col) } - - return(list(palette=paletteShow, posTick=posTick, - labTick=labTick, colTick=colTick)) + # List of results + res = list(palette=paletteShow, posTick=posTick, + labTick=labTick, colTick=colTick) + return(res) } -### 3.3. Palette tester +### 4.3. Palette tester +# Allows to display the current personal palette palette_tester = function (n=256) { + # An arbitrary x vector X = 1:n + # All the same arbitrary y position to create a colorbar Y = rep(0, times=n) + # Recreates a continuous color palette palette = colorRampPalette(palette_perso)(n) + # Open a plot p = ggplot() + - - theme( - plot.background = element_blank(), - panel.grid.major = element_blank(), - panel.grid.minor = element_blank(), - panel.border = element_blank(), - panel.background = element_blank(), - axis.title.x = element_blank(), - axis.title.y = element_blank(), - axis.text.x = element_blank(), - axis.text.y = element_blank(), - axis.ticks = element_blank(), - axis.line = element_blank() - ) + - + # Make the theme blank + theme( + plot.background = element_blank(), + panel.grid.major = element_blank(), + panel.grid.minor = element_blank(), + panel.border = element_blank(), + panel.background = element_blank(), + axis.title.x = element_blank(), + axis.title.y = element_blank(), + axis.text.x = element_blank(), + axis.text.y = element_blank(), + axis.ticks = element_blank(), + axis.line = element_blank() + ) + + # Plot the palette geom_line(aes(x=X, y=Y), color=palette[X], size=60) + scale_y_continuous(expand=c(0, 0)) - + + # Saves the plot ggsave(plot=p, filename=paste('palette_test', '.pdf', sep=''), width=10, height=10, units='cm', dpi=100) } -## 4. OTHER -### 4.1. Number formatting +## 5. OTHER TOOLS +### 5.1. Number formatting +# Returns the power of ten of the scientific expression of a value get_power = function (value) { + + # Do not care about the sign + value = abs(value) - if (value > 1) { + # If the value is greater than one + if (value >= 1) { + # The magnitude is the number of character of integer part + # of the value minus one power = nchar(as.character(as.integer(value))) - 1 + # If the value is less than one } else { + # Extract the decimal part dec = gsub('0.', '', as.character(value), fixed=TRUE) + # Number of decimal with zero ndec = nchar(dec) + # Number of decimal without zero nnum = nchar(as.character(as.numeric(dec))) + # Compute the power of ten associated power = -(ndec - nnum + 1) } - return(power) } -### 4.2. Pourcentage of variable +### 5.2. Pourcentage of variable +# Returns the value corresponding of a certain percentage of a +# data serie gpct = function (pct, L, ref=NULL, shift=FALSE) { - + + # If no reference for the serie is given if (is.null(ref)) { + # The minimum of the serie is computed minL = min(L, na.rm=TRUE) + # If a reference is specified } else { + # The reference is the minimum minL = ref } - + + # Gets the max maxL = max(L, na.rm=TRUE) + # And the span spanL = maxL - minL - + # Computes the value corresponding to the percentage xL = pct/100 * as.numeric(spanL) - + + # If the value needs to be shift by its reference if (shift) { xL = xL + minL } diff --git a/plotting/map.R b/plotting/map.R index 185cde0..564635a 100644 --- a/plotting/map.R +++ b/plotting/map.R @@ -26,6 +26,7 @@ # +## 1. MAP PANEL map_panel = function (list_df2plot, df_meta, df_shapefile, idPer=1, outdirTmp='', codeLight=NULL, margin=NULL, showSea=TRUE, verbose=TRUE) { diff --git a/plotting/matrix.R b/plotting/matrix.R index a76c9e8..69468df 100644 --- a/plotting/matrix.R +++ b/plotting/matrix.R @@ -26,6 +26,7 @@ # +## 1. MATRIX PANEL matrix_panel = function (list_df2plot, df_meta, trend_period, mean_period, slice=NULL, outdirTmp='', outnameTmp='matrix', title=NULL, A3=FALSE) { nbp = length(list_df2plot) @@ -770,17 +771,12 @@ matrix_panel = function (list_df2plot, df_meta, trend_period, mean_period, slice height = 21 dpi = 100 } - - # print('ggg') ggsave(plot=mat, path=outdirTmp, filename=paste(outnameTmp, '_', fL, imat, '.pdf', sep=''), - width=width, height=height, units='cm', dpi=dpi) - - # print('hhh') - + width=width, height=height, units='cm', dpi=dpi) } } } diff --git a/plotting/other.R b/plotting/other.R deleted file mode 100644 index 491513b..0000000 --- a/plotting/other.R +++ /dev/null @@ -1,239 +0,0 @@ -# \\\ -# Copyright 2021-2022 Louis Héraut*1 -# -# *1 INRAE, France -# louis.heraut@inrae.fr -# -# This file is part of ash R toolbox. -# -# ash R toolbox is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or (at -# your option) any later version. -# -# ash R toolbox is distributed in the hope that it will be useful, but -# WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -# General Public License for more details. -# -# You should have received a copy of the GNU General Public License -# along with ash R toolbox. If not, see <https://www.gnu.org/licenses/>. -# /// -# -# -# plotting/other.R -# -# - - -hydrogramme = function (df_data_code, margin=NULL) { - - monthData = as.numeric(format(df_data_code$Date, "%m")) - - monthMean = c() - for (i in 1:12) { - data = df_data_code$Qm3s[monthData == i] - monthMean[i] = mean(data, na.rm=TRUE) - } - - monthNum = 1:12 - monthName = c("J", "F", "M", "A", "M", "J", - "J", "A", "S", "O", "N", "D") - # monthName = factor(monthName, levels=monthName) - - hyd = ggplot() + theme_ash + - - theme( - # plot.background=element_rect(fill=NA, color="#EC4899"), - panel.border=element_blank(), - axis.text.x=element_text(margin=unit(c(0, 0, 0, 0), "mm"), - vjust=1, hjust=0.5), - axis.ticks.x=element_blank(), - axis.line.y=element_line(color='grey85', size=0.3), - plot.title=element_text(size=8, vjust=-1.5, - hjust=-1E-3, color='grey40')) + - - ggtitle(bquote(bold('Q')~~'['*m^{3}*'.'*s^{-1}*']')) - - if (is.null(margin)) { - hyd = hyd + - theme(plot.margin=margin(t=0, r=0, b=0, l=0, unit="mm")) - } else { - hyd = hyd + - theme(plot.margin=margin) - } - - hyd = hyd + - - geom_bar(aes(x=monthNum, y=monthMean), - stat='identity', - fill="grey70", - width=0.75, size=0.2) + - - scale_x_continuous(breaks=monthNum, - labels=monthName, - limits=c(0, max(monthNum)+0.5), - expand=c(0, 0)) + - - scale_y_continuous(limits=c(0, max(monthMean)), - expand=c(0, 0)) - - return (hyd) -} - - - - - -histogram = function (data_bin, df_meta, figdir='', filedir_opt='') { - - # Get all different stations code - Code = levels(factor(df_meta$code)) - nCode = length(Code) - - # If there is not a dedicated figure directory it creats one - outdir = file.path(figdir, filedir_opt, sep='') - if (!(file.exists(outdir))) { - dir.create(outdir) - } - - datebreak = 10 - dateminbreak = 1 - - res_hist = hist(data_bin, breaks='years', plot=FALSE) - counts = res_hist$counts - counts_pct = counts/nCode * 100 - breaks = as.Date(res_hist$breaks) - mids = as.Date(res_hist$mids) - - p = ggplot() + theme_ash + - - theme(panel.grid.major.y=element_line(color='grey85', size=0.15), - axis.title.y=element_blank()) + - - geom_bar(aes(x=mids, y=counts_pct), - stat='identity', - fill="#00A3A8") + - - scale_x_date(date_breaks=paste(as.character(datebreak), - 'year', sep=' '), - date_minor_breaks=paste(as.character(dateminbreak), - 'year', sep=' '), - guide='axis_minor', - date_labels="%Y", - limits=c(min(data_bin)-years(0), - max(data_bin)+years(0)), - expand=c(0, 0)) + - - scale_y_continuous(limits=c(0, - max(counts_pct)*1.1), - expand=c(0, 0)) - - ggsave(plot=p, - path=outdir, - filename=paste('hist_break_date', '.pdf', sep=''), - width=10, height=10, units='cm', dpi=100) -} - - -cumulative = function (data_bin, df_meta, dyear=10, figdir='', filedir_opt='') { - - # Get all different stations code - Code = levels(factor(df_meta$code)) - nCode = length(Code) - - # If there is not a dedicated figure directory it creats one - outdir = file.path(figdir, filedir_opt, sep='') - if (!(file.exists(outdir))) { - dir.create(outdir) - } - - datebreak = 10 - dateminbreak = 1 - - res_hist = hist(data_bin, breaks='years', plot=FALSE) - counts = res_hist$counts - cumul = cumsum(counts) - cumul_pct = cumul/nCode * 100 - breaks = as.Date(res_hist$breaks) - mids = as.Date(res_hist$mids) - - mids = c(mids[1] - years(dyear), mids[1] - years(1), - mids, - mids[length(mids)] + years(dyear)) - cumul_pct = c(0, 0, - cumul_pct, - cumul_pct[length(cumul_pct)]) - - mids = mids + months(6) - - breaks = breaks + 1 - breaks = breaks[-length(breaks)] - - DB = c() - for (i in 1:length(breaks)) { - DB = c(DB, rep(breaks[i], times=counts[i])) - } - q50 = as.Date(quantile(DB, probs=0.5)) + years(1) - - print(paste('mediane :', q50)) - - p = ggplot() + theme_ash + - - theme(panel.grid.major.y=element_line(color='grey85', size=0.15), - axis.title.y=element_blank()) + - - geom_line(aes(x=mids, y=cumul_pct), - color="#00A3A8") + - - geom_line(aes(x=c(q50, q50), y=c(0, 100)), - color="wheat", - lty='dashed') + - - scale_x_date(date_breaks=paste(as.character(datebreak), - 'year', sep=' '), - date_minor_breaks=paste(as.character(dateminbreak), - 'year', sep=' '), - guide='axis_minor', - date_labels="%Y", - limits=c(min(mids)-years(0), - max(mids)+years(0)), - expand=c(0, 0)) + - - scale_y_continuous(limits=c(-1, 101), - expand=c(0, 0)) - - ggsave(plot=p, - path=outdir, - filename=paste('cumul_break_date', '.pdf', sep=''), - width=10, height=10, units='cm', dpi=100) -} - - - - -void = ggplot() + geom_blank(aes(1,1)) + - theme( - plot.background = element_blank(), - panel.grid.major = element_blank(), - panel.grid.minor = element_blank(), - panel.border = element_blank(), - panel.background = element_blank(), - axis.title.x = element_blank(), - axis.title.y = element_blank(), - axis.text.x = element_blank(), - axis.text.y = element_blank(), - axis.ticks = element_blank(), - axis.line = element_blank() - ) - - -gg_circle = function(r, xc, yc, color="black", fill=NA, ...) { - x = xc + r*cos(seq(0, pi, length.out=100)) - ymax = yc + r*sin(seq(0, pi, length.out=100)) - ymin = yc + r*sin(seq(0, -pi, length.out=100)) - annotate("ribbon", x=x, ymin=ymin, ymax=ymax, color=color, fill=fill, ...) -} - - - diff --git a/plotting/time.R b/plotting/time.R deleted file mode 100644 index 3cbe8f2..0000000 --- a/plotting/time.R +++ /dev/null @@ -1,949 +0,0 @@ -# \\\ -# Copyright 2021-2022 Louis Héraut*1 -# -# *1 INRAE, France -# louis.heraut@inrae.fr -# -# This file is part of ash R toolbox. -# -# ash R toolbox is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or (at -# your option) any later version. -# -# ash R toolbox is distributed in the hope that it will be useful, but -# WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -# General Public License for more details. -# -# You should have received a copy of the GNU General Public License -# along with ash R toolbox. If not, see <https://www.gnu.org/licenses/>. -# /// -# -# -# plotting/datasheet.R -# -# - - -time_panel = function (list_df2plot, df_meta, trend_period, info_header, time_header, layout_matrix, info_ratio, time_ratio, var_ratio, outdirTmp) { - - - # Number of type/variable - nbp = length(list_df2plot) - - # Get all different stations code - Code = levels(factor(df_meta$code)) - nCode = length(Code) - - - df_trend = list_df2plot[[1]]$trend - - # Convert 'trend_period' to list - trend_period = as.list(trend_period) - # Number of trend period - nPeriod_trend = length(trend_period) - - nPeriod_max = 0 - for (code in Code) { - - df_trend_code = df_trend[df_trend$code == code,] - - Start = df_trend_code$period_start - UStart = levels(factor(Start)) - - End = df_trend_code$period_end - UEnd = levels(factor(End)) - - nPeriod = max(length(UStart), length(UEnd)) - - if (nPeriod > nPeriod_max) { - nPeriod_max = nPeriod - } - } - - - Start_code = vector(mode='list', length=nCode) - End_code = vector(mode='list', length=nCode) - Code_code = vector(mode='list', length=nCode) - Periods_code = vector(mode='list', length=nCode) - - for (j in 1:nCode) { - - code = Code[j] - - df_trend_code = df_trend[df_trend$code == code,] - - Start = df_trend_code$period_start - UStart = levels(factor(Start)) - - End = df_trend_code$period_end - UEnd = levels(factor(End)) - - nPeriod = max(length(UStart), length(UEnd)) - - Periods = c() - - for (i in 1:nPeriod_trend) { - Periods = append(Periods, - paste(Start[i], - End[i], - sep=' / ')) - } - - Start_code[[j]] = Start - End_code[[j]] = End - Code_code[[j]] = code - Periods_code[[j]] = Periods - - } - - TrendMean_code = array(rep(1, nPeriod_trend*nbp*nCode), - dim=c(nPeriod_trend, nbp, nCode)) - - for (j in 1:nPeriod_max) { - - for (k in 1:nCode) { - - code = Code[k] - - for (i in 1:nbp) { - - df_data = list_df2plot[[i]]$data - df_trend = list_df2plot[[i]]$trend - p_threshold = list_df2plot[[i]]$p_threshold - - df_data_code = df_data[df_data$code == code,] - df_trend_code = df_trend[df_trend$code == code,] - - Start = Start_code[Code_code == code][[1]][j] - End = End_code[Code_code == code][[1]][j] - Periods = Periods_code[Code_code == code][[1]][j] - - df_data_code_per = - df_data_code[df_data_code$Date >= Start - & df_data_code$Date <= End,] - - df_trend_code_per = - df_trend_code[df_trend_code$period_start == Start - & df_trend_code$period_end == End,] - - Ntrend = nrow(df_trend_code_per) - if (Ntrend > 1) { - df_trend_code_per = df_trend_code_per[1,] - } - - dataMean = mean(df_data_code_per$Qm3s, na.rm=TRUE) - trendMean = df_trend_code_per$trend / dataMean - - if (df_trend_code_per$p <= p_threshold){ - TrendMean_code[j, i, k] = trendMean - } else { - TrendMean_code[j, i, k] = NA - } - } - } - } - - minTrendMean = apply(TrendMean_code, c(1, 2), min, na.rm=TRUE) - maxTrendMean = apply(TrendMean_code, c(1, 2), max, na.rm=TRUE) - - for (code in Code) { - - # Print code of the station for the current plotting - print(paste("Datasheet for station :", code)) - - nbh = as.numeric(info_header) + as.numeric(!is.null(time_header)) - nbg = nbp + nbh - - P = vector(mode='list', length=nbg) - - if (info_header) { - time_header_code = time_header[time_header$code == code,] - - Hinfo = info_panel(list_df2plot, - df_meta, - df_shapefile=df_shapefile, - codeLight=code, - df_data_code=time_header_code) - P[[1]] = Hinfo - # P[[1]] = void - } - - if (!is.null(time_header)) { - - time_header_code = time_header[time_header$code == code,] - axis_xlim = c(min(time_header_code$Date), - max(time_header_code$Date)) - - Htime = time_panel_alone(time_header_code, df_trend_code=NULL, - trend_period=trend_period, missRect=TRUE, - unit2day=365.25, type='Q', grid=TRUE, first=FALSE) - - P[[2]] = Htime - } - - # map = map_panel() - - nbcol = ncol(as.matrix(layout_matrix)) - for (i in 1:nbp) { - df_data = list_df2plot[[i]]$data - df_trend = list_df2plot[[i]]$trend - p_threshold = list_df2plot[[i]]$p_threshold - unit2day = list_df2plot[[i]]$unit2day - missRect = list_df2plot[[i]]$missRect - type = list_df2plot[[i]]$type - - df_data_code = df_data[df_data$code == code,] - df_trend_code = df_trend[df_trend$code == code,] - - color = c() - # for (j in 1:nrow(df_trend_code)) { - grey = 85 - for (j in 1:nPeriod_max) { - if (df_trend_code$p[j] <= p_threshold){ - # color_res = get_color(df_trend_code$trend[j], - # minTrend[i], - # maxTrend[i], - # palette_name='perso', - # reverse=TRUE) - - Start = Start_code[Code_code == code][[1]][j] - End = End_code[Code_code == code][[1]][j] - Periods = Periods_code[Code_code == code][[1]][j] - - df_data_code_per = - df_data_code[df_data_code$Date >= Start - & df_data_code$Date <= End,] - - df_trend_code_per = - df_trend_code[df_trend_code$period_start == Start - & df_trend_code$period_end == End,] - - Ntrend = nrow(df_trend_code_per) - if (Ntrend > 1) { - df_trend_code_per = df_trend_code_per[1,] - } - - dataMean = mean(df_data_code$Qm3s, na.rm=TRUE) - trendMean = df_trend_code_per$trend / dataMean - - color_res = get_color(trendMean, - minTrendMean[j, i], - maxTrendMean[j, i], - palette_name='perso', - reverse=TRUE) - - colortmp = color_res - } else { - colortmp = paste('grey', grey, sep='') - grey = grey - 10 - } - - color = append(color, colortmp) - } - - p = time_panel_alone(df_data_code, df_trend_code, type=type, - p_threshold=p_threshold, missRect=missRect, - trend_period=trend_period, - mean_period=mean_period, axis_xlim=axis_xlim, - unit2day=unit2day, grid=FALSE, last=(i > nbp-nbcol), - color=color) - - P[[i+nbh]] = p - } - - layout_matrix = as.matrix(layout_matrix) - nel = nrow(layout_matrix)*ncol(layout_matrix) - - idNA = which(is.na(layout_matrix), arr.ind=TRUE) - - layout_matrix[idNA] = seq(max(layout_matrix, na.rm=TRUE) + 1, - max(layout_matrix, na.rm=TRUE) + 1 + - nel) - - layout_matrix_H = layout_matrix + nbh - - - info_ratio_scale = info_ratio - time_ratio_scale = time_ratio - var_ratio_scale = var_ratio - - ndec_info = 0 - ndec_time = 0 - ndec_var = 0 - - if (info_ratio_scale != round(info_ratio_scale)) { - ndec_info = nchar(gsub('^[0-9]+.', '', - as.character(info_ratio_scale))) - } - - if (time_ratio_scale != round(time_ratio_scale)) { - ndec_time = nchar(gsub('^[0-9]+.', '', - as.character(time_ratio_scale))) - } - - if (var_ratio_scale != round(var_ratio_scale)) { - ndec_var = nchar(gsub('^[0-9]+.', '', - as.character(var_ratio_scale))) - } - - ndec = max(c(ndec_info, ndec_time, ndec_var)) - - info_ratio_scale = info_ratio_scale * 10^ndec - time_ratio_scale = time_ratio_scale * 10^ndec - var_ratio_scale = var_ratio_scale * 10^ndec - - LM = c() - LMcol = ncol(layout_matrix_H) - LMrow = nrow(layout_matrix_H) - for (i in 1:(LMrow+nbh)) { - - if (info_header & i == 1) { - # LM = rbind(LM, rep(i, times=LMcol)) - LM = rbind(LM, - matrix(rep(rep(i, times=LMcol), - times=info_ratio_scale), - ncol=LMcol, byrow=TRUE)) - - } else if (!is.null(time_header) & i == 2) { - LM = rbind(LM, - matrix(rep(rep(i, times=LMcol), - times=time_ratio_scale), - ncol=LMcol, byrow=TRUE)) - - } else { - LM = rbind(LM, - matrix(rep(layout_matrix_H[i-nbh,], - times=var_ratio_scale), - ncol=LMcol, byrow=TRUE)) - }} - - plot = grid.arrange(grobs=P, layout_matrix=LM) - - # plot = grid.arrange(rbind(cbind(ggplotGrob(P[[2]]), ggplotGrob(P[[2]])), cbind(ggplotGrob(P[[3]]), ggplotGrob(P[[3]]))), heights=c(1/3, 2/3)) - - - # Saving - ggsave(plot=plot, - path=outdirTmp, - filename=paste(as.character(code), '.pdf', sep=''), - width=21, height=29.7, units='cm', dpi=100) - - } -} - - -time_panel_alone = function (df_data_code, df_trend_code, type, p_threshold=0.1, missRect=FALSE, unit2day=365.25, trend_period=NULL, mean_period=NULL, axis_xlim=NULL, grid=TRUE, last=FALSE, first=FALSE, color=NULL) { - - # If 'type' is square root apply it to data - if (type == 'sqrt(Q)') { - df_data_code$Qm3s = sqrt(df_data_code$Qm3s) - } - - # Compute max of flow - maxQ = max(df_data_code$Qm3s, na.rm=TRUE) - # Get the magnitude of the max of flow - power = get_power(maxQ) - - # Normalize the max flow by it's magnitude - maxQtmp = maxQ/10^power - # Compute the spacing between y ticks - if (maxQtmp >= 5) { - dbrk = 1.0 - } else if (maxQtmp < 5 & maxQtmp >= 3) { - dbrk = 0.5 - } else if (maxQtmp < 3 & maxQtmp >= 2) { - dbrk = 0.4 - } else if (maxQtmp < 2 & maxQtmp >= 1) { - dbrk = 0.2 - } else if (maxQtmp < 1) { - dbrk = 0.1 - } - - # Get the spacing in the right magnitude - dbrk = dbrk * 10^power - # Fix the accuracy for label - accuracy = NULL - - # Time span in the unit of time - dDate = as.numeric(df_data_code$Date[length(df_data_code$Date)] - - df_data_code$Date[1]) / unit2day - # Compute the spacing between x ticks - if (dDate >= 100) { - datebreak = 25 - dateminbreak = 5 - } else if (dDate < 100 & dDate >= 50) { - datebreak = 10 - dateminbreak = 1 - } else if (dDate < 50) { - datebreak = 5 - dateminbreak = 1 - } - - # Open new plot - p = ggplot() + theme_ash - - # If it is the lats plot of the pages or not - if (last) { - if (first) { - p = p + - theme(plot.margin=margin(5, 5, 5, 5, unit="mm")) - } else { - p = p + - theme(plot.margin=margin(0, 5, 5, 5, unit="mm")) - } - - # If it is the first plot of the pages or not - } else { - if (first) { - p = p + - theme(plot.margin=margin(5, 5, 0, 5, unit="mm")) - } else { - p = p + - theme(plot.margin=margin(0, 5, 0, 5, unit="mm")) - } - } - - ## Sub period background ## - if (!is.null(trend_period)) { - - # trend_period = as.list(trend_period) - # Imin = 10^99 - # for (per in trend_period) { - # I = interval(per[1], per[2]) - # if (I < Imin) { - # Imin = I - # trend_period_min = as.Date(per) - # } - # } - # p = p + - # geom_rect(aes(xmin=min(df_data_code$Date), - # ymin=0, - # xmax=trend_period_min[1], - # ymax= maxQ*1.1), - # linetype=0, fill='grey97') + - - # geom_rect(aes(xmin=trend_period_min[2], - # ymin=0, - # xmax=max(df_data_code$Date), - # ymax= maxQ*1.1), - # linetype=0, fill='grey97') - - # Convert trend period to list if it is not - trend_period = as.list(trend_period) - # Fix a disproportionate minimum for period - Imin = 10^99 - # For all the sub period of analysis in 'trend_period' - for (per in trend_period) { - # Compute time interval of period - I = interval(per[1], per[2]) - # If it is the smallest interval - if (I < Imin) { - # Store it - Imin = I - # Fix min period of analysis - trend_period_min = as.Date(per) - } - } - # Search for the index of the closest existing date - # to the start of the min period of analysis - idMinPer = which.min(abs(df_data_code$Date - trend_period_min[1])) - # Same for the end of the min period of analysis - idMaxPer = which.min(abs(df_data_code$Date - trend_period_min[2])) - # Get the start and end date associated - minPer = df_data_code$Date[idMinPer] - maxPer = df_data_code$Date[idMaxPer] - - # If it is not a flow or sqrt of flow time serie - if (type != 'sqrt(Q)' & type != 'Q') { - # If there is an 'axis_lim' - if (!is.null(axis_xlim)) { - # If the temporary start of period is smaller - # than the fix start of x axis limit - if (minPer < axis_xlim[1]) { - # Set the start of the period to the start of - # the x axis limit - minPer = axis_xlim[1] - } - } - } - # If it is not a flow or sqrt of flow time serie - if (type != 'sqrt(Q)' & type != 'Q') { - # If there is an 'axis_lim' - if (!is.null(axis_xlim)) { - # If the temporary end of period plus one year - # is smaller than the fix end of x axis limit - if (maxPer + years(1) < axis_xlim[2]) { - # Add one year the the temporary end of period - maxPer = maxPer + years(1) - } else { - # Set the start of the period to the start of - # the x axis limit - maxPer = axis_xlim[2] - } - # Add one year the the temporary end of period - # if there is no 'axis_lim' - } else { - maxPer = maxPer + years(1) - } - } - - # Draw rectangle to delimiting the sub period - p = p + - geom_rect(aes(xmin=minPer, - ymin=0, - xmax=maxPer, - ymax= maxQ*1.1), - linetype=0, fill='grey97') - } - - ## Mean step ## - # If there is a 'mean_period' - if (!is.null(mean_period)) { - # Convert 'mean_period' to list - mean_period = as.list(mean_period) - # Number of mean period - nPeriod_mean = length(mean_period) - - # Blank tibble to store variable in order to plot - # rectangle for mean period - plot_mean = tibble() - # Blank tibble to store variable in order to plot - # upper limit of rectangle for mean period - plot_line = tibble() - # For all mean period - for (j in 1:nPeriod_mean) { - # Get the current start and end of the sub period - Start_mean = mean_period[[j]][1] - End_mean = mean_period[[j]][2] - - # Extract the data corresponding to this sub period - df_data_code_per = - df_data_code[df_data_code$Date >= Start_mean - & df_data_code$Date <= End_mean,] - - # Min for the sub period - xmin = min(df_data_code_per$Date) - # If the min over the sub period is greater - # than the min of the entier period and - # it is not the first sub period - if (xmin > min(df_data_code$Date) & j != 1) { - # Substract 6 months to be in the middle of - # the previous year - xmin = xmin - months(6) - } - # If it is not a flow or sqrt of flow time serie and - # it is the first period - if (type != 'sqrt(Q)' & type != 'Q' & j == 1) { - # If there is an x axis limit - if (!is.null(axis_xlim)) { - # If the min of the period is before the x axis min - if (xmin < axis_xlim[1]) { - # The min for the sub period is the x axis - xmin = axis_xlim[1] - } - } - } - # Max for the sub period - xmax = max(df_data_code_per$Date) - # If the max over the sub period is smaller - # than the max of the entier period and - # it is not the last sub period - if (xmax < max(df_data_code$Date) & j != nPeriod_mean) { - # Add 6 months to be in the middle of - # the following year - xmax = xmax + months(6) - } - # If it is not a flow or sqrt of flow time serie and - # it is the last period - if (type != 'sqrt(Q)' & type != 'Q' & j == nPeriod_mean) { - # If there is an x axis limit - if (!is.null(axis_xlim)) { - # If the max of the period plus 1 year - # is smaller thant the max of the x axis limit - if (xmax + years(1) < axis_xlim[2]) { - # Add one year to the max to include - # the entire last year graphically - xmax = xmax + years(1) - } else { - # The max of this sub period is the max - # of the x axis limit - xmax = axis_xlim[2] - } - # If there is no axis limit - } else { - # Add one year to the max to include - # the entire last year graphically - xmax = xmax + years(1) - } - } - # Mean of the flow over the sub period - ymax = mean(df_data_code_per$Qm3s, na.rm=TRUE) - - # Create temporary tibble with variable - # to create rectangle for mean step - plot_meantmp = tibble(xmin=xmin, xmax=xmax, - ymin=0, ymax=ymax, period=j) - # Bind it to the main tibble to store it with other period - plot_mean = bind_rows(plot_mean, plot_meantmp) - - # Create vector for the upper limit of the rectangle - abs = c(xmin, xmax) - ord = c(ymax, ymax) - - # Create temporary tibble with variable - # to create upper limit for rectangle - plot_linetmp = tibble(abs=abs, ord=ord, period=j) - # Bind it to the main tibble to store it with other period - plot_line = bind_rows(plot_line, plot_linetmp) - } - # Plot rectangles - p = p + - geom_rect(data=plot_mean, - aes(xmin=xmin, ymin=ymin, - xmax=xmax, ymax=ymax), - linetype=0, fill='grey93') - # Plot upper line for rectangle - p = p + - geom_line(data=plot_line, - aes(x=abs, y=ord, group=period), - color='grey85', - size=0.15) - - # for all the sub periods except the last one - for (i in 1:(nPeriod_mean-1)) { - # The y limit of rectangle is the max of - # the two neighboring mean step rectangle - yLim = max(c(plot_mean$ymax[i], plot_mean$ymax[i+1])) - # The x limit is the x max of the ith rectangle - xLim = plot_mean$xmax[i] - # Make a tibble to store data - plot_lim = tibble(x=c(xLim, xLim), y=c(0, yLim)) - # Plot the limit of rectangles - p = p + - geom_line(data=plot_lim, aes(x=x, y=y), - linetype='dashed', size=0.15, color='grey85') - } - - } - - ### Grid ### - if (grid) { - # If there is no axis limit - if (is.null(axis_xlim)) { - # The min and the max is set by - # the min and the max of the date data - xmin = min(df_data_code$Date) - xmax = max(df_data_code$Date) - } else { - # Min and max is set with the limit axis parameter - xmin = axis_xlim[1] - xmax = axis_xlim[2] - } - # Create a vector for all the y grid position - ygrid = seq(0, maxQ*10, dbrk) - # Blank vector to store position - ord = c() - abs = c() - # For all the grid element - for (i in 1:length(ygrid)) { - # Store grid position - ord = c(ord, rep(ygrid[i], times=2)) - abs = c(abs, xmin, xmax) - } - # Create a tibble to store all the position - plot_grid = tibble(abs=as.Date(abs), ord=ord) - # Plot the y grid - p = p + - geom_line(data=plot_grid, - aes(x=abs, y=ord, group=ord), - color='grey85', - size=0.15) - } - - ### Data ### - # If it is a square root flow or flow - if (type == 'sqrt(Q)' | type == 'Q') { - # Plot the data as line - p = p + - geom_line(aes(x=df_data_code$Date, y=df_data_code$Qm3s), - color='grey20', - size=0.3, - lineend="round") - } else { - # Plot the data as point - p = p + - geom_point(aes(x=df_data_code$Date, y=df_data_code$Qm3s), - shape=21, color='grey50', fill='grey97', size=1) - } - - ### Missing data ### - # If the option is TRUE - if (missRect) { - # Remove NA data - NAdate = df_data_code$Date[is.na(df_data_code$Qm3s)] - # Get the difference between each point of date data without NA - dNAdate = diff(NAdate) - # If difference of day is not 1 then - # it is TRUE for the beginning of each missing data period - NAdate_Down = NAdate[append(Inf, dNAdate) != 1] - # If difference of day is not 1 then - # it is TRUE for the ending of each missing data period - NAdate_Up = NAdate[append(dNAdate, Inf) != 1] - - # Plot the missing data period - p = p + - geom_rect(aes(xmin=NAdate_Down, - ymin=0, - xmax=NAdate_Up, - ymax=maxQ*1.1), - linetype=0, fill='Wheat', alpha=0.4) - } - - ### Trend ### - # If there is trends - if (!is.null(df_trend_code)) { - - # Extract starting period of trends - Start = df_trend_code$period_start - # Get the name of the different period - UStart = levels(factor(Start)) - # Same for ending - End = df_trend_code$period_end - UEnd = levels(factor(End)) - - # Compute the max of different start and end - # so the number of different period - nPeriod_trend = max(length(UStart), length(UEnd)) - - # Blank tibble to store trend data and legend data - plot_trend = tibble() - leg_trend = tibble() - # For all the different period - for (i in 1:nPeriod_trend) { - - # Get the trend associated to the first period - df_trend_code_per = - df_trend_code[df_trend_code$period_start == Start[i] - & df_trend_code$period_end == End[i],] - - # Number of trend selected - Ntrend = nrow(df_trend_code_per) - # If the number of trend is greater than a unique one - if (Ntrend > 1) { - # Extract only the first hence it is the same period - df_trend_code_per = df_trend_code_per[1,] - } - - # Search for the index of the closest existing date - # to the start of the trend period of analysis - iStart = which.min(abs(df_data_code$Date - Start[i])) - # Same for the end - iEnd = which.min(abs(df_data_code$Date - End[i])) - - # Get the start and end date associated - xmin = df_data_code$Date[iStart] - xmax = df_data_code$Date[iEnd] - - # If there is a x axis limit - if (!is.null(axis_xlim)) { - # If the min of the current period - # is smaller than the min of the x axis limit - if (xmin < axis_xlim[1]) { - # The min of the period is the min - # of the x axis limit - xmin = axis_xlim[1] - } - # Same for end - if (xmax > axis_xlim[2]) { - xmax = axis_xlim[2] - } - } - - # Create vector to store x data - abs = c(xmin, xmax) - # Convert the number of day to the unit of the period - abs_num = as.numeric(abs) / unit2day - # Compute the y of the trend - ord = abs_num * df_trend_code_per$trend + - df_trend_code_per$intercept - - # Create temporary tibble with variable to plot trend - # for each period - plot_trendtmp = tibble(abs=abs, ord=ord, period=i) - # Bind it to the main tibble to store it with other period - plot_trend = bind_rows(plot_trend, plot_trendtmp) - - # If there is a x axis limit - if (!is.null(axis_xlim)) { - # The x axis limit is selected - codeDate = axis_xlim - } else { - # The entire date data is selected - codeDate = df_data_code$Date - } - # The flow data is extract - codeQ = df_data_code$Qm3s - - # Position of the x beginning and end of the legend symbol - x = gpct(2, codeDate, shift=TRUE) - xend = x + gpct(3, codeDate) - - # Position of the y beginning and end of the legend symbol - dy = gpct(7, codeQ, ref=0) - y = gpct(100, codeQ, ref=0) - (i-1)*dy - yend = y - - # Position of x for the beginning of the associated text - xt = xend + gpct(1, codeDate) - - # Position of the background rectangle of the legend - xminR = x - gpct(1, codeDate) - yminR = y - gpct(4, codeQ, ref=0) - xmaxR = x + gpct(24, codeDate) - ymaxR = y + gpct(5, codeQ, ref=0) - - # Get the tendance analyse - trend = df_trend_code_per$trend - # Compute the magnitude of the trend - power = get_power(trend) - # Convert it to character - powerC = as.character(power) - # Get the power of ten of magnitude - brk = 10^power - # Convert trend to character for sientific expression - trendC = as.character(round(trend / brk, 2)) - - # Create temporary tibble with variable to plot legend - leg_trendtmp = tibble(x=x, xend=xend, - y=y, yend=yend, - xt=xt, - trendC=trendC, - powerC=powerC, - xminR=xminR, yminR=yminR, - xmaxR=xmaxR, ymaxR=ymaxR, - period=i) - # Bind it to the main tibble to store it with other period - leg_trend = bind_rows(leg_trend, leg_trendtmp) - } - - # For all periods - for (i in 1:nPeriod_trend) { - # Extract the trend of the current sub period - leg_trend_per = leg_trend[leg_trend$period == i,] - - # Plot the background for legend - p = p + - geom_rect(data=leg_trend_per, - aes(xmin=xminR, - ymin=yminR, - xmax=xmaxR, - ymax=ymaxR), - linetype=0, fill='white', alpha=0.5) - } - - # For all periods - for (i in 1:nPeriod_trend) { - # Extract the trend of the current sub period - leg_trend_per = leg_trend[leg_trend$period == i,] - - # Get the character variable for naming the trend - trendC = leg_trend_per$trendC - powerC = leg_trend_per$powerC - - # Create the name of the trend - label = bquote(bold(.(trendC)~'x'~'10'^{.(powerC)})~'['*m^{3}*'.'*s^{-1}*'.'*an^{-1}*']') - - # Plot the trend symbole and value of the legend - p = p + - annotate("segment", - x=leg_trend_per$x, xend=leg_trend_per$xend, - y=leg_trend_per$y, yend=leg_trend_per$yend, - color=color[i], - linetype='solid', - lwd=1) + - - annotate("text", - label=label, size=3, - x=leg_trend_per$xt, y=leg_trend_per$y, - hjust=0, vjust=0.4, - color=color[i]) - } - - # For all periods - for (i in 1:nPeriod_trend) { - # Extract the trend of the current sub period - plot_trend_per = plot_trend[plot_trend$period == i,] - - # Plot the line of white background of each trend - p = p + - geom_line(data=plot_trend_per, - aes(x=abs, y=ord), - color='white', - linetype='solid', - size=1.5, - lineend="round") - } - - # For all periods - for (i in 1:nPeriod_trend) { - # Extract the trend of the current sub period - plot_trend_per = plot_trend[plot_trend$period == i,] - - # Plot the line of trend - p = p + - geom_line(data=plot_trend_per, - aes(x=abs, y=ord), - color=color[i], - linetype='solid', - size=0.75, - lineend="round") - } - } - - # Title - p = p + - ggtitle(bquote(bold(.(type))~~'['*m^{3}*'.'*s^{-1}*']')) - - # If the is no x axis limit - if (is.null(axis_xlim)) { - # Parameters of the x axis contain the limit of the date data - p = p + - scale_x_date(date_breaks=paste( - as.character(datebreak), - 'year', sep=' '), - date_minor_breaks=paste( - as.character(dateminbreak), - 'year', sep=' '), - guide='axis_minor', - date_labels="%Y", - limits=c(min(df_data_code$Date), - max(df_data_code$Date)), - expand=c(0, 0)) - } else { - # Parameters of the x axis contain the x axis limit - p = p + - scale_x_date(date_breaks=paste( - as.character(datebreak), - 'year', sep=' '), - date_minor_breaks=paste( - as.character(dateminbreak), - 'year', sep=' '), - guide='axis_minor', - date_labels="%Y", - limits=axis_xlim, - expand=c(0, 0)) - } - - # Parameters of the y axis - p = p + - scale_y_continuous(breaks=seq(0, maxQ*10, dbrk), - limits=c(0, maxQ*1.1), - expand=c(0, 0), - labels=label_number(accuracy=accuracy)) - - return(p) -} - - -- GitLab