Package 'nhdplusTools'

Description

Given a river network with required base attributes, adds the NHDPlus network attributes: hydrosequence, levelpath, terminalpath, pathlength, down levelpath, down hydroseq, total drainage area, and terminalflag. The function implements two parallelization schemes for small and large basins respectively. If a number of cores is specified, parallel execution will be used.

Usage

add_plus_network_attributes(
  net,
  override = 5,
  cores = NULL,
  split_temp = NULL,
  status = TRUE
)

Arguments

Value

data.frame with added attributes

Examples

source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))

test_flowline <- prepare_nhdplus(walker_flowline, 0, 0, FALSE)

test_flowline <- data.frame(
  comid = test_flowline$COMID,
  tocomid = test_flowline$toCOMID,
  nameID = walker_flowline$GNIS_ID,
  lengthkm = test_flowline$LENGTHKM,
  areasqkm = walker_flowline$AreaSqKM)

add_plus_network_attributes(test_flowline)

Description

this function takes any NHDPlus dataset and aligns the attribute names with those used in nhdplusTools.

Usage

align_nhdplus_names(x)

Arguments

Value

data.frame renamed sf object

Examples

source(system.file("extdata/new_hope_data.R", package = "nhdplusTools"))

names(new_hope_flowline)

names(new_hope_flowline) <- tolower(names(new_hope_flowline))

new_hope_flowline <- align_nhdplus_names(new_hope_flowline)

names(new_hope_flowline)

Description

Calculates arbolate sum given a dendritic network and incremental lengths. Arbolate sum is the total length of all upstream flowlines.

Usage

calculate_arbolate_sum(x)

Arguments

Value

numeric with arbolate sum.

Examples

library(dplyr)
source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))
catchment_length <- select(walker_flowline, COMID, AreaSqKM) %>%
  right_join(prepare_nhdplus(walker_flowline, 0, 0,
                            purge_non_dendritic = FALSE, warn = FALSE), by = "COMID") %>%
  select(ID = COMID, toID = toCOMID, length = LENGTHKM)

arb_sum <- calculate_arbolate_sum(catchment_length)

catchment_length$arb_sum <- arb_sum
catchment_length$nhd_arb_sum <- walker_flowline$ArbolateSu

mean(abs(catchment_length$arb_sum - catchment_length$nhd_arb_sum))
max(abs(catchment_length$arb_sum - catchment_length$nhd_arb_sum))

Description

Calculates total drainage area given a dendritic network and incremental areas.

Usage

calculate_total_drainage_area(x)

Arguments

Value

numeric with total area.

Examples

library(dplyr)
source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))
catchment_area <- select(walker_flowline, COMID, AreaSqKM) %>%
  right_join(prepare_nhdplus(walker_flowline, 0, 0,
                            purge_non_dendritic = FALSE, warn = FALSE), by = "COMID") %>%
  select(ID = COMID, toID = toCOMID, area = AreaSqKM)

new_da <- calculate_total_drainage_area(catchment_area)

catchment_area$totda <- new_da
catchment_area$nhdptotda <- walker_flowline$TotDASqKM

mean(abs(catchment_area$totda - catchment_area$nhdptotda))
max(abs(catchment_area$totda - catchment_area$nhdptotda))

Description

Given a set of flowline indexes and numeric or ascii criteria, return closest match. If numeric criteria are used, the minimum difference in the numeric attribute is used for disambiguation. If ascii criteria are used, the adist function is used with the following algorithm: '1 - adist_score / max_string_length'. Comparisons ignore case.

Usage

disambiguate_flowline_indexes(indexes, flowpath, hydro_location)

Arguments

Value

data.frame indexes deduplicated according to the minimum difference between the values in the metric columns. If two or more result in the same "minimum" value, duplicates will be returned.

Examples

source(system.file("extdata", "sample_flines.R", package = "nhdplusTools"))

hydro_location <- sf::st_sf(id = c(1, 2, 3),
                            geom = sf::st_sfc(list(sf::st_point(c(-76.86934, 39.49328)),
                                                   sf::st_point(c(-76.91711, 39.40884)),
                                                   sf::st_point(c(-76.88081, 39.36354))),
                                              crs = 4326),
                            totda = c(23.6, 7.3, 427.9),
                            nameid = c("Patapsco", "", "Falls Run River"))

flowpath <- dplyr::select(sample_flines,
                          comid = COMID,
                          totda = TotDASqKM,
                          nameid = GNIS_NAME,
                          REACHCODE,
                          ToMeas,
                          FromMeas)

indexes <- get_flowline_index(flowpath,
                              hydro_location,
                              search_radius = 0.2,
                              max_matches = 10)

disambiguate_flowline_indexes(indexes,
                              dplyr::select(flowpath, comid, totda),
                              dplyr::select(hydro_location, id, totda))

result <- disambiguate_flowline_indexes(indexes,
                                        dplyr::select(flowpath, comid, nameid),
                                        dplyr::select(hydro_location, id, nameid))

result[result$id == 1, ]

result[result$id == 2, ]

result[result$id == 3, ]

Description

Queries the geoconnex.us reference feature server for available layers and attributes.

Usage

discover_geoconnex_reference()

Value

data.frame containing layers available and fields that are available to query.

Examples

discover_geoconnex_reference()

Description

Multipurpose function to find a COMID of interest.

Usage

discover_nhdplus_id(point = NULL, nldi_feature = NULL, raindrop = FALSE)

Arguments

Value

integer COMID or list containing COMID and raindrop trace.

Examples

point <- sf::st_sfc(sf::st_point(c(-76.874, 39.482)), crs = 4326)
discover_nhdplus_id(point)

discover_nhdplus_id(point, raindrop = TRUE)

nldi_nwis <- list(featureSource = "nwissite", featureID = "USGS-08279500")
discover_nhdplus_id(nldi_feature = nldi_nwis)

Description

This functionality is deprecated and will be removed in the near future.

Please use get_characteristics_metadata instead.

Usage

discover_nldi_characteristics(type = "all")

Arguments

Value

data.frame containing available characteristics

Description

Download NHD

Usage

download_nhd(nhd_dir, hu_list, download_files = TRUE)

Arguments

Value

character Paths to geodatabases created.

Examples

hu <- get_huc(sf::st_sfc(sf::st_point(c(-73, 42)), crs = 4326),
                            type = "huc08")

(hu <- substr(hu$huc8, 1, 2))

download_nhd(tempdir(), c(hu, "0203"), download_files = FALSE)

Description

Download NHDPlus HiRes

Usage

download_nhdplushr(nhd_dir, hu_list, download_files = TRUE, archive = FALSE)

Arguments

Value

character Paths to geodatabases created.

Examples

hu <- get_huc(sf::st_sfc(sf::st_point(c(-73, 42)), crs = 4326),
                            type = "huc08")

(hu <- substr(hu$huc8, 1, 2))

download_nhdplushr(tempdir(), c(hu, "0203"), download_files = FALSE)

download_nhdplushr(tempdir(), c(hu, "0203"), download_files = FALSE, archive = TRUE)

Description

This function downloads and decompresses staged seamless NHDPlusV2 data. The following requirements are needed: p7zip (MacOS), 7zip (windows) Please see: https://www.epa.gov/waterdata/get-nhdplus-national-hydrography-dataset-plus-data for more information and metadata about this data.

Default downloads lower-48 only. See examples for islands. No Alaska data are available.

Usage

download_nhdplusv2(
  outdir,
  url = paste0("https://dmap-data-commons-ow.s3.amazonaws.com/NHDPlusV21/",
    "Data/NationalData/NHDPlusV21_NationalData_Seamless", "_Geodatabase_Lower48_07.7z"),
  progress = TRUE
)

Arguments

Value

character path to the local geodatabase

Examples

## Not run: 
  download_nhdplusV2("./data/nhd/")

  download_nhdplusv2(outdir = "./inst/",
      url = paste0("https://dmap-data-commons-ow.s3.amazonaws.com/NHDPlusV21/",
                   "Data/NationalData/NHDPlusV21_NationalData_Seamless",
                   "_Geodatabase_HI_PR_VI_PI_03.7z"))

## End(Not run)

Description

This function downloads and decompresses staged RF1 data. See: https://water.usgs.gov/GIS/metadata/usgswrd/XML/erf1_2.xml for metadata.

Usage

download_rf1(
  outdir,
  url = "https://water.usgs.gov/GIS/dsdl/erf1_2.e00.gz",
  progress = TRUE
)

Arguments

Value

character path to the local e00 file

Examples

## Not run: 
  download_wbd("./data/rf1/")

## End(Not run)

Description

downloads and caches NHDPlusVAA data on your computer

Usage

