Get NHDPlusHR Data

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

get_nhdphr(
  AOI = NULL,
  ids = NULL,
  type = NULL,
  reachcode = NULL,
  t_srs = NULL,
  buffer = 0.5,
  page_size = 2000
)

Arguments

AOI

sf (MULTI)POINT or (MULTI)POLYGON. An 'area of interest' can be provided as either a location (sf POINT) or area (sf POLYGON) in any Spatial Reference System.

ids

character vector of nhdplusid ids

type

character. Type of feature to return e.g. c("networknhdflowline", nonnetworknhdflowline", nhdwaterbody", "nhdpluscatchment"). If NULL (default) a data.frame of available types is returned

reachcode

character vector of reachcodes NOTE: performance of this query is currently very poor, spatial queries are the primary use of this function.

t_srs

character (PROJ string or EPSG code) or numeric (EPSG code). A user specified - target -Spatial Reference System (SRS/CRS) for returned objects. Will default to the CRS of the input AOI if provided, and to 4326 for ID requests.

buffer

numeric. The amount (in meters) to buffer a POINT AOI by for an extended search. Default = 0.5

page_size

numeric default number of features to request at a time. Reducing may help if 500 errors are experienced.

Value

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

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

Examples

# \donttest{
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_nhdphr(AOI = AOI, type = "networknhdflowline")
point <- get_nhdphr(AOI = AOI, type = "nhdpoint")
waterbody <- get_nhdphr(AOI = AOI, type = "nhdwaterbody")

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


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

get_nhdphr(reachcode = "13020101021927", type = "networknhdflowline")
#> Warning: Failed to get JSON from https://hydro.nationalmap.gov/arcgis/rest/services/NHDPlus_HR/MapServer/3/query: HTTP 504 Gateway Timeout.
#> Warning: No networknhdflowline features found in area of interest.
#> NULL
# }