nhdplusTools: NHDPlus Tools

test_that("rescale", {
  skip_on_cran()
  skip_on_os("mac")

  old_opts <- options(arrow.unsafe_metadata = TRUE)

  vars <- data.frame(characteristic_id = c("CAT_EWT", "CAT_EWT", "CAT_EWT", "CAT_EWT", "CAT_BASIN_AREA"),
                     summary_statistic = c("area_weighted_mean", "min", "sum", "max", "sum"))

  # file_name <- "refactor_02.gpkg"
  # file_path <- file.path(tempdir(), file_name)
  # gpkg_file <- sbtools::item_file_download("61fbfdced34e622189cb1b0a",
  #                                          destinations = file_path,
  #                                          names = file_name,
  #                                          overwrite_file = TRUE)
  #
  # lookup_table <- sf::st_read(gpkg_file, layer= "lookup_table") |>
  #   dplyr::rename(comid = NHDPlusV2_COMID, member_comid = member_COMID, id = reconciled_ID) |>
  #   # subset lookup table for this example
  #   dplyr::filter(id %in% c(10012268, 10012979, 10024047, 10024048, 10024049, 10024050))
  # lookup_table
  #
  # split_divides <- sf::st_read(gpkg_file, layer = "split_divides") |>
  #   dplyr::rename("featureid" = "FEATUREID") |>
  #   dplyr::filter(as.integer(featureid) %in% lookup_table$comid)
  #
  #
  # catchment_characteristic <- get_catchment_characteristics(varname = vars$characteristic_id,
  #                                                           ids = unique(lookup_table$comid))
  #
  # catchment_areas <- nhdplusTools:::get_catchment_areas(unique(lookup_table$member_comid), split_divides)
  #
  # saveRDS(list(lookup_table = lookup_table,
  #              split_divides = split_divides,
  #              catchment_characteristic = catchment_characteristic,
  #              catchment_areas = catchment_areas),
  #         "tests/testthat/data/rescale_data.rds")

  d <- readRDS(list.files(pattern = "rescale_data.rds", recursive = TRUE, full.names = TRUE))

  suppressWarnings(
    rescale <- rescale_catchment_characteristics(vars, d$lookup_table, d$split_divides))

  suppressWarnings(
    rescale_2 <- rescale_catchment_characteristics(vars, d$lookup_table, d$split_divides,
                                                   d$catchment_characteristic, d$catchment_areas))
  expect_true(is.data.frame(rescale))

  expect_equal(length(unique(d$lookup_table$id)), nrow(rescale))

  expect_equal(rescale, rescale_2)

  expect_equal(round(rescale$areasqkm_sum,0), round(rescale$CAT_BASIN_AREA_sum,0))

  test_id1 <- 10012268
  comids1 <- filter(d$lookup_table, id == test_id1)
  vars_comids1 <- left_join(x = left_join(x = comids1,
                                          y = filter(d$catchment_characteristic,
                                                     comid %in% comids1$comid,
                                                     characteristic_id == "CAT_EWT"),
                                          by = "comid"),
                            y = select(d$catchment_areas, c("member_comid","split_catchment_areasqkm","split_area_prop")),
                            by = "member_comid")
  vars_comids1 <- mutate(vars_comids1, area_rescaled = split_catchment_areasqkm*split_area_prop)
  rescale_test1 <- filter(rescale, id == test_id1)
  expect_equal(min(vars_comids1$characteristic_value),rescale_test1$CAT_EWT_min)
  expect_equal(sum(vars_comids1$characteristic_value),rescale_test1$CAT_EWT_sum)
  expect_equal(max(vars_comids1$characteristic_value),rescale_test1$CAT_EWT_max)
  expect_equal(weighted.mean(vars_comids1$characteristic_value, vars_comids1$area_rescaled),
               rescale_test1$CAT_EWT_area_wtd)

  test_id2 <- 10024048
  comids2 <- filter(d$lookup_table, id == test_id2)
  vars_comids2 <- left_join(x = left_join(x = comids2,
                                          y = filter(d$catchment_characteristic,
                                                     comid %in% comids2$comid,
                                                     characteristic_id == "CAT_EWT"),
                                          by = "comid"),
                            y = select(d$catchment_areas, c("member_comid","split_catchment_areasqkm","split_area_prop")),
                            by = "member_comid")
  vars_comids2 <- mutate(vars_comids2, area_rescaled = split_catchment_areasqkm*split_area_prop)
  rescale_test2 <- filter(rescale, id == test_id2)
  expect_equal(min(vars_comids2$characteristic_value),rescale_test2$CAT_EWT_min)
  expect_equal(sum(vars_comids2$characteristic_value*vars_comids2$split_area_prop),rescale_test2$CAT_EWT_sum)
  expect_equal(max(vars_comids2$characteristic_value),rescale_test2$CAT_EWT_max)
  expect_equal(weighted.mean(vars_comids2$characteristic_value*vars_comids2$split_area_prop, vars_comids2$area_rescaled),
               rescale_test2$CAT_EWT_area_wtd)

  borked <- dplyr::rename(d$lookup_table, borked = "member_comid")

  expect_error(rescale_catchment_characteristics(vars, borked, d$split_divides),
               "Check that lookup_table")

  borked <- dplyr::rename(vars, borked = "characteristic_id")

  expect_error(rescale_catchment_characteristics(borked, d$lookup_table, d$split_divides),
               "Check that vars contains")

  borked <- vars
  borked$summary_statistic[1] <- "borked"

  expect_error(rescale_catchment_characteristics(borked, d$lookup_table, d$split_divides),
               "Check that all entries in")

  expect_error(rescale_catchment_characteristics(vars, d$lookup_table, NULL,
                                                     d$catchment_characteristic),
               "refactored_areas required when")

  mini <- d$lookup_table[1:2,]

  rescale <- rescale_catchment_characteristics(vars, mini, d$split_divides,
                                               d$catchment_characteristic)

  expect_equal(nrow(rescale), 1)

  options(old_opts)

})