download_vaa(
  path = get_vaa_path(updated_network),
  force = FALSE,
  updated_network = FALSE
)

Arguments

Details

The VAA data is a aggregate table of information from the NHDPlusV2 elevslope.dbf(s), PlusFlowlineVAA.dbf(s); and NHDFlowlines. All data originates from the EPA NHDPlus Homepage here. To see the location of cached data on your machine use get_vaa_path. To view aggregate data and documentation, see here

Value

character path to cached data

Description

This function downloads and decompresses staged seamless WBD data. Please see: https://prd-tnm.s3.amazonaws.com/StagedProducts/Hydrography/WBD/National/GDB/WBD_National_GDB.xml for metadata.

Usage

download_wbd(
  outdir,
  url = paste0("https://prd-tnm.s3.amazonaws.com/StagedProducts/",
    "Hydrography/WBD/National/GDB/WBD_National_GDB.zip"),
  progress = TRUE
)

Arguments

Value

character path to the local geodatabase

Examples

## Not run: 
  download_wbd("./data/wbd/")

## End(Not run)

Description

Calls the 3DHP_all web service and returns sf data.frames for the selected layers. See https://hydro.nationalmap.gov/arcgis/rest/services/3DHP_all/MapServer for source data documentation.

Usage

get_3dhp(
  AOI = NULL,
  ids = NULL,
  type = NULL,
  universalreferenceid = NULL,
  t_srs = NULL,
  buffer = 0.5,
  page_size = 2000
)

Arguments

Details

The returned object(s) will have the same Spatial Reference System (SRS) as the input AOI. If a individual or set of IDs are used to query, then the default CRS of EPSG:4269 is preserved. In all cases, a user-defined SRS can be passed to t_srs which will override all previous SRS (either input or default). All buffer and distance operations are handled internally using in EPSG:5070 Albers Equal Area projection

Value

a simple features (sf) object or valid types if no type supplied

Examples

AOI <- sf::st_as_sfc(sf::st_bbox(c(xmin = -89.56684, ymin = 42.99816,
                                   xmax = -89.24681, ymax = 43.17192),
                                 crs = "+proj=longlat +datum=WGS84 +no_defs"))

# get flowlines and hydrolocations
flowlines <- get_3dhp(AOI = AOI, type = "flowline")
hydrolocation <- get_3dhp(AOI = AOI, type = "hydrolocation")
waterbody <- get_3dhp(AOI = AOI, type = "waterbody")

if(!is.null(waterbody) & !is.null(flowlines) & !is.null(hydrolocation)) {
plot(sf::st_geometry(waterbody), col = "lightblue", border = "lightgrey")
plot(sf::st_geometry(flowlines), col = "blue", add = TRUE)
plot(sf::st_geometry(hydrolocation), col = "grey", pch = "+", add = TRUE) }

# given mainstem ids from any source, can query for them in ids.

CO <- get_3dhp(ids = "https://geoconnex.us/ref/mainstems/29559",
               type = "flowline")

if(!is.null(CO))
  plot(sf::st_geometry(CO), col = "blue")

# get all the waterbodies along the CO river
CO_wb <- get_3dhp(ids = unique(CO$waterbodyid3dhp), type = "waterbody")

if(!is.null(CO_wb)) {
plot(sf::st_geometry(CO_wb[grepl("Powell", CO_wb$gnisidlabel),]),
     col = "blue", border = "NA") }

# given universalreferenceid (reachcodes), can query for them but only
# for hydrolocations. This is useful for looking up mainstem ids.

get_3dhp(universalreferenceid = unique(hydrolocation$universalreferenceid),
         type = "hydrolocation")

Description

Return RPU or VPU boundaries

Usage

get_boundaries(type = "vpu")

Arguments

Value

An object of class "sf"

Description

Downloads (subsets of) catchment characteristics from a cloud data store. See get_characteristics_metadata for available characteristics.

Source: Wieczorek, M.E., Jackson, S.E., and Schwarz, G.E., 2018, Select Attributes for NHDPlus Version 2.1 Reach Catchments and Modified Network Routed Upstream Watersheds for the Conterminous United States (ver. 3.0, January 2021): U.S. Geological Survey data release, doi:10.5066/F7765D7V.

Usage

get_catchment_characteristics(varname, ids, reference_fabric = "nhdplusv2")

Arguments

Examples

get_catchment_characteristics("CAT_BFI", c(5329343, 5329427))

Description

Download and cache table of catchment characteristics.

Wieczorek, M.E., Jackson, S.E., and Schwarz, G.E., 2018, Select Attributes for NHDPlus Version 2.1 Reach Catchments and Modified Network Routed Upstream Watersheds for the Conterminous United States (ver. 3.0, January 2021): U.S. Geological Survey data release, doi:10.5066/F7765D7V.

Usage

get_characteristics_metadata(search, cache = TRUE)

Arguments

Examples

get_characteristics_metadata()

Description

Traverse NHDPlus network downstream with diversions NOTE: This algorithm may not scale well in large watersheds. For reference, the lower Mississippi will take over a minute.

Usage

get_DD(network, comid, distance = NULL)

Arguments

Value

integer vector of all COMIDs downstream of the starting COMID

Examples

library(sf)
start_COMID <- 11688818

source(system.file("extdata", "sample_flines.R", package = "nhdplusTools"))

DD_COMIDs <- get_DD(sample_flines, start_COMID, distance = 4)
plot(dplyr::filter(sample_flines, COMID %in% DD_COMIDs)$geom,
     col = "red", lwd = 2)

DM_COMIDs <- get_DM(sample_flines, start_COMID, distance = 4)
plot(dplyr::filter(sample_flines, COMID %in% DM_COMIDs)$geom,
     col = "blue", add = TRUE, lwd = 2)

Description

Traverse NHDPlus network downstream main stem

Usage

get_DM(network, comid, distance = NULL, sort = FALSE, include = TRUE)

Arguments

Value

integer vector of all COMIDs downstream of the starting COMID along the mainstem

Examples

library(sf)

source(system.file("extdata", "sample_flines.R", package = "nhdplusTools"))

plot(sample_flines$geom)
start_COMID <- 11690092
DM_COMIDs <- get_DM(sample_flines, start_COMID)
plot(dplyr::filter(sample_flines, COMID %in% DM_COMIDs)$geom,
     col = "red", add = TRUE, lwd = 3)

DM_COMIDs <- get_DM(sample_flines, start_COMID, distance = 40)
plot(dplyr::filter(sample_flines, COMID %in% DM_COMIDs)$geom,
     col = "blue", add = TRUE, lwd = 2)

Description

Uses a cross section retrieval web services to retrieve elevation along a path.

Usage

get_elev_along_path(points, num_pts, res = 1, status = TRUE)

Arguments

Value

sf data.frame containing points retrieved. Names include "id", "distance_m", "elevation_m", "spatial_ref", "geometry", and ".group". .group tracks which input point each set of output points belongs to.

Examples

point1 <- sf::st_sfc(sf::st_point(x = c(-105.9667, 36.17602)), crs = 4326)
point2 <- sf::st_sfc(sf::st_point(x = c(-105.97768, 36.17526)), crs = 4326)
point3 <- sf::st_sfc(sf::st_point(x = c(-105.98869, 36.17450)), crs = 4326)

points <- sf::st_as_sf(c(point1, point2, point3))

(xs <- get_elev_along_path(points, 100))

if(inherits(xs, "sf")) {

bbox <- sf::st_bbox(xs) + c(-0.005, -0.005, 0.005, 0.005)

nhdplusTools::plot_nhdplus(bbox = bbox, cache_data = FALSE)

plot(sf::st_transform(sf::st_geometry(xs), 3857), pch = ".", add = TRUE, col = "red")
plot(sf::st_transform(sf::st_sfc(point1, crs = 4326), 3857), add = TRUE)
plot(sf::st_transform(sf::st_sfc(point2, crs = 4326), 3857), add = TRUE)
plot(sf::st_transform(sf::st_sfc(point3, crs = 4326), 3857), add = TRUE)

plot(xs$distance_m, xs$elevation_m)
}

Description

given an sf point geometry column, return COMID, reachcode, and measure for each.

Usage

get_flowline_index(
  flines,
  points,
  search_radius = NULL,
  precision = NA,
  max_matches = 1
)

Arguments

Details

Note 1: Inputs are cast into LINESTRINGS. Because of this, the measure output of inputs that are true multipart lines may be in error.

Note 2: This algorithm finds the nearest node in the input flowlines to identify which flowline the point should belong to. As a second pass, it can calculate the measure to greater precision than the nearest flowline geometry node.

Note 3: Offset is returned in units consistent with the projection of the input points.

