Skip to content

Commit

Permalink
Revision R1
Browse files Browse the repository at this point in the history
  • Loading branch information
@PMassicotte
PMassicotte committed Nov 7, 2018
1 parent c663307 commit 842793e
Show file tree
Hide file tree
Showing 75 changed files with 69,705 additions and 15,771 deletions.
8 changes: 8 additions & 0 deletions .vscode/settings.json
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,8 @@
{
"spellright.language": [],
"spellright.documentTypes": [
"markdown",
"plaintext",
"rmd"
]
}
8 changes: 4 additions & 4 deletions R/effect_scattering_par.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,12 +27,12 @@ read_hydrolight <- function(file, variable) {
return(df)
}

ed_without_fluo <- read_hydrolight("~/Desktop/phil_sans_raman.xlsx", "Ed")
lu_without_fluo <- read_hydrolight("~/Desktop/phil_sans_raman.xlsx", "Lu")
ed_without_fluo <- read_hydrolight("data/raw/hydrolight/phil_sans_raman.xlsx", "Ed")
lu_without_fluo <- read_hydrolight("data/raw/hydrolight/phil_sans_raman.xlsx", "Lu")
without_fluo <- inner_join(ed_without_fluo, lu_without_fluo)

ed_with_fluo <- read_hydrolight("~/Desktop/phil_avec_raman.xlsx", "Ed")
lu_with_fluo <- read_hydrolight("~/Desktop/phil_avec_raman.xlsx", "Lu")
ed_with_fluo <- read_hydrolight("data/raw/hydrolight/phil_avec_raman.xlsx", "Ed")
lu_with_fluo <- read_hydrolight("data/raw/hydrolight/phil_avec_raman.xlsx", "Lu")
with_fluo <- inner_join(ed_with_fluo, lu_with_fluo)

df <- bind_rows(without_fluo, with_fluo) %>%
Expand Down
91 changes: 82 additions & 9 deletions R/fig1.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,9 @@

rm(list = ls())


# Plot 1 ------------------------------------------------------------------

