Nothing
#' Latitude circle
#'
#' Create a circular polygon using a latitude value in a map projection. The
#' longitude range can be modified from global to give a portion of a circle.
#'
#' The argument `nverts` controls the total number of vertices of the circle
#' within a linearly within the `lonlim` range of longitudes at `latitude`
#'
#' This is for use on classic polar projections centred on the north or the
#' south pole, particularly Polar Stereographic and Lambert Azimuthal Equal Area
#' but will also work with some caveats on other families and situations. We
#' have not explored this more general use. Feel free to contact the maintainer
#' if you have interest in less typical usage or find problems.
#'
#' @param latitude latitude value for the boundary, defaults to 0
#' @param crs map projection to use, defaults to southern Polar Stereographic true scale at -71S
#' @param lonlim the range of longitude to use, defaults to entire globe
#' @param nverts total number of vertices to use, see Details
#' @return SpatialPolygonsDataFrame
#' @export
#' @importFrom reproj reproj
#' @examples
#' latitudecircle(seq(0, -65, by = -5))
#' library(raster)
#' plot(ice)
#' circ <- latitudecircle(-71, crs = projection(ice))
#' plot(circ, add = TRUE)
latitudecircle <- function(latitude = 0, crs = "+proj=stere +lon_0=0 +lat_0=-90 +lat_ts=-71 +ellps=WGS84",
lonlim = c(-180, 180),
nverts = 1800) {
stopifnot(is.numeric(latitude))
if (length(latitude) > 1) warning("ignoring multiple latitude values, using first")
latitude <- latitude[1L]
raster::spPolygons(
reproj::reproj(cbind(seq(lonlim[1], lonlim[2], length = nverts), latitude), crs,
source = "+proj=longlat +datum=WGS84 +no_defs")[,1:2],
crs = sp::CRS(crs, doCheckCRSArgs = FALSE), attr = data.frame(latitude = latitude))
}
#' Latitude mask for polar raster
#'
#' Mask out values based on latitude for a raster. This works by finding all cells at
#' latitudes less than `latitude` and setting them to missing. If `southern = FALSE`
#' the inequality is reversed, and all cells at latitudes greater than `latitude` are
#' masked out.
#'
#' The `trim` option allows for the result to be reduced to the common bounding box
#' within which any row or column has a non-missing value.
#' @param x a raster layer
#' @param latitude maximum latitude (effectively a minimum latitude if `southern = FALSe`)
#' @param southern flag for whether south-polar context is used, default is `TRUE`
#' @param trim if `TRUE` runs `raster::trim` on the result, to remove `NA` margin
#' @param ... ignored currently
#'
#' @seealso [raster::trim], [latitudecircle]
#' @return RasterLayer
#' @export
#'
#' @examples
#' library(raster)
#' plot(latmask(ice, -60))
#' plot(latmask(ice, -60, trim = TRUE))
#' ice[!ice > 0] <- NA
#' plot(ice)
#' plot(latmask(ice, -55, trim = TRUE))
latmask <- function(x, latitude = 0, southern = TRUE, trim = FALSE, ...) {
if (length(latitude) > 1) warning("ignoring multiple latitude values, using first")
latitude <- latitude[1L]
if (raster::isLonLat(x)) {
xy <- sp::coordinates(x)
} else {
xy <- reproj::reproj(sp::coordinates(x), "+proj=longlat +datum=WGS84", source = raster::projection(x))
}
if (southern) x[xy[,2] > (latitude )] <- NA else x[xy[,2] < (latitude)] <- NA
if (trim) x <- raster::trim(x)
x
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.