Note 4: See 'dfMaxLength' input to sf::st_segmentize() for details of handling of precision parameter.

Note 5: "from" is downstream – 0 is the outlet "to" is upstream – 100 is the inlet

Value

data.frame with five columns, id, COMID, REACHCODE, REACH_meas, and offset. id is the row or list element in the point input.

Examples

source(system.file("extdata", "sample_flines.R", package = "nhdplusTools"))

point <- sf::st_sfc(sf::st_point(c(-76.87479, 39.48233)),
                    crs = 4326)

get_flowline_index(sample_flines, point)

point <- sf::st_transform(point, 5070)

get_flowline_index(sample_flines, point,
                   search_radius = units::set_units(200, "m"))

get_flowline_index("download_nhdplusv2", point)

get_flowline_index(sample_flines, point, precision = 30)

get_flowline_index(sample_flines,
                   sf::st_sfc(list(sf::st_point(c(-76.86934, 39.49328)),
                                   sf::st_point(c(-76.91711, 39.40884)),
                                   sf::st_point(c(-76.88081, 39.36354))),
                              crs = 4326),
                   search_radius = units::set_units(0.2, "degrees"),
                   max_matches = 10)

Description

Subsets the gagesII dataset by location (POINT), area (POLYGON), or set of IDs. See <doi:10.5066/P96CPHOT> for documentation of source data.

Usage

get_gagesII(AOI = NULL, id = NULL, t_srs = NULL, buffer = 0.5, basin = FALSE)

Arguments

Details

The returned object(s) will have the same Spatial Reference System (SRS) as the input AOI. If a individual or set of IDs are used to query, then the default geoserver CRS of EPSG:4326 is preserved. In all cases, a user-defined SRS can be passed to t_srs which will override all previous SRS (either input or default). All buffer and distance operations are handled internally using in EPSG:5070 Albers Equal Area projection

Value

a simple features (sf) object

Description

Queries the geoconnex reference feature server for features of interest.

Usage

get_geoconnex_reference(
  AOI,
  type = NULL,
  t_srs = NULL,
  buffer = 0.5,
  status = TRUE
)

Arguments

Value

sf data.frame containing requested reference features

Examples

dplyr::distinct(discover_geoconnex_reference()[c("id", "title")])

AOI <- sf::st_as_sfc(sf::st_bbox(c(xmin = -89.56684, ymin = 42.99816,
                                   xmax = -89.24681, ymax = 43.17192),
                                 crs = "+proj=longlat +datum=WGS84 +no_defs"))

get_geoconnex_reference(AOI, type = "hu04")

get_geoconnex_reference("https://geoconnex.us/ref/mainstems/315626", type = "hu04", )

AOI <- sf::st_sfc(sf::st_point(c(-89.56684, 42.99816)),
                  crs = "+proj=longlat +datum=WGS84 +no_defs")

get_geoconnex_reference(AOI, type = "hu04", buffer = 100000)

Description

Use to remove unwanted detail NHDPlusHR data See get_nhdplushr for examples.

Usage

get_hr_data(
  gdb,
  layer = NULL,
  min_size_sqkm = NULL,
  simp = NULL,
  proj = NULL,
  rename = TRUE
)

Arguments

Value

sf data.frame containing requested data

Description

Subsets WBD features by location (POINT), area (POLYGON), or set of HUC IDs.

Usage

get_huc(AOI = NULL, id = NULL, t_srs = NULL, buffer = 0.5, type = "huc12")

Arguments

Details

The returned object(s) will have the same Spatial Reference System (SRS) as the input AOI. If a individual or set of IDs are used to query, then the default geoserver CRS of EPSG:4326 is preserved. In all cases, a user-defined SRS can be passed to t_srs which will override all previous SRS (either input or default). All buffer and distance operations are handled internally using in EPSG:5070 Albers Equal Area projection

Value

a simple features (sf) object

Description

Calculates level paths using the stream-leveling approach of NHD and NHDPlus. In addition to a levelpath identifier, a topological sort and levelpath outlet identifier is provided in output. If arbolate sum is provided in the weight column, this will match the behavior of NHDPlus. Any numeric value can be included in this column and the largest value will be followed when no nameID is available.

Usage

get_levelpaths(x, override_factor = NULL, status = FALSE, cores = NULL)

Arguments

Details

1. levelpath provides an identifier for the collection of flowlines that make up the single mainstem flowpath of a total upstream aggregate catchment.

2. outletID is the catchment ID (COMID in the case of NHDPlus) for the catchment at the outlet of the levelpath the catchment is part of.

3. topo_sort is similar to Hydroseq in NHDPlus in that large topo_sort values are upstream of small topo_sort values. Note that there are many valid topological sort orders of a directed graph.

Value

data.frame with ID, outletID, topo_sort, and levelpath columns. See details for more info.

Examples

source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))

test_flowline <- prepare_nhdplus(walker_flowline, 0, 0, FALSE)

test_flowline <- data.frame(
  ID = test_flowline$COMID,
  toID = test_flowline$toCOMID,
  nameID = walker_flowline$GNIS_ID,
  weight = walker_flowline$ArbolateSu,
  stringsAsFactors = FALSE)

get_levelpaths(test_flowline)

Description

Subsets NHDPlusV2 Area features by location (POINT), area (POLYGON), or set of IDs. See download_nhdplusv2 for source data documentation.

Usage

get_nhdarea(AOI = NULL, id = NULL, t_srs = NULL, buffer = 0.5)

Arguments

Details

The returned object(s) will have the same Spatial Reference System (SRS) as the input AOI. If a individual or set of IDs are used to query, then the default geoserver CRS of EPSG:4326 is preserved. In all cases, a user-defined SRS can be passed to t_srs which will override all previous SRS (either input or default). All buffer and distance operations are handled internally using in EPSG:5070 Albers Equal Area projection

Value

a simple features (sf) object

Description

Subsets NHDPlusV2 features by location (POINT), area (POLYGON), or set of COMIDs. Multi realizations are supported allowing you to query for flowlines, catchments, or outlets.

Usage

get_nhdplus(
  AOI = NULL,
  comid = NULL,
  nwis = NULL,
  realization = "flowline",
  streamorder = NULL,
  t_srs = NULL
)

Arguments

Details

The returned object(s) will have the same Spatial Reference System (SRS) as the input AOI. If a individual or set of IDs are used to query, then the default geoserver CRS of EPSG:4326 is preserved. In all cases, a user-defined SRS can be passed to t_srs which will override all previous SRS (either input or default). All buffer and distance operations are handled internally using in EPSG:5070 Albers Equal Area projection

Value

sfc a single, or list, of simple feature objects

Examples

point <- sf::st_sfc(sf::st_point(c(-119.845, 34.4146)), crs = 4326)
 get_nhdplus(point)
 get_nhdplus(point, realization = "catchment")
 get_nhdplus(point, realization = "all")
 get_nhdplus(comid = 101)
 get_nhdplus(nwis  = c(11120000, 11120500))
 area <- sf::st_as_sfc(sf::st_bbox(c(xmin = -119.8851, xmax =-119.8361,
 ymax = 34.42439, ymin = 34.40473), crs = 4326))
 get_nhdplus(area)
 get_nhdplus(area, realization = "flowline", streamorder = 3)

Description

Get NHDPlus HiRes

Usage

get_nhdplushr(
  hr_dir,
  out_gpkg = NULL,
  layers = c("NHDFlowline", "NHDPlusCatchment"),
  pattern = ".*GDB.gdb$",
  check_terminals = TRUE,
  overwrite = FALSE,
  keep_cols = NULL,
  ...
)

Arguments

Details

NHDFlowline is joined to value added attributes prior to being returned. Names are not modified from the NHDPlusHR geodatabase. Set layers to "NULL" to get all layers.

Value

sf data.frames containing output that may also be written to a geopackage for later use.

Examples

## Not run: 
# Note this will download a lot of data to a temp directory.
# Change 'temp_dir' to your directory of choice.
temp_dir <- file.path(nhdplusTools_data_dir(), "temp_hr_cache")

download_dir <- download_nhdplushr(temp_dir, c("0302", "0303"))

get_nhdplushr(download_dir, file.path(download_dir, "nhdplus_0302-03.gpkg"))

get_nhdplushr(download_dir,
              file.path(download_dir, "nhdplus_0302-03.gpkg"),
              layers = NULL, overwrite = TRUE)

get_nhdplushr(download_dir,
              file.path(download_dir, "nhdplus_0302-03.gpkg"),
              layers = "NHDFlowline", overwrite = TRUE,
              min_size_sqkm = 10, simp = 10, proj = "+init=epsg:5070")

# Cleanup
unlink(temp_dir, recursive = TRUE)