circle <- function(radius, angle = seq(0, 2 * pi, length.out = 200)) {
tibble(
x = radius * cos(angle),
Expand All @@ -29,25 +32,95 @@ configuration <- map2(sqrt(25 / c(0.25, 0.20, 0.15, 0.10, 0.05, 0.01)), seq(0, p
bind_rows(.id = "radius") %>%
mutate(radius = parse_number(radius))

p <- ggplot() +
frame <- tibble(xmin = -50, xmax = 50, ymin = -50, ymax = 50)

p1 <- ggplot() +
geom_rect(data = frame, aes(xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax), fill = "transparent", color = "#57e5d6") +
geom_path(data = sampling_circle, aes(x = x, y = y, color = factor(radius)), size = 0.25) +
geom_polygon(data = melt_pond, aes(x = x, y = y, fill = "Melt pond")) +
scale_x_continuous(limits = c(-50, 50)) +
scale_y_continuous(limits = c(-50, 50)) +
scale_fill_manual(values = c("Melt pond" = "#36454F")) +
scale_fill_manual(values = c("Melt pond" = "#6097ce")) +
# scale_color_brewer(palette = "BuPu") +
coord_equal() +
labs(fill = "") +
labs(color = "Sampling\ndistance (m)") +
xlab("x-distance (meters)") +
ylab("y-distance (meters)") +
xlab("x-distance (m)") +
ylab("y-distance (m)") +
geom_segment(data = configuration, aes(x = 0, y = 0, xend = x, yend = y, group = radius), size = 0.25) +
geom_label(data = configuration, aes(x = x, y = y, label = sprintf("%2.2f%% melt pond cover", radius * 100)), hjust = 0, label.size = 0, nudge_x = 1, size = 2, label.padding = unit(0.05, "lines")) +
geom_label(data = configuration, aes(x = x, y = y, label = sprintf("%2.0f%% melt pond cover", radius * 100)), hjust = 0, label.size = 0, nudge_x = 1, size = 2, label.padding = unit(0.05, "lines")) +
guides(color = guide_legend(
keywidth = 0.15,
keyheight = 0.15,
default.unit = "inch",
ncol = 1
))
nrow = 4
)) +
theme(legend.title = element_text(size = 8),
legend.text = element_text(size = 8)) +
theme(legend.position = "top", legend.direction = "horizontal")


# Plot 2 ------------------------------------------------------------------

frame <- tibble(xmin = -120, xmax = 120, ymin = -150, ymax = 0)
segment <- tibble(x = rep(-120, 6), y = seq(-25, 0, by = 5), xend = rep(120, 6), yend = seq(-25, 0, by = 5))
ice <- tibble(xmin = -120, xmax = 120.25, ymin = 0, ymax = 5)
mp <- tibble(xmin = -5, xmax = 5, ymin = 0, ymax = 5)

p2 <- ggplot() +
geom_rect(data = frame, aes(xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax), fill = "transparent", color = "black") +
geom_segment(data = segment, aes(x = x, y = y, xend = xend, yend = yend), lty = 2, size = 0.25) +
geom_rect(data = ice, aes(xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax), fill = "#36454F") +
geom_rect(data = mp, aes(xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax), fill = "#6097ce", color = "transparent") +
annotate("segment", x = -50, xend = -50, y = 0, yend = -25, lty = 1, size = 0.75, color = "#57e5d6") +
annotate("segment", x = 50, xend = 50, y = 0, yend = -25, lty = 1, size = 0.75, color = "#57e5d6") +
annotate("segment", x = 50, xend = -50, y = -25, yend = -25, lty = 1, size = 0.75, color = "#57e5d6") +

annotate("text", x = 35, y = -50, label = str_wrap("3D volume for which simulated data was extracted", width = 30), size = 3, family = "IBM Plex Sans Light") +
annotate("segment", x = 25, xend = 25, y = -40, yend = -28, lty = 1, size = 0.25, arrow = arrow(length = unit(1.5, "mm"))) +

annotate("text", x = -75, y = -50, label = str_wrap("2D detectors placed every 0.5 m depth", width = 20), size = 3, family = "IBM Plex Sans Light") +
annotate("segment", x = -75, xend = -75, y = -40, yend = -28, lty = 1, size = 0.25, arrow = arrow(length = unit(1.5, "mm"))) +

# annotate("label", x = 28, y = -5, label = "Melt pond", size = 3, label.size = NA) +
# annotate("segment", x = 10, xend = 5, y = -5, yend = -2, lty = 1, size = 0.25, arrow = arrow(length = unit(1.5, "mm"))) +

annotate("segment", x = 0, xend = 0, y = 0, yend = -15, lty = 1, size = 1, arrow = arrow(length = unit(1.5, "mm")), color = "orange") +
annotate("segment", x = 0, xend = 10, y = 0, yend = -10, lty = 1, size = 1, arrow = arrow(length = unit(1.5, "mm")), color = "orange") +
annotate("segment", x = 0, xend = -10, y = 0, yend = -10, lty = 1, size = 1, arrow = arrow(length = unit(1.5, "mm")), color = "orange") +

annotate("segment", x = 75, xend = 75, y = 0, yend = -15, lty = 1, size = 0.25, arrow = arrow(length = unit(1.5, "mm")), color = "orange") +
annotate("segment", x = 75, xend = 85, y = 0, yend = -10, lty = 1, size = 0.25, arrow = arrow(length = unit(1.5, "mm")), color = "orange") +
annotate("segment", x = 75, xend = 65, y = 0, yend = -10, lty = 1, size = 0.25, arrow = arrow(length = unit(1.5, "mm")), color = "orange") +

annotate("segment", x = -75, xend = -75, y = 0, yend = -15, lty = 1, size = 0.25, arrow = arrow(length = unit(1.5, "mm")), color = "orange") +
annotate("segment", x = -75, xend = -85, y = 0, yend = -10, lty = 1, size = 0.25, arrow = arrow(length = unit(1.5, "mm")), color = "orange") +
annotate("segment", x = -75, xend = -65, y = 0, yend = -10, lty = 1, size = 0.25, arrow = arrow(length = unit(1.5, "mm")), color = "orange") +

annotate("text", x = -110, y = -100, label = "Homogeneous water column:", size = 3.5, hjust = 0, fontface = 2, family = "IBM Plex Sans Light") +
annotate("text", x = -110, y = -110, label = expression("a = b = 0.05"*m^{-1}), size = 3, hjust = 0, family = "IBM Plex Sans Light") +
annotate("text", x = -110, y = -120, label = "VSF: Fourrier-Forand 3%", size = 3, hjust = 0, family = "IBM Plex Sans Light") +

theme(panel.grid = element_blank()) +
theme(panel.border = element_blank()) +
scale_x_continuous(expand = c(0, 0), limits = c(-120, 120.25)) +
scale_y_continuous(expand = c(0, 0), breaks = seq(-150, 0, by = 20)) +
xlab("Horizontal distance (m)") +
ylab("Depth (m)")


# Save --------------------------------------------------------------------

p <- cowplot::plot_grid(p1, p2, ncol = 1, labels = "AUTO", rel_heights = c(1, 0.7))

ggsave("graphs/fig1.pdf", device = cairo_pdf, width = 4, height = 8)



# Test --------------------------------------------------------------------

ggsave("graphs/fig1.pdf", device = cairo_pdf, width = 5, height = 5)
# p2 +
# geom_path(data = filter(sampling_circle, y > 0), aes(x = x, y = y, color = factor(radius)), size = 0.25) +
# geom_segment(data = configuration, aes(x = 0, y = 0, xend = x, yend = y, group = radius), size = 0.25) +
# geom_label(data = configuration, aes(x = x, y = y, label = sprintf("%2.0f%% melt pond cover", radius * 100)), hjust = 0, label.size = 0, nudge_x = 1, size = 2, label.padding = unit(0.05, "lines"))
#
# ggsave("graphs/fig1.pdf", device = cairo_pdf, width = 12, height = 8)
90 changes: 26 additions & 64 deletions R/fig2.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3,75 +3,37 @@
#
# DESCRIPTION:
#
# Fig 1 of the article.
# Angular distribution (Simon's data)
# <><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>

rm(list = ls())

source("https://gist.githubusercontent.com/friendly/67a7df339aa999e2bcfcfec88311abfc/raw/761a7688fba3668a84b2dfe42a655a1b246ca193/wavelength_to_rgb.R")
df1 <- readxl::read_excel("data/raw/angular_distribution.xlsx", n_max = 2, na = "NaN") %>%
janitor::remove_empty() %>%
gather(angle, value) %>%
set_names(c("angle", "Measured under-ice\ndownward radiance distribution")) %>%
mutate_all(parse_number)

cops <- read_feather("data/clean/cops_wavelength_interpolated.feather")
df2 <- readxl::read_excel("data/raw/angular_distribution.xlsx", skip = 4, n_max = 2) %>%
janitor::remove_empty() %>%
gather(angle, value) %>%
set_names(c("angle", "Emitting source chosen\nfor the simulation")) %>%
mutate_all(parse_number)

df <- cops %>%
filter(profile_filename == "GE2016.ICMP_ICEP_160620_CAST_006") %>%
dplyr::select(depth, wavelength, edz, luz) %>%
gather(type, light, -depth, -wavelength) %>%
filter(wavelength %in% seq(400, 700, by = 20))
df <- bind_rows(df1, df2) %>%
gather(source, value, -angle) %>%
drop_na()

color <- lapply(unique(df$wavelength), wavelength_to_rgb) %>% unlist()
color <- setNames(color, unique(df$wavelength))

mylabel <- c(
edz = "Downwelling irradiance (Ed)",
luz = "Upwelling radiance (Lu)"
)

p1 <- df %>%
drop_na() %>%
filter(type == "edz") %>%
ggplot(aes(light, depth, color = factor(wavelength))) +
geom_path() +
scale_y_reverse() +
facet_wrap(~type, scales = "free", labeller = labeller(type = mylabel)) +
ylab("Depth (m)") +
xlab(bquote(Ed~"("*mu*W%*%cm^{-2}%*%nm^{-1}*")")) +
labs(color = "Wavelengths (nm)") +
theme(legend.position = "none") +
guides(color = guide_legend(
keywidth = 0.15,
keyheight = 0.15,
default.unit = "inch",
ncol = 2
)) +
scale_color_manual(values = color)

p2 <- df %>%
drop_na() %>%
filter(type == "luz") %>%
ggplot(aes(light, depth, color = factor(wavelength))) +
df %>%
ggplot(aes(x = angle, y = value, color = source)) +
geom_path() +
scale_y_reverse() +
facet_wrap(~type, scales = "free", labeller = labeller(type = mylabel)) +
ylab("Depth (m)") +
xlab(bquote(Lu~"("*mu*W%*%cm^{-2}%*%nm^{-1}%*%sr^{-1}*")")) +
labs(color = "Wavelengths (nm)") +
theme(
legend.position = c(0.99, 0.01),
legend.justification = c(1, 0)
) +
theme(
legend.title = element_text(size = 10),
legend.text = element_text(size = 8)
) +
guides(color = guide_legend(
keywidth = 0.15,
keyheight = 0.15,
default.unit = "inch",
ncol = 3
)) +
scale_color_manual(values = color)

p <- cowplot::plot_grid(p1, p2, ncol = 2, labels = "AUTO", align = "hv")

ggsave("graphs/fig2.pdf", plot = p, device = cairo_pdf, height = 3, width = 7)

xlab(expression("Azimuthal angle"~"("^degree*")")) +
ylab("Normalized raidance") +
theme(legend.title = element_blank()) +
theme(legend.justification = c(0, 0)) +
theme(legend.position = c(0.01, 0.01)) +
scale_x_continuous(limits = c(0, 90), breaks = seq(0, 90, by = 10)) +
scale_y_continuous(limits = c(0, 1), breaks = seq(0, 1, by = 0.1)) +
theme(legend.key.height = unit(1, "cm"))

ggsave("graphs/fig2.pdf", device = cairo_pdf, width = 3 * 1.61803398875, height = 3)
99 changes: 42 additions & 57 deletions R/fig3.R
View file
Original file line number Diff line number Diff line change
@@ -1,40 +1,57 @@
# <><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>
# <><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>
# AUTHOR: Philippe Massicotte
#
# DESCRIPTION:
# DESCRIPTION:
#
# Fig. 2 showing that Ed and Lu light profiles follow roughly the same patterns
# after 10 m.
# Fig 3 of the article.
# <><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>

rm(list = ls())

source("https://gist.githubusercontent.com/friendly/67a7df339aa999e2bcfcfec88311abfc/raw/761a7688fba3668a84b2dfe42a655a1b246ca193/wavelength_to_rgb.R")

# https://stackoverflow.com/questions/27435453/how-to-normalise-subgroups-from-a-grouped-data-frame-in-r
cops <- read_feather("data/clean/cops_wavelength_interpolated.feather")

df <- read_feather("data/clean/cops_wavelength_interpolated.feather") %>%
filter(profile_filename == "GE2016.ICMP_ICEP_160620_CAST_006") %>%
filter(between(wavelength, 400, 700)) %>%
filter(depth >= 10) %>%
group_by(wavelength) %>%
mutate_at(.vars = vars(edz, luz), funs(. / .[depth == 10])) %>%
drop_na(depth, edz, luz)
df <- cops %>%
filter(profile_filename == "GE2016.ICMP_ICEP_160620_CAST_006") %>%
dplyr::select(depth, wavelength, edz, luz) %>%
gather(type, light, -depth, -wavelength) %>%
filter(wavelength %in% seq(400, 700, by = 20))

color <- lapply(unique(df$wavelength), wavelength_to_rgb) %>% unlist()
color <- setNames(color, unique(df$wavelength))

p <- df %>%
filter(wavelength %in% c(seq(400, 700, by = 20))) %>%
ggplot(aes(y = depth, color = factor(wavelength))) +
geom_path(aes(x = edz, linetype = "Downwelling irradiance (Ed)")) +
geom_path(aes(x = luz, linetype = "Upwelling radiance (Lu)")) +
scale_y_reverse(limits = c(NA, 0)) +
scale_x_continuous(labels = scales::percent) +
p1 <- df %>%
drop_na() %>%
filter(type == "edz") %>%
ggplot(aes(light, depth, color = factor(wavelength))) +
geom_path() +
scale_y_reverse() +
ylab("Depth (m)") +
xlab(bquote(E[d]~"("*mu*W~cm^{-2}*")")) +
labs(color = "Wavelengths (nm)") +
theme(legend.position = "none") +
guides(color = guide_legend(
keywidth = 0.15,
keyheight = 0.15,
default.unit = "inch",
ncol = 2
)) +
scale_color_manual(values = color)

p2 <- df %>%
drop_na() %>%
filter(type == "luz") %>%
ggplot(aes(light, depth, color = factor(wavelength))) +
geom_path() +
scale_y_reverse() +
ylab("Depth (m)") +
xlab(bquote(L[u]~"("*mu*W~cm^{-2}~sr^{-1}*")")) +
labs(color = "Wavelengths (nm)") +
theme(
legend.position = c(0.99, -0.05),
legend.justification = c(1, 0)
) +
legend.position = c(0.99, 0.01),
legend.justification = c(1, 0)
) +
theme(
legend.title = element_text(size = 10),
legend.text = element_text(size = 8)
Expand All @@ -45,41 +62,9 @@ p <- df %>%
default.unit = "inch",
ncol = 3
)) +
labs(
color = "Wavelengths (nm)",
linetype = "Light type"
) +
scale_color_manual(values = color) +
facet_wrap(~wavelength, scales = "free", ncol = 3) +
xlab("Normalized light") +
ylab("Depth (m)") +
theme(legend.box = "horizontal")
scale_color_manual(values = color)

