1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233
|
#' Map projections.
#'
#' This coordinate system provides the full range of map projections available
#' in the mapproj package.
#'
#' This is still experimental, and if you have any advice to offer regarding
#' a better (or more correct) way to do this, please let me know
#'
#' @export
#' @param projection projection to use, see
#' \code{\link[mapproj]{mapproject}} for list
#' @param ... other arguments passed on to
#' \code{\link[mapproj]{mapproject}}
#' @param orientation projection orientation, which defaults to
#' \code{c(90, 0, mean(range(x)))}. This is not optimal for many
#' projections, so you will have to supply your own. See
#' \code{\link[mapproj]{mapproject}} for more information.
#' @param xlim manually specific x limits (in degrees of lontitude)
#' @param ylim manually specific y limits (in degrees of latitude)
#' @export
#' @examples
#' if (require("maps")) {
#' # Create a lat-long dataframe from the maps package
#' nz <- map_data("nz")
#' nzmap <- ggplot(nz, aes(x=long, y=lat, group=group)) +
#' geom_polygon(fill="white", colour="black")
#'
#' # Use cartesian coordinates
#' nzmap
#' # With default mercator projection
#' nzmap + coord_map()
#' # Other projections
#' nzmap + coord_map("cylindrical")
#' nzmap + coord_map("azequalarea",orientation=c(-36.92,174.6,0))
#'
#' states <- map_data("state")
#' usamap <- ggplot(states, aes(x=long, y=lat, group=group)) +
#' geom_polygon(fill="white", colour="black")
#'
#' # Use cartesian coordinates
#' usamap
#' # With mercator projection
#' usamap + coord_map()
#' # See ?mapproject for coordinate systems and their parameters
#' usamap + coord_map("gilbert")
#' usamap + coord_map("lagrange")
#'
#' # For most projections, you'll need to set the orientation yourself
#' # as the automatic selection done by mapproject is not available to
#' # ggplot
#' usamap + coord_map("orthographic")
#' usamap + coord_map("stereographic")
#' usamap + coord_map("conic", lat0 = 30)
#' usamap + coord_map("bonne", lat0 = 50)
#'
#' # World map, using geom_path instead of geom_polygon
#' world <- map_data("world")
#' worldmap <- ggplot(world, aes(x=long, y=lat, group=group)) +
#' geom_path() +
#' scale_y_continuous(breaks=(-2:2) * 30) +
#' scale_x_continuous(breaks=(-4:4) * 45)
#'
#' # Orthographic projection with default orientation (looking down at North pole)
#' worldmap + coord_map("ortho")
#' # Looking up up at South Pole
#' worldmap + coord_map("ortho", orientation=c(-90, 0, 0))
#' # Centered on New York (currently has issues with closing polygons)
#' worldmap + coord_map("ortho", orientation=c(41, -74, 0))
#' }
coord_map <- function(projection="mercator", ..., orientation = NULL, xlim = NULL, ylim = NULL) {
try_require("mapproj")
coord(
projection = projection,
orientation = orientation,
limits = list(x = xlim, y = ylim),
params = list(...),
subclass = "map"
)
}
#' @export
coord_transform.map <- function(coord, data, details) {
trans <- mproject(coord, data$x, data$y, details$orientation)
out <- cunion(trans[c("x", "y")], data)
out$x <- rescale(out$x, 0:1, details$x.proj)
out$y <- rescale(out$y, 0:1, details$y.proj)
out
}
mproject <- function(coord, x, y, orientation) {
suppressWarnings(mapproject(x, y,
projection = coord$projection,
parameters = coord$params,
orientation = orientation
))
}
#' @export
coord_distance.map <- function(coord, x, y, details) {
max_dist <- dist_central_angle(details$x.range, details$y.range)
dist_central_angle(x, y) / max_dist
}
#' @export
coord_aspect.map <- function(coord, ranges) {
diff(ranges$y.proj) / diff(ranges$x.proj)
}
#' @export
coord_train.map <- function(coord, scales) {
# range in scale
ranges <- list()
for (n in c("x", "y")) {
scale <- scales[[n]]
limits <- coord$limits[[n]]
if (is.null(limits)) {
expand <- coord_expand_defaults(coord, scale, n)
range <- scale_dimension(scale, expand)
} else {
range <- range(scale_transform(scale, limits))
}
ranges[[n]] <- range
}
orientation <- coord$orientation %||% c(90, 0, mean(ranges$x))
# Increase chances of creating valid boundary region
grid <- expand.grid(
x = seq(ranges$x[1], ranges$x[2], length = 50),
y = seq(ranges$y[1], ranges$y[2], length = 50)
)
ret <- list(x = list(), y = list())
# range in map
proj <- mproject(coord, grid$x, grid$y, orientation)$range
ret$x$proj <- proj[1:2]
ret$y$proj <- proj[3:4]
for (n in c("x", "y")) {
out <- scale_break_info(scales[[n]], ranges[[n]])
ret[[n]]$range <- out$range
ret[[n]]$major <- out$major_source
ret[[n]]$minor <- out$minor_source
ret[[n]]$labels <- out$labels
}
details <- list(
orientation = orientation,
x.range = ret$x$range, y.range = ret$y$range,
x.proj = ret$x$proj, y.proj = ret$y$proj,
x.major = ret$x$major, x.minor = ret$x$minor, x.labels = ret$x$labels,
y.major = ret$y$major, y.minor = ret$y$minor, y.labels = ret$y$labels
)
details
}
#' @export
coord_render_bg.map <- function(coord, details, theme) {
xrange <- expand_range(details$x.range, 0.2)
yrange <- expand_range(details$y.range, 0.2)
# Limit ranges so that lines don't wrap around globe
xmid <- mean(xrange)
ymid <- mean(yrange)
xrange[xrange < xmid - 180] <- xmid - 180
xrange[xrange > xmid + 180] <- xmid + 180
yrange[yrange < ymid - 90] <- ymid - 90
yrange[yrange > ymid + 90] <- ymid + 90
xgrid <- with(details, expand.grid(
y = c(seq(yrange[1], yrange[2], len = 50), NA),
x = x.major
))
ygrid <- with(details, expand.grid(
x = c(seq(xrange[1], xrange[2], len = 50), NA),
y = y.major
))
xlines <- coord_transform(coord, xgrid, details)
ylines <- coord_transform(coord, ygrid, details)
if (nrow(xlines) > 0) {
grob.xlines <- element_render(
theme, "panel.grid.major.x",
xlines$x, xlines$y, default.units = "native"
)
} else {
grob.xlines <- zeroGrob()
}
if (nrow(ylines) > 0) {
grob.ylines <- element_render(
theme, "panel.grid.major.y",
ylines$x, ylines$y, default.units = "native"
)
} else {
grob.ylines <- zeroGrob()
}
ggname("grill", grobTree(
element_render(theme, "panel.background"),
grob.xlines, grob.ylines
))
}
#' @export
coord_render_axis_h.map <- function(coord, details, theme) {
if (is.null(details$x.major)) return(zeroGrob())
x_intercept <- with(details, data.frame(
x = x.major,
y = y.range[1]
))
pos <- coord_transform(coord, x_intercept, details)
guide_axis(pos$x, details$x.labels, "bottom", theme)
}
#' @export
coord_render_axis_v.map <- function(coord, details, theme) {
if (is.null(details$y.major)) return(zeroGrob())
x_intercept <- with(details, data.frame(
x = x.range[1],
y = y.major
))
pos <- coord_transform(coord, x_intercept, details)
guide_axis(pos$y, details$y.labels, "left", theme)
}
|