## End(Not run)

Description

Get a basin boundary for a given NLDI feature.

Usage

get_nldi_basin(nldi_feature, simplify = TRUE, split = FALSE)

Arguments

Details

Only resolves to the nearest NHDPlus catchment divide. See: https://waterdata.usgs.gov/blog/nldi-intro/ for more info on the nldi.

Value

sf data.frame with result basin boundary

Examples

library(sf)
library(dplyr)

nldi_nwis <- list(featureSource = "nwissite", featureID = "USGS-05428500")

site <- get_nldi_feature(nldi_nwis)

basin <- get_nldi_basin(nldi_feature = nldi_nwis)

plot(st_geometry(basin))

basin

basin2 <- get_nldi_basin(nldi_feature = nldi_nwis,
                         simplify = FALSE, split = TRUE)

if(inherits(basin, "sf") & inherits(basin2, "sf")) {

length(st_coordinates(basin))
length(st_coordinates(basin2))

plot(st_geometry(st_buffer(site, units::set_units(3000, "m"))), border = NA)

plot(st_geometry(site), add = TRUE)
plot(st_geometry(basin2), add = TRUE)

plot(st_geometry(basin), border = "red", add = TRUE)

}

Description

This functionality id deprecated and will be removed in the near future.

Please use get_catchment_characteristics instead.

Usage

get_nldi_characteristics(nldi_feature, type = "local")

Arguments

Value

data.frame containing requested characteristics

Description

Get a single feature from the NLDI

Usage

get_nldi_feature(nldi_feature)

Arguments

Value

sf data.frame with one feature

Examples

get_nldi_feature(list("featureSource" = "nwissite", featureID = "USGS-05428500"))

Description

uses the Network Linked Data Index to retrieve and estimated network location for the given point. If not within a grid cell of a flowline, will use a raindrop trace service to find the nearest downslope flowline location.

Usage

get_nldi_index(location)

Arguments

Examples

index <- get_nldi_index(c(-89.276, 42.988))

if(inherits(index, "sf")) {

plot_nhdplus(bbox = sf::st_bbox(sf::st_buffer(index[1,], units::set_units(1000, "m"))))
plot(sf::st_geometry(sf::st_transform(index, 3857)), add = TRUE)

}

Description

Returns a POINT feature class of active, stream network, NWIS gages for an Area of Interest. If a POINT feature is used as an AOI, then the returned sites within the requested buffer, are sorted by distance (in meters) from that POINT.

Usage

get_nwis(AOI = NULL, t_srs = NULL, buffer = 20000)

Arguments

Details

The returned object(s) will have the same Spatial Reference System (SRS) as the input AOI. If a individual or set of IDs are used to query, then the default geoserver CRS of EPSG:4326 is preserved. In all cases, a user-defined SRS can be passed to t_srs which will override all previous SRS (either input or default). All buffer and distance operations are handled internally using in EPSG:5070 Albers Equal Area projection

Value

a simple features (sf) object

Description

Given a network and set of IDs, finds path lengths between all identified flowpath outlets. This algorithm finds distance between outlets regardless of flow direction.

Usage

get_path_lengths(outlets, network, cores = 1, status = FALSE)

Arguments

Value

data.frame containing the distance between pairs of network outlets. For a network with one terminal outlet, the data.frame will have 'nrow(network)^2' rows.

Examples

source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))
fline <- walker_flowline

outlets <- c(5329303, 5329357, 5329317, 5329365, 5329435, 5329817)

# Add toCOMID
fline <- nhdplusTools::get_tocomid(fline, add = TRUE)

fl <- dplyr::select(fline, ID = comid, toID = tocomid, lengthkm)

get_path_lengths(outlets, fl)

Description

Given a network and set of IDs, finds paths between all identified flowpath outlets. This algorithm finds members between outlets regardless of flow direction.

Usage

get_path_members(outlets, network, cores = 1, status = FALSE)

Arguments

Value

list of lists containing flowpath identifiers along path that connect outlets.

Examples

source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))
fline <- walker_flowline

outlets <- c(5329303, 5329357, 5329317, 5329365, 5329435, 5329817)

# Add toCOMID
fline <- nhdplusTools::get_tocomid(fline, add = TRUE)

fl <- dplyr::select(fline, ID = comid, toID = tocomid, lengthkm)

get_path_members(outlets, fl)

Description

Generates the main path length to a basin's terminal path.

Usage

get_pathlength(x)

Arguments

Value

data.frame containing pathlength for each ID

Examples

source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))

fl <- dplyr::select(prepare_nhdplus(walker_flowline, 0, 0),
                    ID = COMID, toID = toCOMID, length = LENGTHKM)

get_pathlength(fl)

Description

Determines Pfafstetter codes for a dendritic network with total drainage area, levelpath, and topo_sort attributes.

Usage

get_pfaf(x, max_level = 2, status = FALSE)

Arguments

Value

data.frame with ID and pfaf columns.

Examples

library(dplyr)
source(system.file("extdata/nhdplushr_data.R", package = "nhdplusTools"))
hr_flowline <- align_nhdplus_names(hr_data$NHDFlowline)

fl <-  select(hr_flowline, COMID, AreaSqKM) %>%
  right_join(prepare_nhdplus(hr_flowline, 0, 0,
                             purge_non_dendritic = FALSE,
                             warn = FALSE),
             by = "COMID") %>%
  sf::st_sf() %>%
  select(ID = COMID, toID = toCOMID, area = AreaSqKM)

fl$nameID = ""
fl$totda <- calculate_total_drainage_area(sf::st_set_geometry(fl, NULL))
fl <- left_join(fl, get_levelpaths(rename(sf::st_set_geometry(fl, NULL),
                                   weight = totda)), by = "ID")

pfaf <- get_pfaf(fl, max_level = 3)

fl <- left_join(fl, pfaf, by = "ID")

plot(fl["pf_level_3"], lwd = 2)

pfaf <- get_pfaf(fl, max_level = 4)

hr_catchment <- left_join(hr_data$NHDPlusCatchment, pfaf, by = c("FEATUREID" = "ID"))

colors <- data.frame(pf_level_4 = unique(hr_catchment$pf_level_4),
                     color = sample(terrain.colors(length(unique(hr_catchment$pf_level_4)))),
                    stringsAsFactors = FALSE)
hr_catchment <- left_join(hr_catchment, colors, by = "pf_level_4")
plot(hr_catchment["color"], border = NA, reset = FALSE)
plot(sf::st_geometry(hr_flowline), col = "blue", add = TRUE)

source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))

fl <- select(walker_flowline, COMID, AreaSqKM) %>%
  right_join(prepare_nhdplus(walker_flowline, 0, 0,
                            purge_non_dendritic = FALSE, warn = FALSE),
            by = "COMID") %>%
  sf::st_sf() %>%
  select(ID = COMID, toID = toCOMID, area = AreaSqKM)

fl$nameID = ""
fl$totda <- calculate_total_drainage_area(sf::st_set_geometry(fl, NULL))
fl <- left_join(fl, get_levelpaths(rename(sf::st_set_geometry(fl, NULL),
                                   weight = totda)), by = "ID")

pfaf <- get_pfaf(fl, max_level = 2)

fl <- left_join(fl, pfaf, by = "ID")

plot(fl["pf_level_2"], lwd = 2)

Description

Uses a raindrop trace web service to trace the nhdplus digital elevation model to the nearest downslope flowline.

Usage

get_raindrop_trace(point, direction = "down")

Arguments

Value

sf data.frame containing raindrop trace and requested portion of flowline.

Examples

point <- sf::st_sfc(sf::st_point(x = c(-89.2158, 42.9561)), crs = 4326)

(trace <- get_raindrop_trace(point))

if(inherits(trace, "sf")) {
bbox <- sf::st_bbox(trace) + c(-0.005, -0.005, 0.005, 0.005)

nhdplusTools::plot_nhdplus(bbox = bbox, cache_data = FALSE)

plot(sf::st_transform(sf::st_sfc(point, crs = 4326), 3857), add = TRUE)
plot(sf::st_transform(sf::st_geometry(trace)[1], 3857), add = TRUE, col = "red")
plot(sf::st_transform(sf::st_geometry(trace)[2], 3857), add = TRUE, col = "black")
}

Description

given a tree with an id and and toid in the first and second columns, returns a sorted and potentially split set of output.

Can also be used as a very fast implementation of upstream with tributaries navigation. The full network from each outlet is returned in sorted order.

Usage

get_sorted(x, split = FALSE, outlets = NULL)

Arguments

Value

data.frame containing a topologically sorted version of the requested network and optionally a terminal id.

Examples

