## ----setup, include=FALSE-----------------------------------------------------
knitr::opts_chunk$set(
  collapse  = TRUE,
  comment   = "#>",
  fig.width = 7,
  fig.height = 5,
  dpi       = 120
)

## ----load---------------------------------------------------------------------
library(okBATHTUB)

## ----summary------------------------------------------------------------------
ok_reservoir_summary()

## ----lookup-------------------------------------------------------------------
# Partial match
ok_reservoir("Tenkiller")[, c("lake_name", "surface_area_ha", "mean_depth_m",
                              "eco_l3_name")]

# Filter by ecoregion
ct <- ok_reservoir(ecoregion = "Cross Timbers")
nrow(ct)
head(ct[, c("lake_name", "surface_area_ha", "mean_depth_m")], 5)

## ----arcadia_inputs-----------------------------------------------------------
arcadia <- ok_reservoir("Arcadia Lake", exact = TRUE)

result <- ok_load(
    inflow_m3yr   = 45e6,    # ~36,500 acre-feet/yr
    tp_inflow_ugl = 120,     # flow-weighted mean TP
    tn_inflow_ugl = 1800,    # flow-weighted mean TN
    coefficients  = "oklahoma",
    ecoregion     = "Cross Timbers"
  ) |>
  ok_hydraulics(
    surface_area_ha = arcadia$surface_area_ha,
    mean_depth_m    = arcadia$mean_depth_m
  ) |>
  ok_retention() |>
  ok_inlake()    |>
  ok_tsi()

summary(result)

## ----compare_observed---------------------------------------------------------
# Suppose lake monitoring data summary for Arcadia gives these
# growing-season means:
obs <- ok_tsi_observed(
  tp_ugl   = 56,
  chla_ugl = 18,
  secchi_m = 0.95
)

cat("--- Predicted (okBATHTUB) vs Observed ---\n")
cat(sprintf("TP    : pred = %5.1f, obs = %5.1f ug/L\n",
            result$data$tp_inlake_ugl, 56))
cat(sprintf("Chl-a : pred = %5.2f, obs = %5.2f ug/L\n",
            result$data$chla_ugl, 18))
cat(sprintf("Secchi: pred = %5.2f, obs = %5.2f m\n",
            result$data$secchi_m, 0.95))
cat(sprintf("TSI   : pred = %5.1f, obs = %5.1f\n",
            result$data$tsi_mean, obs$tsi_mean))

## ----multi_segment------------------------------------------------------------
# A simple 3-segment Tenkiller-like reservoir
result <- ok_segment_chain(
  inflow_m3yr   = 1e9,
  tp_inflow_ugl = 110,
  tn_inflow_ugl = 1700,
  segments = list(
    list(label = "riverine",    surface_area_ha = 1200, mean_depth_m = 4),
    list(label = "transition",  surface_area_ha = 1800, mean_depth_m = 9),
    list(label = "lacustrine",  surface_area_ha = 2060, mean_depth_m = 25)
  ),
  coefficients = "oklahoma",
  ecoregion    = "Ouachita Mountains"
)

ok_segment_summary(result)

## ----scenario, fig.height=5, eval=requireNamespace("ggplot2", quietly=TRUE)----
baseline <- ok_load(
    inflow_m3yr   = 45e6,
    tp_inflow_ugl = 120,
    coefficients  = "oklahoma",
    ecoregion     = "Cross Timbers"
  ) |>
  ok_hydraulics(surface_area_ha = 716, mean_depth_m = 4.6)

sweep <- ok_scenario_sweep(
  baseline,
  max_reduction_pct = 60,
  step_pct          = 10
)

# Print key columns
print(sweep[, c("tp_reduction_pct", "tp_inflow_ugl",
                "tp_inlake_ugl", "chla_ugl", "tsi_mean",
                "trophic_state")],
      row.names = FALSE)