panel.R 67.22 KiB
# Usefull library
library(ggplot2)
library(scales)
library(qpdf)
library(gridExtra)
library(gridtext)
library(dplyr)
library(grid)
library(ggh4x)
library(RColorBrewer)
library(rgdal)
display_type = function (type, bold=FALSE) {
    if (type == "QA") {
        if (bold) {
            disp = bquote(Q[A])
        } else {
            disp = bquote(bold(Q[A]))
    } else if (type == "QMNA") {
        if (bold) {
            disp = bquote(Q[MNA])
        } else {
            disp = bquote(bold(Q[MNA]))
    } else if (type == "VCN10") {
        if (bold) {
            disp = bquote(V[CN10])
        } else {
            disp = bquote(bold(V[CN10]))
    return (disp)
# Personal theme
theme_ash =
    theme(
        # White background
        panel.background=element_rect(fill='white'),
        # Font
        text=element_text(family='sans'),
        # Border of plot
        panel.border = element_rect(color="grey85",
                                    fill=NA,
                                    size=0.7),
        # Grid
        panel.grid.major.x=element_blank(),
        panel.grid.major.y=element_blank(),
        # 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'),
        # Ticks length
        axis.ticks.length=unit(1.5, 'mm'),
        # Ticks minor
        ggh4x.axis.ticks.length.minor=rel(0.5),
        # Title
        plot.title=element_text(size=9, vjust=-2, 
                                hjust=-1E-3, color='grey20'), 
        # Axis title
        axis.title.x=element_blank(),
7172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140
axis.title.y=element_blank(), # Axis line axis.line.x=element_blank(), axis.line.y=element_blank(), ) time_panel = 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, 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) {
141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210
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] } }
211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280
} # 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
281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350
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,
351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420
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 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 ###
421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490
# 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
491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560
# 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,
561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630
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, lineend="round") } # For all periods for (i in 1:nPeriod_trend) { # Extract the trend of the current sub period
631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700
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.5, 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) } matrice_panel = function (list_df2plot, df_meta, trend_period, slice=NULL, outdirTmp='', outnameTmp='matrix') { nbp = length(list_df2plot) # Get all different stations code Code = levels(factor(df_meta$code)) nCode = length(Code) if (!is.null(slice)) { # By nMat = as.integer(nCode/slice) + 1
701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770
sublist_df2plot = list_df2plot for (i in 1:nMat) { subdf_meta = df_meta[(slice*(i-1)+1):(slice*i),] subdf_meta = subdf_meta[!is.na(subdf_meta$code),] subCode = subdf_meta$code for (j in 1:nbp) { df_datatmp = list_df2plot[[j]]$data df_trendtmp = list_df2plot[[j]]$trend subdf_data = df_datatmp[(df_datatmp$code %in% subCode),] subdf_trend = df_trendtmp[(df_trendtmp$code %in% subCode),] sublist_df2plot[[j]]$data = subdf_data sublist_df2plot[[j]]$trend = subdf_trend } mat = matrice_panel(sublist_df2plot, subdf_meta, trend_period=trend_period, outdirTmp=outdirTmp, outnameTmp=paste(outnameTmp, '_', i, sep='')) } } else { print(paste('matrix :', outnameTmp)) # nbp = length(list_df2plot) # # Get all different stations code # Code = levels(factor(df_meta$code)) # nCode = length(Code) df_trend = list_df2plot[[1]]$trend 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
771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840
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_max) { Periods = append(Periods, paste(substr(Start[i], 1, 4), substr(End[i], 1, 4), sep=' / ')) } Start_code[[j]] = Start End_code[[j]] = End Code_code[[j]] = code Periods_code[[j]] = Periods } TrendMean_code = array(rep(1, nPeriod_max*nbp*nCode), dim=c(nPeriod_max, 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 TrendMean_code[j, i, k] = trendMean } } } minTrendMean = apply(TrendMean_code, c(1, 2), min, na.rm=TRUE) maxTrendMean = apply(TrendMean_code, c(1, 2), max, na.rm=TRUE) Periods_mat = c()
841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910
NPeriod_mat = c() Type_mat = list() Code_mat = c() Pthresold_mat = c() TrendMean_mat = c() DataMean_mat = c() Fill_mat = c() Color_mat = c() for (j in 1:nPeriod_max) { for (code in Code) { for (i in 1:nbp) { df_data = list_df2plot[[i]]$data df_trend = list_df2plot[[i]]$trend p_threshold = list_df2plot[[i]]$p_threshold type = list_df2plot[[i]]$type 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){ color_res = get_color(trendMean, minTrendMean[j, i], maxTrendMean[j, i], palette_name='perso', reverse=TRUE) fill = color_res$color color = 'white' Pthresold = p_thresold } else { fill = 'white' color = 'grey85' Pthresold = NA } Periods_mat = append(Periods_mat, Periods) NPeriod_mat = append(NPeriod_mat, j) Type_mat = append(Type_mat, type) Code_mat = append(Code_mat, code) Pthresold_mat = append(Pthresold_mat, Pthresold) TrendMean_mat = append(TrendMean_mat, trendMean) DataMean_mat = append(DataMean_mat, dataMean) Fill_mat = append(Fill_mat, fill) Color_mat = append(Color_mat, color) } } }
911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980
height = length(Code) width = nbp * 2 * nPeriod_max + nPeriod_max options(repr.plot.width=width, repr.plot.height=height) mat = ggplot() + theme_ash + theme( panel.border=element_blank(), axis.text.x=element_blank(), axis.text.y=element_blank(), axis.ticks.y=element_blank(), axis.ticks.x=element_blank(), axis.title.y=element_blank(), ) # xt = -1 # yt = height + 1.75 # Title = bquote(bold(Territoire)) # mat = mat + # annotate("text", x=xt, y=yt, # label=Title, # hjust=0, vjust=0.5, # size=6, color="#00A3A8") for (j in 1:nPeriod_max) { Type_mat_per = Type_mat[NPeriod_mat == j] Code_mat_per = Code_mat[NPeriod_mat == j] Pthresold_mat_per = Pthresold_mat[NPeriod_mat == j] TrendMean_mat_per = TrendMean_mat[NPeriod_mat == j] DataMean_mat_per = DataMean_mat[NPeriod_mat == j] Fill_mat_per = Fill_mat[NPeriod_mat == j] Color_mat_per = Color_mat[NPeriod_mat == j] Xtmp = as.integer(factor(as.character(Type_mat_per))) Xc = j + (j - 1)*nbp*2 Xm = Xtmp + (j - 1)*nbp*2 + j X = Xtmp + (j - 1)*nbp*2 + nbp + j Y = as.integer(factor(Code_mat_per)) x = Xc - 0.4 xend = X[length(X)] + 0.25 y = height + 1 yend = height + 1 mat = mat + annotate("segment", x=x, xend=xend, y=y, yend=yend, color="grey40", size=0.35) yt = y + 0.15 Start = trend_period[[j]][1] End = trend_period[[j]][2] periodName = bquote(bold('Priode')~bold(.(as.character(j)))) # bquote(bold(.(Start))~'/'~bold(.(End))) mat = mat + annotate("text", x=x, y=yt, label=periodName, hjust=0, vjust=0.5, size=3, color='grey40') for (i in 1:length(X)) { mat = mat +
981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050
gg_circle(r=0.45, xc=X[i], yc=Y[i], fill=Fill_mat_per[i], color=Color_mat_per[i]) + gg_circle(r=0.45, xc=Xm[i], yc=Y[i], fill='white', color='grey40') + gg_circle(r=0.45, xc=Xc, yc=Y[i], fill='white', color='grey40') } for (i in 1:length(TrendMean_mat_per)) { trendMean = TrendMean_mat_per[i] trendC = round(trendMean*100, 2) if (!is.na(Pthresold_mat_per[i])) { Tcolor = 'white' } else { Tcolor = 'grey85' } dataMean = round(DataMean_mat_per[i], 2) mat = mat + annotate('text', x=X[i], y=Y[i], label=trendC, hjust=0.5, vjust=0.5, size=3, color=Tcolor) + annotate('text', x=Xm[i], y=Y[i], label=dataMean, hjust=0.5, vjust=0.5, size=3, color='grey40') } mat = mat + annotate('text', x=Xc, y=max(Y) + 0.85, label=bquote(bold('Dbut')), hjust=0.5, vjust=0.5, size=3, color='grey20') + annotate('text', x=Xc, y=max(Y) + 0.6, label=bquote(bold('Fin')), hjust=0.5, vjust=0.5, size=3, color='grey20') for (i in 1:nbp) { type = list_df2plot[[i]]$type mat = mat + annotate('text', x=X[i], y=max(Y) + 0.7, label=bquote(.(type)), hjust=0.5, vjust=0.5, size=3.5, color='grey20') + annotate('text', x=Xm[i], y=max(Y) + 0.7, label=bquote(''*.(type)), hjust=0.5, vjust=0.5, size=3.5, color='grey20') } } for (i in 1:length(Code)) { code = Code[i] name = df_meta[df_meta$code == code,]$nom ncharMax = 30 if (nchar(name) > ncharMax) { name = paste(substr(name, 1, ncharMax), '...', sep='') }
1051105210531054105510561057105810591060106110621063106410651066106710681069107010711072107310741075107610771078107910801081108210831084108510861087108810891090109110921093109410951096109710981099110011011102110311041105110611071108110911101111111211131114111511161117111811191120
mat = mat + annotate('text', x=0.3, y=i + 0.14, label=bquote(bold(.(code))), hjust=1, vjust=0.5, size=3.5, color="#00A3A8") + annotate('text', x=0.3, y=i - 0.14, label=name, hjust=1, vjust=0.5, size=3.5, color="#00A3A8") for (j in 1:nPeriod_max) { Xc = j + (j - 1)*nbp*2 label = Periods_code[Code_code == code][[1]][j] periodStart = substr(label, 1, 4) periodEnd = substr(label, 8, 11) mat = mat + annotate('text', x=Xc, y=i + 0.13, label=bquote(bold(.(periodStart))), hjust=0.5, vjust=0.5, size=3, color='grey40') + annotate('text', x=Xc, y=i - 0.13, label=bquote(bold(.(periodEnd))), hjust=0.5, vjust=0.5, size=3, color='grey40') } } mat = mat + coord_fixed() + scale_x_continuous(limits=c(1 - rel(6), width + rel(0.5)), expand=c(0, 0)) + scale_y_continuous(limits=c(1 - rel(0.5), height + rel(1.5)), expand=c(0, 0)) # Saving matrix plot ggsave(plot=mat, path=outdirTmp, filename=paste(outnameTmp, '.pdf', sep=''), width=29.7, height=21, units='cm', dpi=100) } } map_panel = function (list_df2plot, df_meta, computer_data_path, fr_shpdir, fr_shpname, bs_shpdir, bs_shpname, rv_shpdir, rv_shpname, idPer=1, outdirTmp='', codeLight=NULL, margin=NULL, showSea=TRUE) { fr_shppath = file.path(computer_data_path, fr_shpdir, fr_shpname) rv_shppath = file.path(computer_data_path, rv_shpdir, rv_shpname) bs_shppath = file.path(computer_data_path, bs_shpdir, bs_shpname) # France fr_spdf = readOGR(dsn=fr_shppath, verbose=FALSE) proj4string(fr_spdf) = CRS("+proj=longlat +ellps=WGS84") # Trasnformation en Lambert93 france = spTransform(fr_spdf, CRS("+init=epsg:2154")) df_france = tibble(fortify(france)) # Bassin hydrographique bassin = readOGR(dsn=bs_shppath, verbose=FALSE) df_bassin = tibble(fortify(bassin))
1121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181118211831184118511861187118811891190
# Rseau hydrographique # river = readOGR(dsn=rv_shppath, verbose=FALSE) ### trop long ### # river = river[which(river$Classe == 1),] # df_river = tibble(fortify(river)) nbp = length(list_df2plot) # Get all different stations code Code = levels(factor(df_meta$code)) nCode = length(Code) df_trend = list_df2plot[[1]]$trend 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_max) { Periods = append(Periods, paste(substr(Start[i], 1, 4), substr(End[i], 1, 4), sep=' / ')) } Start_code[[j]] = Start End_code[[j]] = End Code_code[[j]] = code Periods_code[[j]] = Periods } TrendMean_code = array(rep(1, nPeriod_max*nbp*nCode), dim=c(nPeriod_max, nbp, nCode)) for (j in 1:nPeriod_max) {
1191119211931194119511961197119811991200120112021203120412051206120712081209121012111212121312141215121612171218121912201221122212231224122512261227122812291230123112321233123412351236123712381239124012411242124312441245124612471248124912501251125212531254125512561257125812591260
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 TrendMean_code[j, i, k] = trendMean } } } minTrendMean = apply(TrendMean_code, c(1, 2), min, na.rm=TRUE) maxTrendMean = apply(TrendMean_code, c(1, 2), max, na.rm=TRUE) ncolor = 256 nbTick = 10 for (i in 1:nbp) { if (i > 1 & !is.null(codeLight)) { break } outname = paste('map_', i, sep='') print(paste('map :', outname)) type = list_df2plot[[i]]$type map = ggplot() + theme_void() + # theme(plot.background=element_rect(fill=NA, # color="#EC4899")) + coord_fixed() + geom_polygon(data=df_france, aes(x=long, y=lat, group=group), color=NA, fill="grey97") + # geom_path(data=df_river, # aes(x=long, y=lat, group=group), # color="grey85", size=0.3) +
1261126212631264126512661267126812691270127112721273127412751276127712781279128012811282128312841285128612871288128912901291129212931294129512961297129812991300130113021303130413051306130713081309131013111312131313141315131613171318131913201321132213231324132513261327132813291330
geom_polygon(data=df_bassin, aes(x=long, y=lat, group=group), color="grey70", fill=NA, size=0.1) + geom_polygon(data=df_france, aes(x=long, y=lat, group=group), color="grey40", fill=NA, size=0.2) if (showSea) { xlim = c(280000, 790000) ylim = c(6110000, 6600000) } else { xlim = c(305000, 790000) ylim = c(6135000, 6600000) } if (is.null(codeLight)) { xmin = gpct(7, xlim, shift=TRUE) xint = c(0, 10*1E3, 50*1E3, 100*1E3) ymin = gpct(5, ylim, shift=TRUE) ymax = ymin + gpct(1, ylim) map = map + geom_line(aes(x=c(xmin, max(xint)+xmin), y=c(ymin, ymin)), color="grey40", size=0.2) + annotate("text", x=max(xint)+xmin+gpct(1, xlim), y=ymin, vjust=0, hjust=0, label="km", color="grey40", size=3) for (x in xint) { map = map + annotate("segment", x=x+xmin, xend=x+xmin, y=ymin, yend=ymax, color="grey40", size=0.2) + annotate("text", x=x+xmin, y=ymax+gpct(0.5, ylim), vjust=0, hjust=0.5, label=x/1E3, color="grey40", size=3) } } map = map + coord_sf(xlim=xlim, ylim=ylim, expand=FALSE) if (is.null(margin)) { map = map + theme(plot.margin=margin(t=0, r=0, b=0, l=0, unit="mm")) } else { map = map + theme(plot.margin=margin) } lon = c() lat = c() fill = c() shape = c() for (code in Code) { 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,]
1331133213331334133513361337133813391340134113421343134413451346134713481349135013511352135313541355135613571358135913601361136213631364136513661367136813691370137113721373137413751376137713781379138013811382138313841385138613871388138913901391139213931394139513961397139813991400
Start = Start_code[Code_code == code][[1]][idPer] End = End_code[Code_code == code][[1]][idPer] 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 color_res = get_color(trendMean, minTrendMean[idPer, i], maxTrendMean[idPer, i], palette_name='perso', reverse=TRUE, ncolor=ncolor, nbTick=nbTick) if (df_trend_code_per$p <= p_threshold){ filltmp = color_res$color palette = color_res$palette if (trendMean >= 0) { shapetmp = 24 } else { shapetmp = 25 } } else { filltmp = color_res$color palette = color_res$palette shapetmp = 21 } lontmp = df_meta[df_meta$code == code,]$L93X_m_BH lattmp = df_meta[df_meta$code == code,]$L93Y_m_BH lon = c(lon, lontmp) lat = c(lat, lattmp) fill = c(fill, filltmp) shape = c(shape, shapetmp) } plot_map = tibble(lon=lon, lat=lat, fill=fill, shape=shape, code=Code) if (is.null(codeLight)) { map = map + geom_point(data=plot_map, aes(x=lon, y=lat), shape=shape, size=5, stroke=1, color='grey50', fill=fill) } else { plot_map_codeNo = plot_map[plot_map$code != codeLight,] plot_map_code = plot_map[plot_map$code == codeLight,] map = map +
1401140214031404140514061407140814091410141114121413141414151416141714181419142014211422142314241425142614271428142914301431143214331434143514361437143814391440144114421443144414451446144714481449145014511452145314541455145614571458145914601461146214631464146514661467146814691470
geom_point(data=plot_map_codeNo, aes(x=lon, y=lat), shape=21, size=0.5, stroke=0.5, color='grey70', fill='grey97') + geom_point(data=plot_map_code, aes(x=lon, y=lat), shape=21, size=1.5, stroke=0.5, color='grey70', fill='grey70') } idTick = color_res$idTick labTick = color_res$labTick colTick = color_res$colTick nbTickmod = length(idTick) valNorm = nbTickmod * 10 ytick = idTick / max(idTick) * valNorm labTick = as.character(round(labTick*100, 2)) xtick = rep(0, times=nbTickmod) plot_palette = tibble(xtick=xtick, ytick=ytick, colTick=colTick, labTick=labTick) title = ggplot() + theme_void() + annotate('text', x=-0.3, y=0.15, label=bquote(bold(.(type))), hjust=0, vjust=0, size=10, color="#00A3A8") + geom_line(aes(x=c(-0.3, 3.3), y=c(0.05, 0.05)), size=0.6, color="#00A3A8") + scale_x_continuous(limits=c(-1, 1 + 3), expand=c(0, 0)) + scale_y_continuous(limits=c(0, 10), expand=c(0, 0)) + theme(plot.margin=margin(t=5, r=5, b=0, l=0, unit="mm")) pal = ggplot() + theme_void() + geom_point(data=plot_palette, aes(x=xtick, y=ytick), shape=21, size=5, stroke=1, color='white', fill=colTick) pal = pal + annotate('text', x=-0.3, y= valNorm + 23, label="Tendance", hjust=0, vjust=0.5, size=6, color='grey40') + annotate('text', x=-0.2, y= valNorm + 13, label=bquote(bold("% par an")), hjust=0, vjust=0.5, size=4, color='grey40') for (j in 1:nbTickmod) {
1471147214731474147514761477147814791480148114821483148414851486148714881489149014911492149314941495149614971498149915001501150215031504150515061507150815091510151115121513151415151516151715181519152015211522152315241525152615271528152915301531153215331534153515361537153815391540
pal = pal + annotate('text', x=xtick[j]+0.3, y=ytick[j], label=bquote(bold(.(labTick[j]))), hjust=0, vjust=0.7, size=3, color='grey40') } pal = pal + geom_point(aes(x=0, y=-20), shape=24, size=4, stroke=1, color='grey50', fill='grey97') + annotate('text', x=0.3, y=-20, label=bquote(bold("Hausse significative 10%")), hjust=0, vjust=0.5, size=3, color='grey40') pal = pal + geom_point(aes(x=0, y=-29), shape=21, size=4, stroke=1, color='grey50', fill='grey97') + annotate('text', x=0.3, y=-29, label=bquote(bold("Non significatif 10%")), hjust=0, vjust=0.7, size=3, color='grey40') pal = pal + geom_point(aes(x=0, y=-40), shape=25, size=4, stroke=1, color='grey50', fill='grey97') + annotate('text', x=0.3, y=-40, label=bquote(bold("Baisse significative 10%")), hjust=0, vjust=0.5, size=3, color='grey40') pal = pal + scale_x_continuous(limits=c(-1, 1 + 3), expand=c(0, 0)) + scale_y_continuous(limits=c(-60, valNorm + 35), expand=c(0, 0)) + theme(plot.margin=margin(t=0, r=5, b=5, l=0, unit="mm")) Map = list(map, title, pal) plot = grid.arrange(grobs=Map, layout_matrix= matrix(c(1, 1, 1, 2, 1, 1, 1, 3), nrow=2, byrow=TRUE)) if (is.null(codeLight)) { # Saving matrix plot ggsave(plot=plot, path=outdirTmp, filename=paste(outname, '.pdf', sep=''), width=29.7, height=21, units='cm', dpi=100) }
1541154215431544154515461547154815491550155115521553155415551556155715581559156015611562156315641565156615671568156915701571157215731574157515761577157815791580158115821583158415851586158715881589159015911592159315941595159615971598159916001601160216031604160516061607160816091610
} return (map) } 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("jan", "fv", "mar", "avr", "mai", "jun", "jul", "ao", "sep", "oct", "nov", "dc") 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"), angle=90, vjust=0.5, hjust=1), axis.ticks.x=element_blank(), axis.line.y=element_line(color='grey85', size=0.3)) 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) } info_panel = function(list_df2plot, df_meta, computer_data_path, fr_shpdir, fr_shpname, bs_shpdir, bs_shpname, rv_shpdir, rv_shpname, 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 } # yearLast = format(databin$Date[nrow(databin)], "%Y") # yearFirst = format(databin$Date[1], "%Y") # Nyear = yearLast - yearFirst + 1 map = map_panel(list_df2plot,
1611161216131614161516161617161816191620162116221623162416251626162716281629163016311632163316341635163616371638163916401641164216431644164516461647164816491650165116521653165416551656165716581659166016611662166316641665166616671668166916701671167216731674167516761677167816791680
df_meta, computer_data_path=computer_data_path, fr_shpdir=fr_shpdir, fr_shpname=fr_shpname, bs_shpdir=bs_shpdir, bs_shpname=bs_shpname, rv_shpdir=rv_shpdir, rv_shpname=rv_shpname, codeLight=codeLight, margin=margin(t=5, r=2, b=0, l=0, unit="mm"), showSea=FALSE) df_meta_code = df_meta[df_meta$code == codeLight,] nom = df_meta_code$nom nom = gsub("-", "-&nbsp;", nom) text1 = paste( "<b>", codeLight, '</b> - ', nom, sep='') text2 = paste( "<b>", "Gestionnaire : ", df_meta_code$gestionnaire, "<br>", "Rgion hydro : ", df_meta_code$region_hydro, "</b>", sep='') text3 = paste( "<b>", "Superficie : ", df_meta_code$surface_km2_BH, " [km<sup>2</sup>] <br>", "X = ", df_meta_code$L93X_m_BH, " [m ; Lambert 93]", "</b>", sep='') text4 = paste( "<b>", "Altitude : ", df_meta_code$altitude_m_BH, " [m]<br>", "Y = ", df_meta_code$L93Y_m_BH, " [m ; Lambert 93]", "</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=0.55, 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"))
1681168216831684168516861687168816891690169116921693169416951696169716981699170017011702170317041705170617071708170917101711171217131714171517161717171817191720172117221723172417251726172717281729173017311732173317341735173617371738173917401741174217431744174517461747174817491750
P = list(gtext1, gtext2, gtext3, gtext4, hyd, map) # P = list(void, void, void, void, void, void) plot = grid.arrange(grobs=P, layout_matrix=matrix(c(1, 1, 1, 6, 2, 2, 5, 6, 2, 2, 5, 6, 3, 4, 5, 6, 3, 4, 5, 6), nrow=5, byrow=TRUE)) return(plot) } 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='') {
1751175217531754175517561757175817591760176117621763176417651766176717681769177017711772177317741775177617771778177917801781178217831784178517861787178817891790179117921793179417951796179717981799180018011802180318041805180618071808180918101811181218131814181518161817181818191820
# 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)
1821182218231824182518261827182818291830183118321833183418351836183718381839184018411842184318441845184618471848184918501851185218531854185518561857185818591860186118621863186418651866186718681869187018711872187318741875187618771878187918801881188218831884188518861887188818891890
} get_color = function (value, min, max, ncolor=256, palette_name='perso', reverse=FALSE, nbTick=10) { if (palette_name == 'perso') { colorList = c('#0f3b57', '#1d7881', '#80c4a9', '#e2dac6', #mid '#fadfad', '#d08363', '#7e392f') } else { colorList = brewer.pal(11, palette_name) } nSample = length(colorList) palette = colorRampPalette(colorList)(ncolor) 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) if (reverse) { palette = rev(palette) palette_hot = rev(palette_hot) palette_cold = rev(palette_cold) } if (value < 0) { idNorm = (value - min) / (0 - min) id = round(idNorm*(ncolor - 1) + 1, 0) color = palette_cold[id] } else { idNorm = (value - 0) / (max - 0) id = round(idNorm*(ncolor - 1) + 1, 0) color = palette_hot[id] } if (min < 0 & max < 0) { paletteShow = palette_cold idTick = c() for (i in 1:nbTick) { id = round((ncolor-1)/(nbTick-1)*(i-1)) + 1 idTick = c(idTick, id) } labTick = seq(min, max, length.out=nbTick) colTick = paletteShow[idTick] } else if (min > 0 & max > 0) { paletteShow = palette_hot idTick = c() for (i in 1:nbTick) { id = round((ncolor-1)/(nbTick-1)*(i-1)) + 1 idTick = c(idTick, id) } labTick = seq(min, max, length.out=nbTick) colTick = paletteShow[idTick] } else { paletteShow = palette nbSemiTick = round(nbTick/2) + 1 idSemiTick = c() for (i in 1:nbSemiTick) { id = round((ncolor-1)/(nbSemiTick-1)*(i-1)) + 1 idSemiTick = c(idSemiTick, id)
1891189218931894189518961897189818991900190119021903190419051906190719081909191019111912191319141915191619171918191919201921192219231924192519261927192819291930193119321933193419351936193719381939194019411942194319441945194619471948194919501951195219531954195519561957195819591960
} labTick_hot = seq(0, max, length.out=nbSemiTick) labTick_cold = seq(min, 0, length.out=nbSemiTick) colTick_hot = palette_hot[idSemiTick] colTick_cold = palette_cold[idSemiTick] idTick = as.integer(seq(1, ncolor, length.out=nbTick+1)) labTick = c(labTick_cold, labTick_hot[-1]) colTick = c(colTick_cold, colTick_hot[-1]) } return(list(color=color, palette=paletteShow, idTick=idTick, labTick=labTick, colTick=colTick)) } 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() ) palette_tester = function (n=256) { X = 1:n Y = rep(0, times=n) palette = colorRampPalette(c( # '#1a4157', # '#00af9d', # '#fbdd7e', # '#fdb147', # '#fd4659' # '#543005', # '#8c510a', # '#bf812d', # '#dfc27d', # '#f6e8c3', # '#f5f5f5', # '#c7eae5', # '#80cdc1', # '#35978f', # '#01665e', # '#003c30' '#0f3b57', '#1d7881', '#80c4a9', '#e2dac6', #mid '#fadfad', '#d08363', '#7e392f' # '#193830', # '#2A6863', # '#449C93', # '#7ACEB9',
1961196219631964196519661967196819691970197119721973197419751976197719781979198019811982198319841985198619871988198919901991199219931994199519961997199819992000200120022003200420052006200720082009201020112012201320142015201620172018201920202021202220232024202520262027202820292030
# '#BCE6DB', # '#EFE0B0', # '#D4B86A', # '#B3762A', # '#7F4A23', # '#452C1A' ))(n) 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() ) + geom_line(aes(x=X, y=Y), color=palette[X], size=60) + scale_y_continuous(expand=c(0, 0)) ggsave(plot=p, filename=paste('palette_test', '.pdf', sep=''), width=10, height=10, units='cm', dpi=100) } # palette_tester() get_power = function (value) { if (value > 1) { power = nchar(as.character(as.integer(value))) - 1 } else { dec = gsub('0.', '', as.character(value), fixed=TRUE) ndec = nchar(dec) nnum = nchar(as.character(as.numeric(dec))) power = -(ndec - nnum + 1) } return(power) } 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, ...) } gpct = function (pct, L, ref=NULL, shift=FALSE) { if (is.null(ref)) { minL = min(L, na.rm=TRUE) } else { minL = ref } maxL = max(L, na.rm=TRUE) spanL = maxL - minL
203120322033203420352036203720382039
xL = pct/100 * as.numeric(spanL) if (shift) { xL = xL + minL } return (xL) }