Skip to contents

Spatially filter an sf layer by any polygon

Source: R/pdok_filter_by.R
pdok_filter_by.Rd

Keeps the features of data that relate to filter_geometry under a spatial predicate. It does what sf::st_filter() does, and additionally reconciles coordinate reference systems for you: filter_geometry is transformed to the CRS of data before filtering.

Usage

pdok_filter_by(data, filter_geometry, predicate = "intersects")

Arguments

data

An sf object to filter (for example a layer loaded with pdok_read()).

filter_geometry

An sf or sfc object whose geometry defines the area of interest.

predicate

The spatial relationship to test, one of "intersects", "within", "contains", "overlaps", "touches", "crosses", "covers", "covered_by", or "disjoint". See the Spatial predicates section below, and sf::geos_binary_pred for the full definitions.

Value

An sf object: the subset of data whose features satisfy predicate with respect to filter_geometry. A zero-row sf is returned (silently) when no feature matches.

Details

filter_geometry can be any polygon: a municipality from pdok_read() on the CBS administrative boundaries, a nature reserve, a water-authority area, a hand-drawn polygon, or another PDOK layer.

The plain-sf equivalent is data[filter_geometry, , op = sf::st_intersects] (after matching CRS); use that if you prefer to drop down to sf.

Spatial predicates

The predicate selects how a feature must relate to filter_geometry to be kept. Each value calls the matching sf predicate function sf::st_<predicate>() (so "within" uses sf::st_within(), "intersects" uses sf::st_intersects(), and so on). The common choices:

  • "intersects" (default) — the feature touches filter_geometry in any way (overlaps, is inside, or shares a boundary). The most permissive; the usual choice for "everything in this area".

  • "within" — the feature lies entirely inside filter_geometry. Use this to exclude features that only stick partly into the area.

  • "contains" — the feature entirely encloses filter_geometry (the inverse of "within").

  • "disjoint" — the feature does not touch filter_geometry at all (everything outside the area).

"overlaps", "touches", "crosses", "covers" and "covered_by" are the remaining, more specialized DE-9IM relationships. For the precise definition of each, see the sf help page sf::geos_binary_pred (the topic that documents st_intersects(), st_within(), and the rest).

See also

sf::geos_binary_pred for the definition of each spatial predicate; sf::st_filter(), the sf filter this wraps; and pdok_read(), whose filter_by argument applies this filter while loading.

Examples

# \donttest{
# All national parks that intersect the province of Utrecht
utrecht <- pdok_read(
  "cbs/gebiedsindelingen", "provincie_gegeneraliseerd",
  datetime = 2024
)
utrecht <- utrecht[utrecht$statnaam == "Utrecht", ]
parks <- pdok_read("rvo/nationale-parken-geharmoniseerd", "protectedsite")
#> ⠙ Downloading PDOK features: 22 fetched
pdok_filter_by(parks, utrecht)
#> Simple feature collection with 1 feature and 14 fields
#> Geometry type: MULTIPOLYGON
#> Dimension:     XY
#> Bounding box:  xmin: 5.147185 ymin: 51.95757 xmax: 5.566058 ymax: 52.22923
#> Geodetic CRS:  WGS 84
#> # A tibble: 1 × 15
#>   id          gml_id language legalfoundationdate legalfoundationdocum…¹ localid
#>   <chr>       <chr>  <chr>    <chr>               <chr>                  <chr>  
#> 1 4d5b664b-7… G_4a2… nld      20131001            prov. Utrecht          L_4a24…
#> # ℹ abbreviated name: ¹​legalfoundationdocument
#> # ℹ 9 more variables: namespace <chr>, namestatus <chr>, nativeness <chr>,
#> #   pronunciation <chr>, script <chr>, siteprotectionclassification <chr>,
#> #   sourceofname <chr>, text <chr>, geometry <MULTIPOLYGON [°]>
# }