source(system.file("extdata/new_hope_data.R", package = "nhdplusTools"))

fpath <- get_tocomid(
  dplyr::select(new_hope_flowline, COMID, FromNode, ToNode, Divergence, FTYPE,
                AreaSqKM, LENGTHKM, GNIS_ID)
)

head(fpath <- get_sorted(fpath, split = TRUE))

fpath['sort_order'] <- 1:nrow(fpath)

plot(fpath['sort_order'])

Description

Uses a catchment splitting web service to retrieve the portion of a catchment upstream of the point provided.

Usage

get_split_catchment(point, upstream = TRUE)

Arguments

Details

This service works within the coterminous US NHDPlusV2 domain. If the point provided falls on an NHDPlusV2 flowline as retrieved from get_raindrop_trace the catchment will be split across the flow line. IF the point is not along the flowline a small sub catchment will typically result. As a result, most users of this function will want to use get_raindrop_trace prior to calls to this function.

An attempt is made to eliminate polygon shards if they exist in the output. However, there is a chance that this function will return a multipolygon data.frame.

Value

sf data.frame containing the local catchment, the split portion and optionally the total drainage basin.

Examples

point <- sf::st_sfc(sf::st_point(x = c(-89.2158, 42.9561)), crs = 4326)

trace <- get_raindrop_trace(point)

if(inherits(trace, "sf")) {

(snap_point <- sf::st_sfc(sf::st_point(trace$intersection_point[[1]]),
                          crs = 4326))

(catchment <- get_split_catchment(snap_point))

bbox <- sf::st_bbox(catchment) + c(-0.005, -0.005, 0.005, 0.005)

nhdplusTools::plot_nhdplus(bbox = bbox, cache_data = FALSE)

plot(sf::st_transform(sf::st_geometry(catchment)[2], 3857), add = TRUE, col = "black")
plot(sf::st_transform(sf::st_geometry(catchment)[1], 3857), add = TRUE, col = "red")
plot(sf::st_transform(sf::st_sfc(point, crs = 4326), 3857), add = TRUE, col = "white")

(catchment <- get_split_catchment(snap_point, upstream = FALSE))

bbox <- sf::st_bbox(catchment) + c(-0.005, -0.005, 0.005, 0.005)

nhdplusTools::plot_nhdplus(bbox = bbox, cache_data = FALSE)

plot(sf::st_transform(sf::st_geometry(catchment)[1], 3857), add = TRUE, col = "red")
plot(sf::st_transform(sf::st_geometry(catchment)[2], 3857), add = TRUE, col = "black")
plot(sf::st_transform(sf::st_sfc(point, crs = 4326), 3857), add = TRUE, col = "white")

pour_point <- sf::st_sfc(sf::st_point(x = c(-89.25619, 42.98646)), crs = 4326)

(catchment <- get_split_catchment(pour_point, upstream = FALSE))

bbox <- sf::st_bbox(catchment) + c(-0.005, -0.005, 0.005, 0.005)

nhdplusTools::plot_nhdplus(bbox = bbox, cache_data = FALSE)

plot(sf::st_transform(sf::st_geometry(catchment)[1], 3857), add = TRUE, col = "red")
plot(sf::st_transform(sf::st_geometry(catchment)[2], 3857), add = TRUE, col = "black")
plot(sf::st_transform(sf::st_sfc(pour_point, crs = 4326), 3857), add = TRUE, col = "white")
}

Description

Applies a topological sort and calculates stream level. Algorithm: Terminal level paths are assigned level 1 (see note 1). Paths that terminate at a level 1 are assigned level 2. This pattern is repeated until no paths remain.

If a TRUE/FALSE coastal attribute is included, coastal terminal paths begin at 1 and internal terminal paths begin at 4 as is implemented by the NHD stream leveling rules.

Usage

get_streamlevel(x)

Arguments

Value

numeric stream order in same order as input

Examples

source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))

test_flowline <- data.frame(
 levelpathi = walker_flowline$LevelPathI,
 dnlevelpat = walker_flowline$DnLevelPat)

 test_flowline$dnlevelpat[1] <- 0

(level <- get_streamlevel(test_flowline))

walker_flowline$level <- level

plot(sf::st_geometry(walker_flowline), lwd = walker_flowline$level, col = "blue")

test_flowline$coastal <- rep(FALSE, nrow(test_flowline))
(level <- get_streamlevel(test_flowline))

test_flowline$coastal[!test_flowline$dnlevelpat %in% test_flowline$levelpathi] <- TRUE
(level <- get_streamlevel(test_flowline))

Description

Applies a topological sort and calculates strahler stream order. Algorithm: If more than one upstream flowpath has an order equal to the maximum upstream order then the downstream flowpath is assigned the maximum upstream order plus one. Otherwise it is assigned the max upstream order.

Usage

get_streamorder(x, status = TRUE)

Arguments

Value

numeric stream order in same order as input

Examples

source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))

test_flowline <- prepare_nhdplus(walker_flowline, 0, 0, FALSE)

test_flowline <- data.frame(
  ID = test_flowline$COMID,
  toID = test_flowline$toCOMID)

(order <- get_streamorder(test_flowline))

walker_flowline$order <- order

plot(sf::st_geometry(walker_flowline), lwd = walker_flowline$order, col = "blue")

Description

Get the ID of the basin outlet for each flowline. This function has been deprecated in favor of get_sorted.

Usage

get_terminal(x, outlets)

Arguments

Value

data.frame containing the terminal ID for each outlet

Examples

source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))

fl <- dplyr::select(prepare_nhdplus(walker_flowline, 0, 0),
                    ID = COMID, toID = toCOMID)

outlet <- fl$ID[which(!fl$toID %in% fl$ID)]

get_terminal(fl, outlet)

Description

Given flowlines with fromnode and tonode attributes, will return a toid attribute that is the result of joining tonode and fromnode attributes. In the case that a terminalpa attribute is included, the join is executed by terminalpa group. This is done grouped by terminalpathID because duplicate node ids have been encountered across basins in some datasets. If 'remove_coastal' is 'TRUE' (the default) either ftype or fcode are required. Uses the add_toids function.

Usage

get_tocomid(
  x,
  return_dendritic = TRUE,
  missing = 0,
  remove_coastal = TRUE,
  add = TRUE
)

Arguments

Value

data.frame containing comid and tocomid attributes or all attributes provided with comid and tocomid in the first and second columns..

Examples

source(system.file("extdata", "sample_flines.R", package = "nhdplusTools"))

tocomid <- get_tocomid(sample_flines)

tocomid <- get_tocomid(sample_flines, return_dendritic = FALSE)

Description

Traverse NHDPlus network upstream main stem

Usage

get_UM(network, comid, distance = NULL, sort = FALSE, include = TRUE)

Arguments

Value

integer vector of all COMIDs upstream of the starting COMID along the mainstem

Examples

library(sf)

source(system.file("extdata", "sample_flines.R", package = "nhdplusTools"))

plot(sample_flines$geom)
start_COMID <- 11690196
UM_COMIDs <- get_UM(sample_flines, start_COMID)
plot(dplyr::filter(sample_flines, COMID %in% UM_COMIDs)$geom,
     col = "red", add = TRUE, lwd = 3)

UM_COMIDs <- get_UM(sample_flines, start_COMID, distance = 50)
plot(dplyr::filter(sample_flines, COMID %in% UM_COMIDs)$geom,
     col = "blue", add = TRUE, lwd = 2)

Description

Traverse NHDPlus network upstream with tributaries

Usage

get_UT(network, comid, distance = NULL)

Arguments

Value

integer vector of all COMIDs upstream with tributaries of the starting COMID.

Examples

library(sf)
source(system.file("extdata", "sample_flines.R", package = "nhdplusTools"))
plot(sample_flines$geom)
start_COMID <- 11690196
UT_COMIDs <- get_UT(sample_flines, start_COMID)
plot(dplyr::filter(sample_flines, COMID %in% UT_COMIDs)$geom,
     col = "red", add = TRUE)

UT_COMIDs <- get_UT(sample_flines, start_COMID, distance = 50)
plot(dplyr::filter(sample_flines, COMID %in% UT_COMIDs)$geom,
     col = "blue", add = TRUE)

Description

Return requested NHDPlusV2 Attributes.

Usage

get_vaa(
  atts = NULL,
  path = get_vaa_path(),
  download = TRUE,
  updated_network = FALSE
)

Arguments

Details

The VAA data is a aggregate table of information from the NHDPlusV2 elevslope.dbf(s), PlusFlowlineVAA.dbf(s); and NHDFlowlines. All data originates from the EPA NHDPlus Homepage here. To see the location of cached data on your machine use get_vaa_path. To view aggregate data and documentation, see here

