Extract lake hypsographic curves from bathymetry rasters

hypso_from_raster.R

get_hypso <- function(rsub){  
  maxdepth <-  abs(cellStats(rsub, "min")) # set to "max" if depths are positive

  # define depth intervals by raster resolution  
  min_res   <- 0.5
  depth_int <- -1 * seq(0, round(maxdepth/min_res) * min_res, by = min_res)

  # calculate area of raster between depth intervals
  # reclassify raster based on depth intervals
  # calculate area of each class
  rc <- cut(rsub, breaks = depth_int) %>%
    as.data.frame()
  # drop depth_int entries for missing layers
  depth_int <- depth_int[1:length(depth_int) %in% unique(rc$layer)]

  rc <- rc %>%
    tidyr::drop_na(layer) %>%
    group_by(layer) %>% tally() %>% cumsum() %>%
    mutate(area_m2 = n * 5 * 5) %>%
    mutate(depth_int = c(rev(# add interval midpoints
      as.numeric(na.omit((depth_int + lag(depth_int))/2)) * -1), 0)) %>%
    mutate(area_percent = scales::rescale(area_m2, to = c(0, 100))) %>%
    mutate(depth_percent = scales::rescale(depth_int, to = c(0, 100)))

  rc
}