ggsave("graphs/fig3.pdf", device = cairo_pdf, height = 8, width = 6)
p <- cowplot::plot_grid(p1, p2, ncol = 2, labels = "AUTO", align = "hv")

ggsave("graphs/fig3.pdf", plot = p, device = cairo_pdf, height = 3, width = 7)

# Correlation coefficients ------------------------------------------------

# df <- read_feather("data/clean/cops_wavelength_interpolated.feather") %>%
# filter(between(wavelength, 400, 700)) %>%
# filter(depth >= 10) %>%
# group_by(profile_filename, wavelength) %>%
# mutate_at(.vars = vars(edz, luz), funs(. / .[depth == 10])) %>%
# drop_na(depth, edz, luz)
#
# res <- df %>%
# group_by(profile_filename, wavelength) %>%
# summarise(r = cor(edz, luz)) %>%
# group_by(wavelength) %>%
# summarise(mean_r = mean(r), sd_r = sd(r), n = n())
#
# res %>%
# ggplot(aes(x = wavelength, y = mean_r)) +
# geom_line() +
# geom_point() +
# geom_ribbon(aes(ymin = mean_r - sd_r, ymax = mean_r + sd_r), alpha = 0.25) +
# xlab("Wavelength (nm)") +
# ylab("Mean Pearson correlation\ncoefficient (r)") +
# scale_x_continuous(breaks = seq(400, 700, by = 50))
#
# ggsave("graphs/fig2b.pdf", device = cairo_pdf, height = 3, width = 4)
Loading

0 comments on commit 842793e

Please sign in to comment.