Value

data.frame containing requested VAA data

Examples

## Not run: 
# This will download the vaa file to the path from get_vaa_path()

get_vaa("slope")
get_vaa(c("slope", "lengthkm"))

get_vaa(updated_network = TRUE)
get_vaa("reachcode", updated_network = TRUE)

#cleanup if desired
unlink(dirname(get_vaa_path()), recursive = TRUE)

## End(Not run)

Description

Find variables available from the NHDPlusV2 attribute data.frame

Usage

get_vaa_names(updated_network = FALSE)

Arguments

Details

The VAA data is a aggregate table of information from the NHDPlusV2 elevslope.dbf(s), PlusFlowlineVAA.dbf(s); and NHDFlowlines. All data originates from the EPA NHDPlus Homepage here. To see the location of cached data on your machine use get_vaa_path. To view aggregate data and documentation, see here

Value

character vector

Examples

## Not run: 
# This will download the vaa file to the path from get_vaa_path()
get_vaa_names()

#cleanup if desired
unlink(dirname(get_vaa_path()), recursive = TRUE)

## End(Not run)

Description

nhdplusTools will download and cache an 'fst' file with NHDPlusV2 attribute data sans geometry. This function returns the file path to the cached file. Will use the user data dir indicated by nhdplusTools_data_dir.

Usage

get_vaa_path(updated_network = FALSE)

Arguments

Details

The VAA data is a aggregate table of information from the NHDPlusV2 elevslope.dbf(s), PlusFlowlineVAA.dbf(s); and NHDFlowlines. All data originates from the EPA NHDPlus Homepage here. To see the location of cached data on your machine use get_vaa_path. To view aggregate data and documentation, see here

Value

character file path

Examples

get_vaa_path()

get_vaa_path(updated_network = TRUE)

Description

Subsets NHDPlusV2 waterbody features by location (POINT), area (POLYGON), or set of IDs. See download_nhdplusv2 for source data documentation.

Usage

get_waterbodies(AOI = NULL, id = NULL, t_srs = NULL, buffer = 0.5)

Arguments

Details

The returned object(s) will have the same Spatial Reference System (SRS) as the input AOI. If a individual or set of IDs are used to query, then the default geoserver CRS of EPSG:4326 is preserved. In all cases, a user-defined SRS can be passed to t_srs which will override all previous SRS (either input or default). All buffer and distance operations are handled internally using in EPSG:5070 Albers Equal Area projection

Value

a simple features (sf) object

Description

given an sf point geometry column, return waterbody id, and COMID of dominant artificial path

Usage

get_waterbody_index(waterbodies, points, flines = NULL, search_radius = NULL)

Arguments

Value

data.frame with two columns, COMID, in_wb_COMID, near_wb_COMID, near_wb_dist, and outlet_fline_COMID. Distance is in units of provided projection.

Examples

source(system.file("extdata/sample_data.R", package = "nhdplusTools"))

waterbodies <- sf::st_transform(
  sf::read_sf(sample_data, "NHDWaterbody"), 5070)

points <- sf::st_transform(
  sf::st_sfc(sf::st_point(c(-89.356086, 43.079943)),
             crs = 4326), 5070)

get_waterbody_index(waterbodies, points,
                    search_radius = units::set_units(500, "m"))

Description

Get Waterbody Outlet

Usage

get_wb_outlet(lake_id, network)

Arguments

Value

sf data.frame with single record of network COMID associated with most-downstream reach in the NHD Waterbody

Examples

source(system.file("extdata/sample_data.R", package = "nhdplusTools"))

fline <- sf::read_sf(sample_data, "NHDFlowline_Network")
wtbdy <- sf::read_sf(sample_data, "NHDWaterbody")

lake_COMID <- wtbdy$COMID[wtbdy$GNIS_NAME=='Lake Mendota 254']

get_wb_outlet(13293262, fline)

Description

Uses a cross section retrieval web services to retrieve a cross section given a point and specified width. Orientation is determined based on direction of a the flowline found near point. This function uses a 10m National Elevation Dataset request on the back end.

Usage

get_xs_point(point, width, num_pts)

Arguments

Value

sf data.frame containing points retrieved.

Examples

point <- sf::st_sfc(sf::st_point(x = c(-105.97218, 36.17592)), crs = 4326)

(xs <- get_xs_point(point, 300, 100))

if(inherits(xs, "sf")) {

bbox <- sf::st_bbox(xs) + c(-0.005, -0.005, 0.005, 0.005)

nhdplusTools::plot_nhdplus(bbox = bbox, cache_data = FALSE)

plot(sf::st_transform(sf::st_geometry(xs), 3857), pch = ".", add = TRUE, col = "red")
plot(sf::st_transform(sf::st_sfc(point, crs = 4326), 3857), add = TRUE)

plot(xs$distance_m, xs$elevation_m)
}

Description

Uses a cross section retrieval web services to retrieve a cross section between two endpoints.

Usage

get_xs_points(point1, point2, num_pts, res = 1)

Arguments

Value

sf data.frame containing points retrieved.

Examples

point1 <- sf::st_sfc(sf::st_point(x = c(-105.9667, 36.17602)), crs = 4326)
point2 <- sf::st_sfc(sf::st_point(x = c(-105.97768, 36.17526)), crs = 4326)

(xs <- get_xs_points(point1, point2, 100))

if(inherits(xs, "sf")) {

bbox <- sf::st_bbox(xs) + c(-0.005, -0.005, 0.005, 0.005)

nhdplusTools::plot_nhdplus(bbox = bbox, cache_data = FALSE)

plot(sf::st_transform(sf::st_geometry(xs), 3857), pch = ".", add = TRUE, col = "red")
plot(sf::st_transform(sf::st_sfc(point1, crs = 4326), 3857), add = TRUE)
plot(sf::st_transform(sf::st_sfc(point2, crs = 4326), 3857), add = TRUE)

plot(xs$distance_m, xs$elevation_m)
}

Description

creates a node topology table from an edge topology

Usage

make_node_topology(x, add_div = NULL, add = TRUE)

Arguments

Value

data.frame containing id, fromnode, and tonode attributes or all attributes provided with id, fromnode and tonode in the first three columns.

Examples

source(system.file("extdata/new_hope_data.R", package = "nhdplusTools"))

x <- dplyr::select(get_tocomid(
  dplyr::select(new_hope_flowline, COMID, FromNode, ToNode, Divergence, FTYPE,
                AreaSqKM, LENGTHKM, GNIS_ID)
), -tonode, -fromnode)

head(y <- make_node_topology(x))

# just the divergences which have unique fromids in x but don't in new hope.
div <- get_tocomid(dplyr::select(new_hope_flowline, COMID, FromNode, ToNode),
                   return_dendritic = FALSE,
                   remove_coastal = FALSE)
div <- div[div$tocomid %in%
             new_hope_flowline$COMID[new_hope_flowline$Divergence == 2],]

y <- make_node_topology(x, div)

Description

Cleans up and prepares NHDPlusHR regional data for use as complete NHDPlus data. The primary modification applied is to ensure that any flowpath that exits the domain is labeled as a terminal path and attributes are propagated upstream such that the domain is independently complete.

Usage

make_standalone(flowlines)

Arguments

Value

sf data.frame containing standalone network

Examples

library(dplyr)
library(sf)
source(system.file("extdata/nhdplushr_data.R", package = "nhdplusTools"))

(outlet <- filter(hr_data$NHDFlowline, Hydroseq == min(Hydroseq)))
nrow(filter(hr_data$NHDFlowline, TerminalPa == outlet$Hydroseq))

hr_data$NHDFlowline <- make_standalone(hr_data$NHDFlowline)

(outlet <- filter(hr_data$NHDFlowline, Hydroseq == min(Hydroseq)))
nrow(filter(hr_data$NHDFlowline, TerminalPa == outlet$Hydroseq))

source(system.file("extdata/nhdplushr_data.R", package = "nhdplusTools"))

# Remove mainstem and non-dendritic stuff.
subset <- filter(hr_data$NHDFlowline,
                        StreamLeve > min(hr_data$NHDFlowline$StreamLeve) &
                          StreamOrde == StreamCalc)

subset <- subset_nhdplus(subset$COMID, nhdplus_data = hr_gpkg)$NHDFlowline

plot(sf::st_geometry(hr_data$NHDFlowline))

flowline_mod <- make_standalone(subset)

terminals <- unique(flowline_mod$TerminalPa)

colors <- sample(hcl.colors(length(terminals), palette = "Zissou 1"))

for(i in 1:length(terminals)) {
  fl <- flowline_mod[flowline_mod$TerminalPa == terminals[i], ]
  plot(st_geometry(fl), col = colors[i], lwd = 2, add = TRUE)
}

ol <- filter(flowline_mod, TerminalFl == 1 & TerminalPa %in% terminals)

plot(st_geometry(ol), lwd = 2, add = TRUE)

Description

Given a list of outlets, get their basin boundaries and network and return a leaflet map in EPSG:4326.

Usage

map_nhdplus(
  outlets = NULL,
  bbox = NULL,
  streamorder = NULL,
  nhdplus_data = NULL,
  gpkg = NULL,
  flowline_only = NULL,
  plot_config = NULL,
  overwrite = TRUE,
  cache_data = NULL,
  return_map = FALSE
)

Arguments

Details

map_nhdplus supports several input specifications. An unexported function "as_outlet" is used to convert the outlet formats as described below.

1. if outlets is omitted, the bbox input is required and all nhdplus data in the bounding box is plotted.

2. If outlets is a list of integers, it is assumed to be NHDPlus IDs (comids) and all upstream tributaries are plotted.

3. if outlets is an integer vector, it is assumed to be all NHDPlus IDs (comids) that should be plotted. Allows custom filtering.

4. If outlets is a character vector, it is assumed to be NWIS site ids.

5. if outlets is a list containing only characters, it is assumed to be a list of nldi features and all upstream tributaries are plotted.

6. if outlets is a data.frame with point geometry, a point in polygon match is performed and upstream with tributaries from the identified catchments is plotted.

See plot_nhdplus for details on plot configuration.

Value

data.frame or leaflet map (see return_map)

Examples

map_nhdplus("05428500")

map_nhdplus("05428500", streamorder = 2)

map_nhdplus(list(13293970, 13293750))

source(system.file("extdata/sample_data.R", package = "nhdplusTools"))

map_nhdplus(list(13293970, 13293750), streamorder = 3, nhdplus_data = sample_data)

#return leaflet object
map_nhdplus("05428500", return_map = TRUE)

Description

Provides a full feature network navigation function that will work with local or web service data. Parameter details provide context.

Usage

navigate_network(
  start,
  mode = "UM",
  network = NULL,
  output = "flowlines",
  distance_km = 10,
  trim_start = FALSE,
  trim_stop = FALSE,
  trim_tolerance = 5
)

Arguments

Examples

navigate_network(list(featureSource = "nwissite", featureID = "USGS-06287800"),
                "UM",
                output = "flowlines",
                trim_start = TRUE)


source(system.file("extdata", "walker_data.R", package = "nhdplusTools"))
hydro_location <- list(comid = 5329339,
                      reachcode = "18050005000078",
                       reach_meas = 30)

hydro_location <- sf::st_sf(
  hydro_location,
  geom = nhdplusTools::get_hydro_location(data.frame(hydro_location),
                                          walker_flowline))

net <- navigate_network(hydro_location,
                       mode = "DM", network = walker_flowline,
                       trim_start = TRUE, distance_km = 20)

plot(sf::st_geometry(walker_flowline))
plot(sf::st_geometry(hydro_location), add = TRUE)
plot(sf::st_geometry(net), add = TRUE, col = "blue", lwd = 2)

Description

Navigate the Network Linked Data Index network.

Usage

navigate_nldi(
  nldi_feature,
  mode = "upstreamMain",
  data_source = "flowlines",
  distance_km = 10
)

Arguments

Value

sf data.frame with result

Examples

library(sf)
library(dplyr)

nldi_nwis <- list(featureSource = "nwissite", featureID = "USGS-05428500")

navigate_nldi(nldi_feature = nldi_nwis,
              mode = "upstreamTributaries")$UT %>%
  st_geometry() %>%
  plot()

navigate_nldi(nldi_feature = nldi_nwis,
              mode = "UM")$UM %>%
  st_geometry() %>%
  plot(col = "blue", add = TRUE)



nwissite <- navigate_nldi(nldi_feature = nldi_nwis,
                          mode = "UT",
                          data_source = "nwissite")$UT_nwissite

st_geometry(nwissite) %>%
  plot(col = "green", add = TRUE)

nwissite

Description

Allows specification of a custom path to a source dataset. Typically this will be the national seamless dataset in geodatabase or geopackage format.

Usage

nhdplus_path(path = NULL, warn = FALSE)

Arguments

Value

0 (invisibly) if set successfully, character path if no input.

Examples

nhdplus_path("/data/NHDPlusV21_National_Seamless.gdb")

nhdplus_path("/data/NHDPlusV21_National_Seamless.gdb", warn=FALSE)

nhdplus_path()

Description

Provides an interface to adjust nhdplusTools 'memoise' cache.

Mode and timeout can also be set using environment variables. 'NHDPLUSTOOLS_MEMOISE_CACHE' and 'NHDPLUSTOOLS_MEMOISE_TIMEOUT' are used unless overriden with this function.

Usage

nhdplusTools_cache_settings(mode = NULL, timeout = NULL)

Arguments

Value

list containing settings at time of calling. If inputs are NULL, current settings. If settings are altered, previous setting values.

Description

if left unset, will return the user data dir as returned by 'tools::R_user_dir' for this package.

Usage

nhdplusTools_data_dir(dir = NULL)

Arguments

Value

character path of data directory (silent when setting)

Examples

nhdplusTools_data_dir()

nhdplusTools_data_dir("demo")

nhdplusTools_data_dir(tools::R_user_dir("nhdplusTools"))

Description

Given a list of outlets, get their basin boundaries and network and return a plot in EPSG:3857 Web Mercator Projection.

Usage

plot_nhdplus(
  outlets = NULL,
  bbox = NULL,
  streamorder = NULL,
  nhdplus_data = NULL,
  gpkg = NULL,
  plot_config = NULL,
  basemap = "Esri.NatGeoWorldMap",
  zoom = NULL,
  add = FALSE,
  actually_plot = TRUE,
  overwrite = TRUE,
  flowline_only = NULL,
  cache_data = NULL
)

Arguments

Details

plot_nhdplus supports several input specifications. An unexported function "as_outlet" is used to convert the outlet formats as described below.

1. if outlets is omitted, the bbox input is required and all nhdplus data in the bounding box is plotted.

2. If outlets is a list of integers, it is assumed to be NHDPlus IDs (comids) and all upstream tributaries are plotted.

3. if outlets is an integer vector, it is assumed to be all NHDPlus IDs (comids) that should be plotted. Allows custom filtering.

4. If outlets is a character vector, it is assumed to be NWIS site ids.

5. if outlets is a list containing only characters, it is assumed to be a list of nldi features and all upstream tributaries are plotted.

6. if outlets is a data.frame with point geometry, a point in polygon match is performed and upstream with tributaries from the identified catchments is plotted.

The plot_config parameter is a list with names "basin", "flowline", "outlets", "network_wtbd", and "off_network_wtbd". The following shows the defaults that can be altered.

1. basin

       list(lwd = 1, col = NA, border = "black")

2. flowline

       list(lwd = 1, col = "blue")

3. outlets

      list(default = list(col = "black", border = NA, pch = 19, cex = 1),
           nwissite = list(col = "grey40", border = NA, pch = 17, cex = 1),
           huc12pp = list(col = "white", border = "black", pch = 22, cex = 1),
           wqp = list(col = "red", border = NA, pch = 20, cex = 1))

4. network_wtbd list(lwd = 1, col = "lightblue", border = "black")

5. off_network_wtbd list(lwd = 1, col = "darkblue", border = "black")

If adding additional layers to the plot, data must be projected to EPSG:3857 with 'sf::st_transform(x, 3857)' prior to adding to the plot.

Value

data.frame plot data is returned invisibly in NAD83 Lat/Lon.

Examples

options("rgdal_show_exportToProj4_warnings"="none")
# Beware plot_nhdplus caches data to the default location.
# If you do not want data in "user space" change the default.
old_dir <- nhdplusTools::nhdplusTools_data_dir()
nhdplusTools_data_dir(tempdir())

plot_nhdplus("05428500")

plot_nhdplus("05428500", streamorder = 2)

plot_nhdplus(list(13293970, 13293750))

source(system.file("extdata/sample_data.R", package = "nhdplusTools"))

plot_nhdplus(list(13293970, 13293750), streamorder = 3, nhdplus_data = sample_data)

plot_nhdplus(list(list("comid", "13293970"),
                  list("nwissite", "USGS-05428500"),
                  list("huc12pp", "070900020603"),
                  list("huc12pp", "070900020602")),
             streamorder = 2,
             nhdplus_data = sample_data)


plot_nhdplus(sf::st_as_sf(data.frame(x = -89.36083,
                                     y = 43.08944),
                          coords = c("x", "y"), crs = 4326),
             streamorder = 2,
             nhdplus_data = sample_data)

