# R/st_ellipse.R In michaeldorman/nngeo: k-Nearest Neighbor Join for Spatial Data

#### Documented in st_ellipse

```#' Calculate ellipse polygon
#'
#' The function calculates ellipse polygons, given centroid locations and sizing on the x and y axes.
#' @param pnt Object of class \code{sf} or \code{sfc} (type \code{"POINT"}) representing centroid locations
#' @param ex Size along x-axis, in CRS units
#' @param ey Size along y-axis, in CRS units
#' @param res Number of points the ellipse polygon consists of (default \code{30})
#' @return Object of class \code{sfc} (type \code{"POLYGON"}) containing ellipse polygons
#' @references Based on StackOverflow answer by user \code{fdetsch}:
#'
#' @export
#'
#' @examples
#' # Sample data
#' dat = data.frame(
#'   x = c(1, 1, -1, 3, 3),
#'   y = c(0, -3, 2, -2, 0),
#'   ex = c(0.5, 2, 2, 0.3, 0.6),
#'   ey = c(0.5, 0.2, 1, 1, 0.3),
#'   stringsAsFactors = FALSE
#' )
#' dat = st_as_sf(dat, coords = c("x", "y"))
#' dat
#'
#' # Plot 1
#' plot(st_geometry(dat), graticule = TRUE, axes = TRUE, main = "Input")
#' text(st_coordinates(dat), as.character(1:nrow(dat)), pos = 2)
#'
#' # Calculate ellipses
#' el = st_ellipse(pnt = dat, ex = dat\$ex, ey = dat\$ey)
#'
#' # Plot 2
#' plot(el, graticule = TRUE, axes = TRUE, main = "Output")
#' plot(st_geometry(dat), pch = 3, add = TRUE)
#' text(st_coordinates(dat), as.character(1:nrow(dat)), pos = 2)

st_ellipse = function(pnt, ex, ey, res = 30) {

# Checks
stopifnot(length(st_geometry(pnt)) == length(ex))
stopifnot(length(ex) == length(ey))
stopifnot(class(ex) == "numeric")
stopifnot(class(ey) == "numeric")
stopifnot(all(!st_is_empty(pnt)))
stopifnot(all(!is.na(ex)))
stopifnot(all(!is.na(ey)))
stopifnot(class(st_geometry(pnt))[1] == "sfc_POINT")

# Point coordinates
coords = st_coordinates(pnt)

# Create ellipses
result = list()
for(i in 1:nrow(coords)) {
theta = seq(0, 2 * pi, length = res)
x = coords[i, 1] - ex[i] * cos(theta)
y = coords[i, 2] - ey[i] * sin(theta)
elp = cbind(x[!is.na(x)], y[!is.na(y)])
pnt = st_multipoint(elp)
pol = st_cast(pnt, "POLYGON")
result[[i]] = pol
}

# Combine
result = st_sfc(result)

# Return result
return(result)

}
```
michaeldorman/nngeo documentation built on June 4, 2022, 12:49 a.m.