---
title: "Metrics"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Metrics}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

## Load Packages

The first step is to load the rTwig package. Real Twig works well when paired with packages from the [Tidyverse](https://www.tidyverse.org/), so we will also load the dplyr, tidyr, and ggplot packages to help with data manipulation and visualization, ggpubr for multi-panel plots, and rgl for point cloud plotting.

```{r, message=FALSE, warning=FALSE}
# Load rTwig
library(rTwig)

# Other useful packages
library(rgl)
library(tidyr)
library(dplyr)
library(ggplot2)
library(ggpubr)
```

## Run Real Twig & Calculate Metrics

Next, let's run Real Twig with `run_rtwig()` and calculate our tree metrics with `tree_metrics()`.

```{r, message=FALSE}
# File path to QSM
file <- system.file("extdata/QSM.mat", package = "rTwig")

# Correct QSM cylinders
qsm <- run_rtwig(file, twig_radius = 4.23)

# Calculate detailed tree metrics
metrics <- tree_metrics(qsm$cylinder)
```

```{r, echo=FALSE, warning=FALSE, message=FALSE}
# Future Package Cleanup
future::plan("sequential")
```

#### Stem Taper

```{r, fig.width=6, fig.height=3, fig.align='center', out.width="100%", echo=FALSE}
metrics$stem_taper %>%
  ggplot(aes(x = height_m, y = diameter_cm)) +
  geom_point() +
  stat_smooth(method = "loess", color = "black", formula = y ~ x) +
  theme_classic() +
  labs(
    # title = "Stem Taper",
    x = "Height (m)",
    y = "Diameter (cm)"
  )
```

#### Tree Metrics

```{r, echo=FALSE}
tree_h <- metrics$tree_height_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = height_class_m,
    y = value,
    color = height_class_m,
    fill = height_class_m
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Height Distributions",
    x = "",
    y = "",
    fill = "Height Class (m)",
    color = "Height Class (m)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
tree_d <- metrics$tree_diameter_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = diameter_class_cm,
    y = value,
    color = diameter_class_cm,
    fill = diameter_class_cm
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Diameter Distributions",
    x = "",
    y = "",
    fill = "Diameter Class (cm)",
    color = "Diameter Class (cm)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
tree_z <- metrics$tree_zenith_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = zenith_class_deg,
    y = value,
    color = zenith_class_deg,
    fill = zenith_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Zenith Distributions",
    x = "",
    y = "",
    fill = "Zenith Class (deg)",
    color = "Zenith Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
tree_az <- metrics$tree_azimuth_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = azimuth_class_deg,
    y = value,
    color = azimuth_class_deg,
    fill = azimuth_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Azimuth Distributions",
    x = "",
    y = "",
    fill = "Azimuth Class (deg)",
    color = "Azimuth Class (deg)"
  ) +
  theme_classic() +
  facet_wrap(~type) +
  theme(legend.position = "right")
```

```{r, fig.width=12, fig.height=6, fig.align='center', out.width="100%", echo=FALSE}
ggarrange(tree_d, tree_h, tree_az, tree_z)
# annotate_figure(tree, top = text_grob("Tree Metrics", size = 16, face = "bold"))
```

#### Branch Metrics

```{r, echo=FALSE}
branch_d <- metrics$branch_diameter_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = diameter_class_cm,
    y = value,
    color = diameter_class_cm,
    fill = diameter_class_cm
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Branch Diameter Distributions",
    x = "",
    y = "",
    fill = "Diameter Class (cm)",
    color = "Diameter Class (cm)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
branch_h <- metrics$branch_height_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = height_class_m,
    y = value,
    color = height_class_m,
    fill = height_class_m
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Branch Height Distributions",
    x = "",
    y = "",
    fill = "Height Class (m)",
    color = "Height Class (m)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
branch_ang <- metrics$branch_angle_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = angle_class_deg,
    y = value,
    color = angle_class_deg,
    fill = angle_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Branch Angle Distributions",
    x = "",
    y = "",
    fill = "Angle Class (deg)",
    color = "Angle Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
branch_z <- metrics$branch_zenith_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = zenith_class_deg,
    y = value,
    color = zenith_class_deg,
    fill = zenith_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Branch Zenith Distributions",
    x = "",
    y = "",
    fill = "Zenith Class (deg)",
    color = "Zenith Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
branch_az <- metrics$branch_azimuth_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = azimuth_class_deg,
    y = value,
    color = azimuth_class_deg,
    fill = azimuth_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Branch Azimuth Distributions",
    x = "",
    y = "",
    fill = "Azimuth Class (deg)",
    color = "Azimuth Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
branch_o <- metrics$branch_order_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = branch_order,
    y = value,
    color = branch_order,
    fill = branch_order
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Branch Order Distributions",
    x = "",
    y = "",
    fill = "Branch Order",
    color = "Branch Order"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, fig.width=12, fig.height=8, out.width="100%", fig.align='center', echo=FALSE}
ggarrange(branch_d, branch_h, branch_ang, branch_az, branch_z, branch_o, ncol = 2, nrow = 3)
# annotate_figure(branch, top = text_grob("Branch Metrics", size = 16, face = "bold"))
```

#### Segment Metrics

```{r, echo=FALSE}
segment_d <- metrics$segment_diameter_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 3:5, names_to = "type") %>%
  ggplot(aes(
    x = diameter_class_cm,
    y = value,
    color = diameter_class_cm,
    fill = diameter_class_cm
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Segment Diameter Distributions",
    x = "",
    y = "",
    fill = "Diameter Class (cm)",
    color = "Diameter Class (cm)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
segment_h <- metrics$segment_height_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = height_class_m,
    y = value,
    color = height_class_m,
    fill = height_class_m
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Segment Height Distributions",
    x = "",
    y = "",
    fill = "Height Class (m)",
    color = "Height Class (m)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
segment_ang <- metrics$segment_angle_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = angle_class_deg,
    y = value,
    color = angle_class_deg,
    fill = angle_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Segment Angle Distributions",
    x = "",
    y = "",
    fill = "Angle Class (deg)",
    color = "Angle Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
segment_z <- metrics$segment_zenith_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = zenith_class_deg,
    y = value,
    color = zenith_class_deg,
    fill = zenith_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Segment Zenith Distributions",
    x = "",
    y = "",
    fill = "Zenith Class (deg)",
    color = "Zenith Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
segment_az <- metrics$segment_azimuth_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = azimuth_class_deg,
    y = value,
    color = azimuth_class_deg,
    fill = azimuth_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Segment Azimuth Distributions",
    x = "",
    y = "",
    fill = "Azimuth Class (deg)",
    color = "Azimuth Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, echo=FALSE}
segment_o <- metrics$segment_order_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 3:5, names_to = "type") %>%
  ggplot(aes(
    x = reverse_order,
    y = value,
    color = reverse_order,
    fill = reverse_order
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    # title = "Reverse Order Distributions",
    x = "",
    y = "",
    fill = "Reverse Order",
    color = "Reverse Order"
  ) +
  facet_wrap(~type) +
  theme_classic()
```

```{r, fig.width=12, fig.height=8, out.width="100%", fig.align='center', echo=FALSE}
ggarrange(segment_d, segment_h, segment_ang, segment_az, segment_z, segment_o, ncol = 2, nrow = 3)
# annotate_figure(segment, top = text_grob("Branch Metrics", size = 16, face = "bold"))
```

#### Spreads

```{r, echo=FALSE, fig.width=6, fig.height=6, out.width="100%"}
spreads <- metrics$spreads

# Generate colors from red to blue
n <- length(unique(spreads$height_class))
colors <- data.frame(r = seq(0, 1, length.out = n), g = 0, b = seq(1, 0, length.out = n))

# Convert colors to hex
color_hex <- apply(colors, 1, function(row) {
  rgb(row[1], row[2], row[3])
})

# Assign colors to height classes
height_class_colors <- data.frame(
  height_class = unique(spreads$height_class),
  color = color_hex
)

spreads %>%
  left_join(height_class_colors, by = "height_class") %>%
  ggplot(aes(x = azimuth_deg, y = spread_m, group = height_class, color = color)) +
  geom_line() +
  scale_color_identity() +
  coord_polar(start = 0) +
  labs(x = "", y = "") +
  theme_minimal() +
  theme(axis.text.y = element_blank())
```

#### Vertical Profile

```{r, fig.width=6, fig.height=3, fig.align='center', out.width="100%", echo=FALSE}
metrics$spreads %>%
  group_by(height_class) %>%
  summarize(
    max = max(spread_m),
    mean = mean(spread_m),
    min = min(spread_m)
  ) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(x = height_class, y = value, color = type)) +
  geom_line() +
  labs(
    x = "Height Class",
    y = "Spread (m)",
    color = ""
  ) +
  theme_classic()
```

#### Point Cloud

We can even plot the simulated point cloud with the rgl library, and look at the cylinder connectivity.

```{r, eval=FALSE}
rgl::plot3d(metrics$cloud, aspect = FALSE, decorate = FALSE)
plot_qsm(qsm$cylinder, skeleton = TRUE)
```

## Plotting Code

```{r, eval=FALSE}
# Stem Taper -------------------------------------------------------------------
metrics$stem_taper %>%
  ggplot(aes(x = height_m, y = diameter_cm)) +
  geom_point() +
  stat_smooth(method = "loess", color = "black", formula = y ~ x) +
  theme_classic() +
  labs(
    title = "Stem Taper",
    x = "Height (m)",
    y = "Diameter (cm)"
  )

# Tree Metrics -----------------------------------------------------------------
# Tree Height Distributions
metrics$tree_height_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = height_class_m,
    y = value,
    color = height_class_m,
    fill = height_class_m
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Height Distributions",
    x = "",
    y = "",
    fill = "Height Class (m)",
    color = "Height Class (m)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Tree Diameter Distributions
metrics$tree_diameter_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = diameter_class_cm,
    y = value,
    color = diameter_class_cm,
    fill = diameter_class_cm
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
     title = "Diameter Distributions",
    x = "",
    y = "",
    fill = "Diameter Class (cm)",
    color = "Diameter Class (cm)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Tree Zenith Distributions
metrics$tree_zenith_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = zenith_class_deg,
    y = value,
    color = zenith_class_deg,
    fill = zenith_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Zenith Distributions",
    x = "",
    y = "",
    fill = "Zenith Class (deg)",
    color = "Zenith Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Tree azimuth distributions
metrics$tree_azimuth_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = azimuth_class_deg,
    y = value,
    color = azimuth_class_deg,
    fill = azimuth_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Azimuth Distributions",
    x = "",
    y = "",
    fill = "Azimuth Class (deg)",
    color = "Azimuth Class (deg)"
  ) +
  theme_classic() +
  facet_wrap(~type) +
  theme(legend.position = "right")

# Branch Metrics ---------------------------------------------------------------
# Branch diameter distributions
metrics$branch_diameter_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = diameter_class_cm,
    y = value,
    color = diameter_class_cm,
    fill = diameter_class_cm
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Branch Diameter Distributions",
    x = "",
    y = "",
    fill = "Diameter Class (cm)",
    color = "Diameter Class (cm)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Branch height distributions
metrics$branch_height_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = height_class_m,
    y = value,
    color = height_class_m,
    fill = height_class_m
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Branch Height Distributions",
    x = "",
    y = "",
    fill = "Height Class (m)",
    color = "Height Class (m)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Branch angle distributions
metrics$branch_angle_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = angle_class_deg,
    y = value,
    color = angle_class_deg,
    fill = angle_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Branch Angle Distributions",
    x = "",
    y = "",
    fill = "Angle Class (deg)",
    color = "Angle Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Branch zenith distributions
metrics$branch_zenith_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = zenith_class_deg,
    y = value,
    color = zenith_class_deg,
    fill = zenith_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Branch Zenith Distributions",
    x = "",
    y = "",
    fill = "Zenith Class (deg)",
    color = "Zenith Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Branch azimuth distributions
metrics$branch_azimuth_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = azimuth_class_deg,
    y = value,
    color = azimuth_class_deg,
    fill = azimuth_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Branch Azimuth Distributions",
    x = "",
    y = "",
    fill = "Azimuth Class (deg)",
    color = "Azimuth Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Branch order distributions
metrics$branch_order_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = branch_order,
    y = value,
    color = branch_order,
    fill = branch_order
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Branch Order Distributions",
    x = "",
    y = "",
    fill = "Branch Order",
    color = "Branch Order"
  ) +
  facet_wrap(~type) +

# Segment Metrics --------------------------------------------------------------
# Segment diameter distributions
metrics$segment_diameter_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 3:5, names_to = "type") %>%
  ggplot(aes(
    x = diameter_class_cm,
    y = value,
    color = diameter_class_cm,
    fill = diameter_class_cm
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Segment Diameter Distributions",
    x = "",
    y = "",
    fill = "Diameter Class (cm)",
    color = "Diameter Class (cm)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Segment height distributions
metrics$segment_height_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = height_class_m,
    y = value,
    color = height_class_m,
    fill = height_class_m
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Segment Height Distributions",
    x = "",
    y = "",
    fill = "Height Class (m)",
    color = "Height Class (m)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Segment angle distributions
metrics$segment_angle_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = angle_class_deg,
    y = value,
    color = angle_class_deg,
    fill = angle_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Segment Angle Distributions",
    x = "",
    y = "",
    fill = "Angle Class (deg)",
    color = "Angle Class (deg)"
  ) +
  facet_wrap(~type) +

# Segment zenith distributions
metrics$segment_zenith_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = zenith_class_deg,
    y = value,
    color = zenith_class_deg,
    fill = zenith_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Segment Zenith Distributions",
    x = "",
    y = "",
    fill = "Zenith Class (deg)",
    color = "Zenith Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Segment azimuth distributions
metrics$segment_azimuth_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(
    x = azimuth_class_deg,
    y = value,
    color = azimuth_class_deg,
    fill = azimuth_class_deg
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Segment Azimuth Distributions",
    x = "",
    y = "",
    fill = "Azimuth Class (deg)",
    color = "Azimuth Class (deg)"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Segment order distributions
metrics$segment_order_dist %>%
  mutate(volume_L = volume_m3 * 1000) %>%
  select(-volume_m3) %>%
  relocate(volume_L, .before = area_m2) %>%
  pivot_longer(cols = 3:5, names_to = "type") %>%
  ggplot(aes(
    x = reverse_order,
    y = value,
    color = reverse_order,
    fill = reverse_order
  )) +
  geom_bar(stat = "identity", position = "dodge2") +
  labs(
    title = "Reverse Order Distributions",
    x = "",
    y = "",
    fill = "Reverse Order",
    color = "Reverse Order"
  ) +
  facet_wrap(~type) +
  theme_classic()

# Spreads ----------------------------------------------------------------------
spreads <- metrics$spreads

# Generate colors from red to blue
n <- length(unique(spreads$height_class))
colors <- data.frame(r = seq(0, 1, length.out = n), g = 0, b = seq(1, 0, length.out = n))

# Convert colors to hex
color_hex <- apply(colors, 1, function(row) {
  rgb(row[1], row[2], row[3])
})

# Assign colors to height classes
height_class_colors <- data.frame(
  height_class = unique(spreads$height_class),
  color = color_hex
)

spreads %>%
  left_join(height_class_colors, by = "height_class") %>%
  ggplot(aes(x = azimuth_deg, y = spread_m, group = height_class, color = color)) +
  geom_line() +
  scale_color_identity() +
  coord_polar(start = 0) +
  labs(x = "", y = "") +
  theme_minimal() +
  theme(axis.text.y = element_blank())

# Vertical profile -------------------------------------------------------------
metrics$spreads %>%
  group_by(height_class) %>%
  summarize(
    max = max(spread_m),
    mean = mean(spread_m),
    min = min(spread_m)
  ) %>%
  pivot_longer(cols = 2:4, names_to = "type") %>%
  ggplot(aes(x = height_class, y = value, color = type)) +
  geom_line() +
  labs(
    x = "Height Class",
    y = "Spread (m)",
    color = ""
  ) +
  theme_classic()
```