plot_nhdplus(list(list("comid", "13293970"),
                  list("nwissite", "USGS-05428500"),
                  list("huc12pp", "070900020603"),
                  list("huc12pp", "070900020602")),
             streamorder = 2,
             nhdplus_data = sample_data,
             plot_config = list(basin = list(lwd = 2),
                                outlets = list(huc12pp = list(cex = 1.5),
                                               comid = list(col = "green"))))

bbox <- sf::st_bbox(c(xmin = -89.43, ymin = 43, xmax = -89.28, ymax = 43.1),
                    crs = "+proj=longlat +datum=WGS84 +no_defs")

fline <- sf::read_sf(sample_data, "NHDFlowline_Network")
comids <- nhdplusTools::get_UT(fline, 13293970)

plot_nhdplus(comids)

#' # With Local Data
plot_nhdplus(bbox = bbox, nhdplus_data = sample_data)

# With downloaded data
plot_nhdplus(bbox = bbox, streamorder = 3)

# Can also plot on top of the previous!
plot_nhdplus(bbox = bbox, nhdplus_data = sample_data,
             plot_config = list(flowline = list(lwd = 0.5)))
plot_nhdplus(comids, nhdplus_data = sample_data, streamorder = 3, add = TRUE,
             plot_config = list(flowline = list(col = "darkblue")))

nhdplusTools::nhdplusTools_data_dir(old_dir)

Description

Function to prep NHDPlus data for use by nhdplusTools functions

Usage

prepare_nhdplus(
  flines,
  min_network_size = 0,
  min_path_length = 0,
  min_path_size = 0,
  purge_non_dendritic = TRUE,
  warn = TRUE,
  error = TRUE,
  skip_toCOMID = FALSE,
  align_names = TRUE
)

Arguments

Value

data.frame ready to be used with the refactor_flowlines function.

Examples

source(system.file("extdata", "sample_flines.R", package = "nhdplusTools"))

prepare_nhdplus(sample_flines,
                min_network_size = 10,
                min_path_length = 1,
                warn = FALSE)

Description

Given catchment characteristics to retrieve or process will aggregate and / or split the characteristics according to a lookup table.

Usage

rescale_catchment_characteristics(
  vars,
  lookup_table,
  refactored_areas = NULL,
  catchment_characteristics = NULL,
  catchment_areas = NULL
)

Arguments

Details

NOTE: Since this algorithm works on catchment characteristics that are spatial averages, when splitting, the average condition is apportioned evenly to each split. In some cases, such as with land cover or elevation, this may not be appropriate and source data should be used to derive new characteristics. In addition, this function handles catchment areas for split catchments but makes no adjustments for the length of flowlines in those catchments. Therefore, requests for length-weighted mean values may not be appropriate when working with split catchments.

Examples

vars <- data.frame(characteristic_id = c("CAT_IMPV11","CAT_BASIN_AREA"),
                   summary_statistic = c("area_weighted_mean","sum"))
lookup_table <- data.frame(id = rep(10012268, 2),
                           comid = c(4146596, 4147382),
                           member_comid = c(4146596, 4147382))
rescale_catchment_characteristics(vars, lookup_table)

vars <- data.frame(characteristic_id = c("CAT_ELEV_MIN","CAT_ELEV_MAX"),
                   summary_statistic = c("min","max"))
lookup_table <- data.frame(id = rep(10012268, 2),
                           comid = c(4146596, 4147382),
                           member_comid = c(4146596, 4147382))
rescale_catchment_characteristics(vars, lookup_table)

vars <- data.frame(characteristic_id = c("CAT_EWT","CAT_TWI", "CAT_BASIN_AREA"),
                   summary_statistic = c("area_weighted_mean", "area_weighted_mean","sum"))
lookup_table <- data.frame(id = c(10012268, 10012268, 10024047, 10024048),
                           comid = c(4146596, 4147382, 4147396, 4147396),
                           member_comid = c("4146596", "4147382", "4147396.1", "4147396.2"))
comid_areas <- data.frame(featureid = c("4146596", "4147382", "4147396.1", "4147396.2"),
                               areasqkm = c(0.9558, 11.9790, 6.513294, 1.439999))
rescale_catchment_characteristics(vars, lookup_table, refactored_areas = comid_areas)

Description

RPU Boundaries Raster Processing Unit boundaries

Usage

rpu_boundaries

Format

An object of class "sf"

Description

Saves a subset of the National Seamless database or other nhdplusTools compatible data based on a specified collection of COMIDs. This function uses get_nhdplus for the "download" data source but returns data consistent with local data subsets in a subset file.

Usage

subset_nhdplus(
  comids = NULL,
  output_file = NULL,
  nhdplus_data = NULL,
  bbox = NULL,
  simplified = TRUE,
  overwrite = FALSE,
  return_data = TRUE,
  status = TRUE,
  flowline_only = NULL,
  streamorder = NULL,
  out_prj = 4269
)

Arguments

Details

This function relies on the National Seamless Geodatabase or Geopackage. It can be downloaded here.

The "download" option of this function should be considered preliminary and subject to revision. It does not include as many layers and may not be available permanently.

Value

character path to the saved subset geopackage

Examples

source(system.file("extdata/sample_data.R", package = "nhdplusTools"))

nhdplus_path(sample_data)

sample_flines <- sf::st_zm(sf::read_sf(nhdplus_path(), "NHDFlowline_Network"))

plot(sf::st_geometry(sample_flines),
     lwd = 3)

start_point <- sf::st_sfc(sf::st_point(c(-89.362239, 43.090266)),
                          crs = 4326)

plot(start_point, cex = 1.5, lwd = 2, col = "red", add = TRUE)

start_comid <- discover_nhdplus_id(start_point)

comids <- get_UT(sample_flines, start_comid)

plot(sf::st_geometry(dplyr::filter(sample_flines, COMID %in% comids)),
     add=TRUE, col = "red", lwd = 2)

output_file <- tempfile(fileext = ".gpkg")

subset_nhdplus(comids = comids,
               output_file = output_file,
               nhdplus_data = sample_data,
               overwrite = TRUE,
               status = TRUE)

sf::st_layers(output_file)

catchment <- sf::read_sf(output_file, "CatchmentSP")

plot(sf::st_geometry(catchment), add = TRUE)

waterbody <- sf::read_sf(output_file, "NHDWaterbody")

plot(sf::st_geometry(waterbody),
     col = rgb(0, 0, 1, alpha = 0.5), add = TRUE)

# Cleanup temp
unlink(output_file)

# Download Option:
subset_nhdplus(comids = comids,
               output_file = output_file,
               nhdplus_data = "download",
               overwrite = TRUE,
               status = TRUE, flowline_only = FALSE)

sf::st_layers(output_file)

# NHDPlusHR
source(system.file("extdata/nhdplushr_data.R", package = "nhdplusTools"))

up_ids <- get_UT(hr_data$NHDFlowline, 15000500028335)

sub_gpkg <- file.path(work_dir, "sub.gpkg")
sub_nhdhr <- subset_nhdplus(up_ids, output_file = sub_gpkg,
                            nhdplus_data = hr_gpkg, overwrite = TRUE)

sf::st_layers(sub_gpkg)
names(sub_nhdhr)

plot(sf::st_geometry(hr_data$NHDFlowline), lwd = 0.5)
plot(sf::st_geometry(sub_nhdhr$NHDFlowline), lwd = 0.6, col = "red", add = TRUE)

unlink(output_file)
unlink(sub_gpkg)

Description

Given flowlines and an rpu_code, performs a network-safe subset such that the result can be used in downstream processing. Has been tested to work against the entire NHDPlusV2 domain and satisfies a number of edge cases.

Usage

subset_rpu(fline, rpu, run_make_standalone = TRUE, strict = FALSE)

Arguments

Value

data.frame containing subset network

Examples

source(system.file("extdata/sample_data.R", package = "nhdplusTools"))

sample_flines <- sf::read_sf(sample_data, "NHDFlowline_Network")

subset_rpu(sample_flines, rpu = "07b")

Description

Calls subset_rpu for all raster processing units for the requested vector processing unit.

Usage

subset_vpu(fline, vpu, include_null_rpuid = TRUE, run_make_standalone = TRUE)

Arguments

Value

data.frame containing subset network

Examples

source(system.file("extdata/sample_data.R", package = "nhdplusTools"))

sample_flines <- sf::read_sf(sample_data, "NHDFlowline_Network")

subset_vpu(sample_flines, "07")

Description

VPU Boundaries Vector Processing Unit boundaries

Usage

vpu_boundaries

Format

An object of class "sf"