Skip to content
GitLab
Commit d15e12b3 authored by Heraut Louis's avatar Heraut Louis
Browse files

Reform code

parent 709cea16
master
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 269 additions and 2951 deletions
+269 -2951
plotting/layout.R
+ 225
− 301
View file @ d15e12b3
...@@ -36,12 +36,65 @@ library(dplyr) ...@@ -36,12 +36,65 @@ library(dplyr)
library(grid) library(grid)
library(ggh4x) library(ggh4x)
library(RColorBrewer) library(RColorBrewer)
library(rgdal)
library(shadowtext)
# Sourcing R file # Sourcing R file
source('plotting/panel.R', encoding='latin1') source('plotting/time.R', encoding='latin1')
source('plotting/info.R', encoding='latin1')
source('plotting/map.R', encoding='latin1')
panels_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) { source('plotting/matrix.R', encoding='latin1')
source('plotting/other.R', encoding='latin1')
## 1. PERSONALISATION
### 1.1. 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(),
axis.title.y=element_blank(),
# Axis line
axis.line.x=element_blank(),
axis.line.y=element_blank(),
)
### 1.2. Color palette
palette_perso = c('#0f3b57',
'#1d7881',
'#80c4a9',
'#e2dac6', #mid
'#fadfad',
'#d08363',
'#7e392f')
## 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) {
outfile = "Panels" outfile = "Panels"
if (filename_opt != '') { if (filename_opt != '') {
...@@ -99,34 +152,6 @@ panels_layout = function (df_data, df_meta, layout_matrix, isplot=c('datasheet', ...@@ -99,34 +152,6 @@ panels_layout = function (df_data, df_meta, layout_matrix, isplot=c('datasheet',
missRect = replicate(nbp, missRect) missRect = replicate(nbp, missRect)
}} }}
# Get all different stations code
Code = levels(factor(df_meta$code))
nCode = length(Code)
# print(df_trend)
df_trendtmp = df_trend[[1]]
# print(df_trendtmp)
nPeriod_max = 0
for (code in Code) {
df_trend_code = df_trendtmp[df_trendtmp$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
}
}
list_df2plot = vector(mode='list', length=nbp) list_df2plot = vector(mode='list', length=nbp)
# minTrend = c() # minTrend = c()
# maxTrend = c() # maxTrend = c()
...@@ -150,279 +175,15 @@ panels_layout = function (df_data, df_meta, layout_matrix, isplot=c('datasheet', ...@@ -150,279 +175,15 @@ panels_layout = function (df_data, df_meta, layout_matrix, isplot=c('datasheet',
if ('datasheet' %in% isplot) { 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)
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_trendtmp[df_trendtmp$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) {
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(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(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)
}
} }
if ('matrix' %in% isplot) { if ('matrix' %in% isplot) {
matrice_panel(list_df2plot, df_meta, trend_period, mean_period, matrix_panel(list_df2plot, df_meta, trend_period, mean_period,
slice=12, outdirTmp=outdirTmp, A3=TRUE) slice=12, outdirTmp=outdirTmp, A3=TRUE)
} }
if ('map' %in% isplot) { if ('map' %in% isplot) {
map_panel(list_df2plot, map_panel(list_df2plot,
df_meta, df_meta,
...@@ -436,5 +197,168 @@ panels_layout = function (df_data, df_meta, layout_matrix, isplot=c('datasheet', ...@@ -436,5 +197,168 @@ panels_layout = function (df_data, df_meta, layout_matrix, isplot=c('datasheet',
pdf_combine(input=file.path(outdirTmp, list.files(outdirTmp)), pdf_combine(input=file.path(outdirTmp, list.files(outdirTmp)),
output=file.path(outdir, outfile)) output=file.path(outdir, outfile))
# unlink(outdirTmp, recursive=TRUE) # unlink(outdirTmp, recursive=TRUE)
} }
## 3. COLOR MANAGEMENT
### 3.1. Color on colorbar
get_color = function (value, min, max, ncolor=256, palette_name='perso', reverse=FALSE) {
if (palette_name == 'perso') {
colorList = palette_perso
} 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 (is.na(value)) {
return (NA)
}
maxAbs = max(abs(max), abs(min))
if (value < 0) {
idNorm = (value + maxAbs) / maxAbs
id = round(idNorm*(ncolor - 1) + 1, 0)
color = palette_cold[id]
} else {
idNorm = value / maxAbs
id = round(idNorm*(ncolor - 1) + 1, 0)
color = palette_hot[id]
}
return(color)
}
### 3.2. Colorbar
get_palette = function (min, max, ncolor=256, palette_name='perso', reverse=FALSE, nbTick=10) {
if (palette_name == 'perso') {
colorList = palette_perso
} 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 (min < 0 & max < 0) {
paletteShow = palette_cold
} else if (min > 0 & max > 0) {
paletteShow = palette_hot
} else {
paletteShow = palette
}
posTick = seq(0, 1, length.out=nbTick)
labTick = c()
colTick = c()
for (i in 1:nbTick) {
lab = (i-1)/(nbTick-1) * (max - min) + min
labTick = c(labTick, lab)
col = get_color(lab, min=min, max=max,
ncolor=ncolor,
palette_name=palette_name,
reverse=reverse)
colTick = c(colTick, col)
}
return(list(palette=paletteShow, posTick=posTick,
labTick=labTick, colTick=colTick))
}
### 3.3. Palette tester
palette_tester = function (n=256) {
X = 1:n
Y = rep(0, times=n)
palette = colorRampPalette(palette_perso)(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)
}
## 4. OTHER
### 4.1. Number formatting
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)
}
### 4.2. Pourcentage of variable
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
xL = pct/100 * as.numeric(spanL)
if (shift) {
xL = xL + minL
}
return (xL)
}
plotting/panel.R deleted 100644 → 0
+ 0
− 2605
View file @ 709cea16
# \\\
# Copyright 2021-2022 Louis Hraut*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/panel.R
#
#
# Usefull library
library(ggplot2)
library(scales)
library(qpdf)
library(gridExtra)
library(gridtext)
library(dplyr)
library(grid)
library(ggh4x)
library(RColorBrewer)
library(rgdal)
library(shadowtext)
palette_perso = c('#0f3b57',
'#1d7881',
'#80c4a9',
'#e2dac6', #mid
'#fadfad',
'#d08363',
'#7e392f')
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(),
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, 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)
}
matrice_panel = function (list_df2plot, df_meta, trend_period, mean_period, slice=NULL, outdirTmp='', outnameTmp='matrix', title=NULL, A3=FALSE) {
nbp = length(list_df2plot)
# Get all different stations code
Code = levels(factor(df_meta$code))
nCode = length(Code)
# 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
}
}
# print(nPeriod_trend)
# print(nPeriod_max)
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_trend) {
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)
Periods_trend = c()
NPeriod_trend = c()
Type_trend = list()
Code_trend = c()
Pthresold_trend = c()
TrendMean_trend = c()
DataMean_trend = c()
Fill_trend = c()
Color_trend = c()
for (j in 1:nPeriod_trend) {
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 = 'white'
Pthresold = p_thresold
} else {
fill = 'white'
color = 'grey85'
Pthresold = NA
}
Periods_trend = append(Periods_trend, Periods)
NPeriod_trend = append(NPeriod_trend, j)
Type_trend = append(Type_trend, type)
Code_trend = append(Code_trend, code)
Pthresold_trend = append(Pthresold_trend, Pthresold)
TrendMean_trend = append(TrendMean_trend, trendMean)
DataMean_trend = append(DataMean_trend, dataMean)
Fill_trend = append(Fill_trend, fill)
Color_trend = append(Color_trend, color)
}
}
}
# If there is a 'mean_period'
if (!is.null(mean_period)) {
Periods_mean = c()
NPeriod_mean = c()
Type_mean = list()
Code_mean = c()
DataMean_mean = c()
BreakMean_mean = c()
# Convert 'mean_period' to list
mean_period = as.list(mean_period)
# Number of mean period
nPeriod_mean = length(mean_period)
BreakMean_code = array(rep(1, nPeriod_mean*nbp*nCode),
dim=c(nPeriod_mean, nbp, nCode))
dataMeantmp = array(rep(NA, nbp*nCode),
dim=c(nbp, nCode))
# For all mean period
for (j in 1:nPeriod_mean) {
for (k in 1:nCode) {
code = Code[k]
for (i in 1:nbp) {
df_data = list_df2plot[[i]]$data
type = list_df2plot[[i]]$type
df_data_code = df_data[df_data$code == code,]
# 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
Datemin = min(df_data_code_per$Date)
# Max for the sub period
Datemax = max(df_data_code_per$Date)
Periods = paste(Datemin, Datemax,
sep=' / ')
# Mean of the flow over the sub period
dataMean = mean(df_data_code_per$Qm3s,
na.rm=TRUE)
if (j > 1) {
Break = dataMean - dataMeantmp[i, k]
} else {
Break = NA
}
BreakMean = Break / dataMeantmp[i, k]
BreakMean_code[j, i, k] = BreakMean
dataMeantmp[i, k] = dataMean
Periods_mean = append(Periods_mean, Periods)
NPeriod_mean = append(NPeriod_mean, j)
Type_mean = append(Type_mean, type)
Code_mean = append(Code_mean, code)
DataMean_mean = append(DataMean_mean, dataMean)
BreakMean_mean = append(BreakMean_mean,
BreakMean)
}
}
}
minBreakMean = apply(BreakMean_code, c(1, 2),
min, na.rm=TRUE)
maxBreakMean = apply(BreakMean_code, c(1, 2),
max, na.rm=TRUE)
Fill_mean = c()
Color_mean = c()
ii = 1
for (j in 1:nPeriod_mean) {
for (k in 1:nCode) {
code = Code[k]
for (i in 1:nbp) {
BreakMean = BreakMean_mean[ii]
color_res = get_color(BreakMean,
minBreakMean[j, i],
maxBreakMean[j, i],
palette_name='perso',
reverse=TRUE)
fill = color_res
color = 'white'
Fill_mean = append(Fill_mean, fill)
Color_mean = append(Color_mean, color)
ii = ii + 1
}
}
}
}
if (is.null(slice)) {
slice = nCode
}
firstLetter = levels(factor(substr(Code, 1, 1)))
for (fL in firstLetter) {
print(paste('Matrix for region :', fL))
# Get only station code with the same first letter
subCodefL = Code[substr(Code, 1, 1) == fL]
nsubCodefL = length(subCodefL)
nMat = as.integer(nsubCodefL/slice) + 1
# print(nsubCodefL)
# print(nMat)
for (imat in 1:nMat) {
subCode = subCodefL[(slice*(imat-1)+1):(slice*imat)]
subCode = subCode[!is.na(subCode)]
nsubCode = length(subCode)
CodefL_trend = Code_trend %in% subCode
subPeriods_trend = Periods_trend[CodefL_trend]
subNPeriod_trend = NPeriod_trend[CodefL_trend]
subType_trend = Type_trend[CodefL_trend]
subCode_trend = Code_trend[CodefL_trend]
subPthresold_trend = Pthresold_trend[CodefL_trend]
subTrendMean_trend = TrendMean_trend[CodefL_trend]
subDataMean_trend = DataMean_trend[CodefL_trend]
subFill_trend = Fill_trend[CodefL_trend]
subColor_trend = Color_trend[CodefL_trend]
CodefL_mean = Code_mean %in% subCode
subPeriods_mean = Periods_mean[CodefL_mean]
subNPeriod_mean = NPeriod_mean[CodefL_mean]
subType_mean = Type_mean[CodefL_mean]
subCode_mean = Code_mean[CodefL_mean]
subDataMean_mean = DataMean_mean[CodefL_mean]
subBreakMean_mean = BreakMean_mean[CodefL_mean]
subFill_mean = Fill_mean[CodefL_mean]
subColor_mean = Color_mean[CodefL_mean]
title = df_meta[df_meta$code == subCode[1],]$region_hydro
### Plot ###
height = nsubCode
width = nbp * 2 * nPeriod_trend + nPeriod_trend + nPeriod_mean * nbp + nPeriod_mean + nbp
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(),
plot.margin=margin(t=5, r=5, b=5, l=5, unit="mm")
)
xt = 1 - 6
yt = height + 2
Title = bquote(bold(.(title)))
mat = mat +
annotate("text", x=xt, y=yt,
label=Title,
hjust=0, vjust=1,
size=6, color="#00A3A8")
### Trend ###
for (j in 1:nPeriod_trend) {
Type_trend_per =
subType_trend[subNPeriod_trend == j]
Code_trend_per =
subCode_trend[subNPeriod_trend == j]
Pthresold_trend_per =
subPthresold_trend[subNPeriod_trend == j]
TrendMean_trend_per =
subTrendMean_trend[subNPeriod_trend == j]
DataMean_trend_per =
subDataMean_trend[subNPeriod_trend == j]
Fill_trend_per =
subFill_trend[subNPeriod_trend == j]
Color_trend_per =
subColor_trend[subNPeriod_trend == j]
Xtmp = as.integer(factor(as.character(Type_trend_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_trend_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))))
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 +
gg_circle(r=0.45, xc=X[i], yc=Y[i],
fill=Fill_trend_per[i],
color=Color_trend_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_trend_per)) {
trendMean = TrendMean_trend_per[i]
trendC = signif(trendMean*100, 2)
if (!is.na(Pthresold_trend_per[i])) {
Tcolor = 'white'
} else {
Tcolor = 'grey85'
}
dataMean = signif(DataMean_trend_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.82,
label=bquote(.(type)),
hjust=0.5, vjust=0.5,
size=3.25, color='grey20') +
annotate('text', x=X[i], y=max(Y) + 0.6,
label=bquote('[%.'*ans^{-1}*']'),
hjust=0.5, vjust=0.5,
size=2, color='grey40') +
annotate('text', x=Xm[i], y=max(Y) + 0.82,
label=bquote(''*.(type)),
hjust=0.5, vjust=0.5,
size=3.25, color='grey20') +
annotate('text', x=Xm[i], y=max(Y) + 0.6,
label=bquote('['*m^3*'.'*s^{-1}*']'),
hjust=0.5, vjust=0.5,
size=2, color='grey40')
}
for (k in 1:nsubCode) {
code = subCode[k]
label = Periods_trend[subNPeriod_trend == j
& subCode_trend == code][1]
periodStart = substr(label, 1, 4)
periodEnd = substr(label, 14, 17)
mat = mat +
annotate('text', x=Xc, y=k + 0.13,
label=bquote(bold(.(periodStart))),
hjust=0.5, vjust=0.5,
size=3, color='grey40') +
annotate('text', x=Xc, y=k - 0.13,
label=bquote(bold(.(periodEnd))),
hjust=0.5, vjust=0.5,
size=3, color='grey40')
}
}
### Mean ###
for (j in 1:nPeriod_mean) {
Type_mean_per =
subType_mean[subNPeriod_mean == j]
Code_mean_per =
subCode_mean[subNPeriod_mean == j]
DataMean_mean_per =
subDataMean_mean[subNPeriod_mean == j]
BreakMean_mean_per =
subBreakMean_mean[subNPeriod_mean == j]
Fill_mean_per =
subFill_mean[subNPeriod_mean == j]
Color_mean_per =
subColor_mean[subNPeriod_mean == j]
Xtmp_mean = as.integer(factor(as.character(Type_mean_per)))
Xc_mean = j + (j - 1)*nbp + X[length(X)]
Xm_mean = Xtmp_mean + (j - 1)*nbp + j + X[length(X)]
Xr_mean = Xtmp_mean + (j - 1)*nbp*2 + j + X[length(X)]
Y_mean = as.integer(factor(Code_mean_per))
x = Xc_mean - 0.4
xend = Xm_mean[length(Xm_mean)] + 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 = mean_period[[j]][1]
End = mean_period[[j]][2]
periodName = bquote(bold('Priode')~bold(.(as.character(j+nPeriod_trend))))
mat = mat +
annotate("text", x=x, y=yt,
label=periodName,
hjust=0, vjust=0.5,
size=3, color='grey40')
if (j > 1) {
x = Xr_mean[1] - 0.4
xend = Xr_mean[length(Xr_mean)] + 0.25
mat = mat +
annotate("segment",
x=x, xend=xend,
y=y, yend=yend,
color="grey40", size=0.35)
breakName = bquote(bold('cart')~bold(.(as.character(j-1+nPeriod_trend)))*bold('-')*bold(.(as.character(j+nPeriod_trend))))
mat = mat +
annotate("text", x=x, y=yt,
label=breakName,
hjust=0, vjust=0.5,
size=3, color='grey40')
}
for (i in 1:length(Xm_mean)) {
mat = mat +
gg_circle(r=0.45, xc=Xm_mean[i], yc=Y[i],
fill='white', color='grey40') +
gg_circle(r=0.45, xc=Xc_mean, yc=Y[i],
fill='white', color='grey40')
if (j > 1) {
mat = mat +
gg_circle(r=0.45, xc=Xr_mean[i], yc=Y[i],
fill=Fill_mean_per[i],
color=Color_mean_per[i])
}
}
for (i in 1:length(DataMean_mean_per)) {
dataMean = signif(DataMean_mean_per[i], 2)
mat = mat +
annotate('text', x=Xm_mean[i], y=Y[i],
label=dataMean,
hjust=0.5, vjust=0.5,
size=3, color='grey40')
if (j > 1) {
BreakMean = BreakMean_mean_per[i]
BreakC = signif(BreakMean*100, 2)
mat = mat +
annotate('text', x=Xr_mean[i], y=Y[i],
label=BreakC,
hjust=0.5, vjust=0.5,
size=3, color='white')
}
}
mat = mat +
annotate('text', x=Xc_mean, y=max(Y) + 0.85,
label=bquote(bold('Dbut')),
hjust=0.5, vjust=0.5,
size=3, color='grey20') +
annotate('text', x=Xc_mean, 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=Xm_mean[i], y=max(Y) + 0.82,
label=bquote(''*.(type)),
hjust=0.5, vjust=0.5,
size=3.25, color='grey20') +
annotate('text', x=Xm_mean[i], y=max(Y) + 0.6,
label=bquote('['*m^3*'.'*s^{-1}*']'),
hjust=0.5, vjust=0.5,
size=2, color='grey40')
if (j > 1) {
mat = mat +
annotate('text', x=Xr_mean[i], y=max(Y) + 0.82,
label=bquote('d'*.(type)),
hjust=0.5, vjust=0.5,
size=3.25, color='grey20') +
annotate('text', x=Xr_mean[i], y=max(Y) + 0.6,
label=bquote('[%]'),
hjust=0.5, vjust=0.5,
size=2, color='grey40')
}
}
for (k in 1:nsubCode) {
code = subCode[k]
label = Periods_mean[subNPeriod_mean == j
& subCode_mean == code][1]
periodStart = substr(label, 1, 4)
periodEnd = substr(label, 14, 17)
mat = mat +
annotate('text', x=Xc_mean, y=k + 0.13,
label=bquote(bold(.(periodStart))),
hjust=0.5, vjust=0.5,
size=3, color='grey40') +
annotate('text', x=Xc_mean, y=k - 0.13,
label=bquote(bold(.(periodEnd))),
hjust=0.5, vjust=0.5,
size=3, color='grey40')
}
}
### Code ###
for (k in 1:nsubCode) {
code = subCode[k]
name = df_meta[df_meta$code == code,]$nom
ncharMax = 38
if (nchar(name) > ncharMax) {
name = paste(substr(name, 1, ncharMax), '...', sep='')
}
mat = mat +
annotate('text', x=0.3, y=k + 0.14,
label=bquote(bold(.(code))),
hjust=1, vjust=0.5,
size=3.5, color="#00A3A8") +
annotate('text', x=0.3, y=k - 0.14,
label=name,
hjust=1, vjust=0.5,
size=3.5, color="#00A3A8")
}
### Environment ###
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(2)),
expand=c(0, 0))
# Saving matrix plot
if (A3) {
width = 42
height = 29.7
dpi = 300
} else {
width = 29.7
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')
}
}
}
map_panel = function (list_df2plot, df_meta, df_shapefile, idPer=1, outdirTmp='', codeLight=NULL, margin=NULL, showSea=TRUE, verbose=TRUE) {
df_france = df_shapefile$france
df_bassin = df_shapefile$bassin
df_river = df_shapefile$river
# 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
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) {
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)
ncolor = 256
nbTick = 10
for (i in 1:nbp) {
if (i > 1 & !is.null(codeLight)) {
break
}
type = list_df2plot[[i]]$type
outname = paste('map_', type, sep='')
if (verbose) {
print(paste('Map for variable :', type))
}
if (is.null(codeLight)) {
sizefr = 0.45
sizebs = 0.4
sizerv = 0.3
} else {
sizefr = 0.35
sizebs = 0.3
sizerv = 0.2
}
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")
if (!is.null(df_river)) {
map = map +
geom_path(data=df_river,
aes(x=long, y=lat, group=group),
color="grey85", size=sizerv)
}
map = map +
geom_polygon(data=df_bassin,
aes(x=long, y=lat, group=group),
color="grey70", fill=NA, size=sizebs) +
geom_polygon(data=df_france,
aes(x=long, y=lat, group=group),
color="grey40", fill=NA, size=sizefr)
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)
size = 3
sizekm = 2.5
} else {
xmin = gpct(5, xlim, shift=TRUE)
xint = c(0, 100*1E3)
ymin = gpct(1, ylim, shift=TRUE)
ymax = ymin + gpct(3, ylim)
size = 2
sizekm = 1.5
}
map = map +
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=sizekm)
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=size)
}
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()
trend = c()
p_threshold_Ok = 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,]
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)
palette_res = get_palette(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
if (trendMean >= 0) {
shapetmp = 24
} else {
shapetmp = 25
}
} else {
filltmp = 'grey97'
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)
trend = c(trend, trendMean)
p_threshold_Ok = c(p_threshold_Ok,
df_trend_code_per$p <= p_threshold)
}
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 +
geom_point(data=plot_map_codeNo,
aes(x=lon, y=lat),
shape=21, size=0.5, stroke=0.5,
color='grey70', fill='grey70') +
geom_point(data=plot_map_code,
aes(x=lon, y=lat),
shape=21, size=1.5, stroke=0.5,
color='grey40', fill='grey40')
}
posTick = palette_res$posTick
labTick = palette_res$labTick
colTick = palette_res$colTick
nbTickmod = length(posTick)
valNorm = nbTickmod * 10
ytick = posTick / max(posTick) * 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() +
geom_line(aes(x=c(-0.3, 3.3), y=c(0.05, 0.05)),
size=0.6, color="#00A3A8") +
geom_shadowtext(data=tibble(x=-0.3, y=0.2,
label=type),
aes(x=x, y=y, label=label),
fontface="bold",
color="#00A3A8",
bg.colour="white",
hjust=0, vjust=0, size=10) +
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) {
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')
yTrend = (trend - minTrendMean[idPer, i]) /
(maxTrendMean[idPer, i] - minTrendMean[idPer, i]) * valNorm
yTrend = yTrend[p_threshold_Ok]
## Random distribution ##
# xTrend = rnorm(length(yTrend), mean=1.75, sd=0.1)
## Histogram distribution ##
res_hist = hist(yTrend, breaks=ytick, plot=FALSE)
counts = res_hist$counts
breaks = res_hist$breaks
mids = res_hist$mids
xTrend = c()
yTrend = c()
start_hist = 1.25
hist_sep = 0.15
for (ii in 1:length(mids)) {
if (counts[ii] != 0) {
xTrend = c(xTrend,
seq(start_hist,
start_hist+(counts[ii]-1)*hist_sep,
by=hist_sep))
}
yTrend = c(yTrend, rep(mids[ii], times=counts[ii]))
}
## No touch distribution ##
# start_hist = 1.25
# min_xsep = 0.15
# min_ysep = 4
# xTrend = rep(start_hist, times=length(yTrend))
# for (ii in 1:length(yTrend)) {
# yTrendtmp = yTrend
# yTrendtmp[ii] = 1E99
# isinf_ysep = abs(yTrendtmp - yTrend[ii]) < min_ysep
# if (any(isinf_ysep) & !all(xTrend[which(isinf_ysep)] > start_hist)) {
# xTrend[ii] = max(xTrend[which(isinf_ysep)]) + min_xsep
# }
# }
plot_trend = tibble(xTrend=xTrend, yTrend=yTrend)
pal = pal +
geom_point(data=plot_trend,
aes(x=xTrend, y=yTrend),
# shape=21, size=1,
# color="grey20", fill="grey20")
alpha=0.4)
pal = pal +
geom_segment(aes(x=2.7, y=valNorm*0.75,
xend=2.7, yend=valNorm*0.25),
color='grey50', size=0.3,
arrow=arrow(length=unit(2, "mm"))) +
annotate('text',
x=2.8, y=valNorm*0.5,
label= "Plus svre",
angle=90,
hjust=0.5, vjust=1,
size=3, color='grey50')
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)
}
}
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("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)
}
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("-", "-&nbsp;", 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>",
"Rgion 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 dbut : ", debut, "<br>",
"Date de fin : ", fin, "<br>",
"Nombre d'annes : ", 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)
}
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)
}
# 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)
# }
# 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))
# }
get_color = function (value, min, max, ncolor=256, palette_name='perso', reverse=FALSE) {
if (palette_name == 'perso') {
colorList = palette_perso
} 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 (is.na(value)) {
return (NA)
}
maxAbs = max(abs(max), abs(min))
if (value < 0) {
idNorm = (value + maxAbs) / maxAbs
id = round(idNorm*(ncolor - 1) + 1, 0)
color = palette_cold[id]
} else {
idNorm = value / maxAbs
id = round(idNorm*(ncolor - 1) + 1, 0)
color = palette_hot[id]
}
return(color)
}
get_palette = function (min, max, ncolor=256, palette_name='perso', reverse=FALSE, nbTick=10) {
if (palette_name == 'perso') {
colorList = palette_perso
} 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 (min < 0 & max < 0) {
paletteShow = palette_cold
} else if (min > 0 & max > 0) {
paletteShow = palette_hot
} else {
paletteShow = palette
}
posTick = seq(0, 1, length.out=nbTick)
labTick = c()
colTick = c()
for (i in 1:nbTick) {
lab = (i-1)/(nbTick-1) * (max - min) + min
labTick = c(labTick, lab)
col = get_color(lab, min=min, max=max,
ncolor=ncolor,
palette_name=palette_name,
reverse=reverse)
colTick = c(colTick, col)
}
return(list(palette=paletteShow, posTick=posTick,
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(palette_perso)(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
xL = pct/100 * as.numeric(spanL)
if (shift) {
xL = xL + minL
}
return (xL)
}
processing/format.R
+ 1
− 1
View file @ d15e12b3
...@@ -92,7 +92,7 @@ reprepare = function(df_XEx, df_Xlist, colnamegroup=NULL) { ...@@ -92,7 +92,7 @@ reprepare = function(df_XEx, df_Xlist, colnamegroup=NULL) {
df_XEx$Date = paste(df_XEx$Date, '01', '01', sep='-') df_XEx$Date = paste(df_XEx$Date, '01', '01', sep='-')
# If there is more than 2 dashes # If there is more than 2 dashes
} else if (nbt != 2) { } else if (nbt != 2) {
This is not a classical date # This is not a classical date
stop('erreur of date format') stop('erreur of date format')
} }
# Recreates the outing of the 'extract.Var' function nicer # Recreates the outing of the 'extract.Var' function nicer
......
script.R
+ 43
− 44
View file @ d15e12b3
...@@ -57,17 +57,17 @@ filedir = ...@@ -57,17 +57,17 @@ filedir =
filename = filename =
# "" # ""
c( # c(
"S2235610_HYDRO_QJM.txt", # "S2235610_HYDRO_QJM.txt",
"P1712910_HYDRO_QJM.txt", # "P1712910_HYDRO_QJM.txt",
"P0885010_HYDRO_QJM.txt", # "P0885010_HYDRO_QJM.txt",
"O5055010_HYDRO_QJM.txt", # "O5055010_HYDRO_QJM.txt",
"O0384010_HYDRO_QJM.txt" # "O0384010_HYDRO_QJM.txt"
) # )
# c("S4214010_HYDRO_QJM.txt", c("S4214010_HYDRO_QJM.txt",
# "O0384010_HYDRO_QJM.txt", "O0384010_HYDRO_QJM.txt",
# "Q7002910_HYDRO_QJM.txt") "Q7002910_HYDRO_QJM.txt")
...@@ -130,7 +130,6 @@ setwd(computer_work_path) ...@@ -130,7 +130,6 @@ setwd(computer_work_path)
source('processing/extract.R', encoding='latin1') source('processing/extract.R', encoding='latin1')
source('processing/format.R', encoding='latin1') source('processing/format.R', encoding='latin1')
source('processing/analyse.R', encoding='latin1') source('processing/analyse.R', encoding='latin1')
source('plotting/panel.R', encoding='latin1')
source('plotting/layout.R', encoding='latin1') source('plotting/layout.R', encoding='latin1')
# Result directory # Result directory
...@@ -256,39 +255,39 @@ df_shapefile = ini_shapefile(computer_data_path, fr_shpdir, fr_shpname, bs_shpdi ...@@ -256,39 +255,39 @@ df_shapefile = ini_shapefile(computer_data_path, fr_shpdir, fr_shpname, bs_shpdi
### 4.1. Simple time panel to criticize station data ### 4.1. Simple time panel to criticize station data
# Plot time panel of debit by stations # Plot time panel of debit by stations
# panels_layout(list(df_data, df_data), # datasheet_layout(list(df_data, df_data),
# layout_matrix=c(1, 2), # layout_matrix=c(1, 2),
# df_meta=df_meta, # df_meta=df_meta,
# missRect=list(TRUE, TRUE), # missRect=list(TRUE, TRUE),
# type=list('Q', 'sqrt(Q)'), # type=list('Q', 'sqrt(Q)'),
# info_header=TRUE, # info_header=TRUE,
# time_header=NULL, # time_header=NULL,
# var_ratio=3, # var_ratio=3,
# figdir=figdir, # figdir=figdir,
# filename_opt='time') # filename_opt='time')
### 4.2. Analysis layout ### 4.2. Analysis layout
panels_layout(isplot=c( datasheet_layout(isplot=c(
'datasheet', 'datasheet',
'matrix', 'matrix',
'map' 'map'
), ),
df_data=list(res_QAtrend$data, res_QMNAtrend$data, df_data=list(res_QAtrend$data, res_QMNAtrend$data,
res_VCN10trend$data), res_VCN10trend$data),
layout_matrix=c(1, 2, 3), layout_matrix=c(1, 2, 3),
df_meta=df_meta, df_meta=df_meta,
df_trend=list(res_QAtrend$trend, res_QMNAtrend$trend, df_trend=list(res_QAtrend$trend, res_QMNAtrend$trend,
res_VCN10trend$trend), res_VCN10trend$trend),
type=list('QA', 'QMNA', 'VCN10'), type=list('QA', 'QMNA', 'VCN10'),
missRect=list(TRUE, TRUE, TRUE), missRect=list(TRUE, TRUE, TRUE),
trend_period=trend_period, trend_period=trend_period,
mean_period=mean_period, mean_period=mean_period,
info_header=TRUE, info_header=TRUE,
time_header=df_data, time_header=df_data,
info_ratio=2, info_ratio=2,
time_ratio=2, time_ratio=2,
var_ratio=3, var_ratio=3,
df_shapefile=df_shapefile, df_shapefile=df_shapefile,
figdir=figdir, figdir=figdir,
filename_opt='') filename_opt='')
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment