diff --git a/DESCRIPTION b/DESCRIPTION
index 349f905e59..af8d1371c2 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -79,7 +79,7 @@ Config/testthat/edition: 3
 Encoding: UTF-8
 LazyData: true
 Roxygen: list(markdown = TRUE)
-RoxygenNote: 7.2.3
+RoxygenNote: 7.3.1
 Collate: 
     'ggproto.R'
     'ggplot-global.R'
@@ -280,5 +280,6 @@ Collate:
     'utilities-patterns.R'
     'utilities-resolution.R'
     'utilities-tidy-eval.R'
+    'utilities-unit.R'
     'zxx.R'
     'zzz.R'
diff --git a/NAMESPACE b/NAMESPACE
index 967573b174..340ad4ad10 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -138,6 +138,11 @@ S3method(vec_cast,mapped_discrete.factor)
 S3method(vec_cast,mapped_discrete.integer)
 S3method(vec_cast,mapped_discrete.logical)
 S3method(vec_cast,mapped_discrete.mapped_discrete)
+S3method(vec_cast,simpleUnit.unit)
+S3method(vec_cast,unit.simpleUnit)
+S3method(vec_cast,unit.unit)
+S3method(vec_proxy,simpleUnit)
+S3method(vec_proxy,unit)
 S3method(vec_ptype2,character.mapped_discrete)
 S3method(vec_ptype2,double.mapped_discrete)
 S3method(vec_ptype2,factor.mapped_discrete)
@@ -147,6 +152,10 @@ S3method(vec_ptype2,mapped_discrete.double)
 S3method(vec_ptype2,mapped_discrete.factor)
 S3method(vec_ptype2,mapped_discrete.integer)
 S3method(vec_ptype2,mapped_discrete.mapped_discrete)
+S3method(vec_ptype2,simpleUnit.unit)
+S3method(vec_ptype2,unit.simpleUnit)
+S3method(vec_ptype2,unit.unit)
+S3method(vec_restore,unit)
 S3method(widthDetails,titleGrob)
 S3method(widthDetails,zeroGrob)
 export("%+%")
@@ -286,6 +295,7 @@ export(aes_all)
 export(aes_auto)
 export(aes_q)
 export(aes_string)
+export(after_coord)
 export(after_scale)
 export(after_stat)
 export(alpha)
@@ -719,7 +729,10 @@ import(scales)
 import(vctrs)
 importFrom(glue,glue)
 importFrom(glue,glue_collapse)
+importFrom(grid,arrow)
+importFrom(grid,unit)
 importFrom(lifecycle,deprecated)
+importFrom(scales,alpha)
 importFrom(stats,setNames)
 importFrom(tibble,tibble)
 importFrom(utils,.DollarNames)
diff --git a/R/aes-evaluation.R b/R/aes-evaluation.R
index e128fd2c15..5e6e82a101 100644
--- a/R/aes-evaluation.R
+++ b/R/aes-evaluation.R
@@ -26,6 +26,8 @@
 #'   expression using variables calculated by the stat.
 #' @param after_scale <[`data-masking`][rlang::topic-data-mask]> An aesthetic
 #'   expression using layer aesthetics.
+#' @param after_coord <[`data-masking`][rlang::topic-data-mask]> An aesthetic
+#'   expression using variables calculated by the coord.
 #'
 #' @details
 #' # Staging
@@ -194,15 +196,35 @@ after_scale <- function(x) {
 }
 #' @rdname aes_eval
 #' @export
-stage <- function(start = NULL, after_stat = NULL, after_scale = NULL) {
+after_coord <- function(x) {
+  # Need to put a non-unit() value here as a placeholder until the after_coord
+  # stage of the pipeline, because geoms/scales will not work on unit()s. We
+  # need something that (1) won't affect training of scales (so it can't be an
+  # arbitrary finite number like 0, 1, etc); (2) won't be removed (so it can't
+  # be `NA`); and (3) won't raise errors about missing required aesthetics (so
+  # it can't be `NULL`). The only value satisfying all these properties is `Inf`
+  # (or `-Inf`). If property (3) is relaxed, this could also be `NULL`, though
+  # that would mean users could not use after_coord() on required aesthetics
+  # and would have to do something like stage(Inf, after_coord = ...).
+  rep.int(Inf, length(x))
+}
+after_coord_eval <- function(x) {
+  x
+}
+#' @rdname aes_eval
+#' @export
+stage <- function(start = NULL, after_stat = NULL, after_scale = NULL, after_coord = NULL) {
   start
 }
-stage_calculated <- function(start = NULL, after_stat = NULL, after_scale = NULL) {
+stage_calculated <- function(start = NULL, after_stat = NULL, after_scale = NULL, after_coord = NULL) {
   after_stat
 }
-stage_scaled <- function(start = NULL, after_stat = NULL, after_scale = NULL) {
+stage_scaled <- function(start = NULL, after_stat = NULL, after_scale = NULL, after_coord = NULL) {
   after_scale
 }
+stage_coord <- function(start = NULL, after_stat = NULL, after_scale = NULL, after_coord = NULL) {
+  after_coord
+}
 
 # Regex to determine if an identifier refers to a calculated aesthetic
 match_calculated_aes <- "^\\.\\.([a-zA-Z._]+)\\.\\.$"
@@ -218,6 +240,9 @@ is_calculated_aes <- function(aesthetics, warn = FALSE) {
 is_scaled_aes <- function(aesthetics) {
   vapply(aesthetics, is_scaled, logical(1), USE.NAMES = FALSE)
 }
+is_coord_aes <- function(aesthetics) {
+  vapply(aesthetics, is_coord_stage, logical(1), USE.NAMES = FALSE)
+}
 is_staged_aes <- function(aesthetics) {
   vapply(aesthetics, is_staged, logical(1), USE.NAMES = FALSE)
 }
@@ -260,10 +285,60 @@ is_calculated <- function(x, warn = FALSE) {
 is_scaled <- function(x) {
   is_call(get_expr(x), "after_scale")
 }
+is_coord_stage <- function(x) {
+  is_call(get_expr(x), "after_coord")
+}
 is_staged <- function(x) {
   is_call(get_expr(x), "stage")
 }
 
+#' Compute aesthetic mappings for the after_scale or after_coord stages
+#' @param data data frame of layer data
+#' @param mapping aesthetic mappings containing calls to `stage()`,
+#' `after_scale()`, or `after_coord()`
+#' @param stage one of `"after_scale"` or `"after_coord"`: the stage to apply
+#' @returns modified version of `data` with mappings corresponding to the
+#' given `stage` applied.
+#' @noRd
+compute_staged_aes <- function(data, mapping, stage = "after_scale", call = caller_env()) {
+  if (length(mapping) == 0) return(data)
+
+  # Set up evaluation environment and mask so they return the correct expressions
+  switch(stage,
+    after_scale = {
+      stage_mask <- child_env(emptyenv(), stage = stage_scaled, after_scale = after_scale)
+    },
+    after_coord = {
+      stage_mask <- child_env(emptyenv(), stage = stage_coord, after_coord = after_coord_eval)
+    }
+  )
+  mask <- new_data_mask(as_environment(data, stage_mask), stage_mask)
+  mask$.data <- as_data_pronoun(mask)
+  modified_aes <- lapply(substitute_aes(mapping), eval_tidy, mask, baseenv())
+
+  # Check that all output are valid data
+  nondata_modified <- check_nondata_cols(modified_aes)
+  if (length(nondata_modified) > 0) {
+    issues <- paste0("{.code ", nondata_modified, " = ", as_label(mapping[[nondata_modified]]), "}")
+    names(issues) <- rep("x", length(issues))
+    cli::cli_abort(
+      c(
+        "Aesthetic modifiers returned invalid values",
+        "x" = "The following mappings are invalid",
+        issues,
+        "i" = "Did you map the modifier in the wrong layer?"
+      ),
+      call = call
+    )
+  }
+
+  modified_aes <- vec_recycle_common(!!!modified_aes, .size = nrow(data))
+  names(modified_aes) <- names(rename_aes(mapping))
+  modified_aes <- data_frame0(!!!compact(modified_aes))
+
+  cunion(modified_aes, data)
+}
+
 # Strip dots from expressions
 strip_dots <- function(expr, env, strip_pronoun = FALSE) {
   if (is.null(expr) || is.atomic(expr)) {
diff --git a/R/aes.R b/R/aes.R
index 87870bccb1..924edaaa3a 100644
--- a/R/aes.R
+++ b/R/aes.R
@@ -386,7 +386,11 @@ mapped_aesthetics <- function(x) {
   }
 
   is_null <- vapply(x, is.null, logical(1))
-  names(x)[!is_null]
+  # Ignore mappings using after_coord here, because they may use aesthetics that
+  # are not aesthetics externally supported by the geom, but rather aesthetics
+  # generated internally just prior to the after_coord stage. Thus, an aesthetic
+  # using after_coord() should not generate errors about unknown aesthetics.
+  names(x)[!is_null & !is_coord_aes(x)]
 }
 
 
diff --git a/R/coord-.R b/R/coord-.R
index 8c4313baf7..9f3eccea04 100644
--- a/R/coord-.R
+++ b/R/coord-.R
@@ -170,7 +170,12 @@ Coord <- ggproto("Coord",
     panel_params
   },
 
-  transform = function(data, range) NULL,
+  transform = function(self, data, panel_params) {
+    data <- self$transform_numeric(data, panel_params)
+    compute_staged_aes(data, panel_params$coord_mapping, stage = "after_coord")
+  },
+
+  transform_numeric = function(data, panel_params) NULL,
 
   distance = function(x, y, panel_params) NULL,
 
diff --git a/R/coord-cartesian-.R b/R/coord-cartesian-.R
index 74f46433db..688e9f63b8 100644
--- a/R/coord-cartesian-.R
+++ b/R/coord-cartesian-.R
@@ -93,7 +93,7 @@ CoordCartesian <- ggproto("CoordCartesian", Coord,
     self$range(panel_params)
   },
 
-  transform = function(data, panel_params) {
+  transform_numeric = function(data, panel_params) {
     data <- transform_position(data, panel_params$x$rescale, panel_params$y$rescale)
     transform_position(data, squish_infinite, squish_infinite)
   },
diff --git a/R/coord-flip.R b/R/coord-flip.R
index 1f3848fb8a..46159bbc07 100644
--- a/R/coord-flip.R
+++ b/R/coord-flip.R
@@ -59,7 +59,7 @@ coord_flip <- function(xlim = NULL, ylim = NULL, expand = TRUE, clip = "on") {
 #' @export
 CoordFlip <- ggproto("CoordFlip", CoordCartesian,
 
-  transform = function(data, panel_params) {
+  transform_numeric = function(data, panel_params) {
     data <- flip_axis_labels(data)
     CoordCartesian$transform(data, panel_params)
   },
diff --git a/R/coord-map.R b/R/coord-map.R
index ee0f6ad139..4d1bdd0640 100644
--- a/R/coord-map.R
+++ b/R/coord-map.R
@@ -153,7 +153,7 @@ coord_map <- function(projection="mercator", ..., parameters = NULL, orientation
 #' @export
 CoordMap <- ggproto("CoordMap", Coord,
 
-  transform = function(self, data, panel_params) {
+  transform_numeric = function(self, data, panel_params) {
     trans <- mproject(self, data$x, data$y, panel_params$orientation)
     out <- cunion(trans[c("x", "y")], data)
 
diff --git a/R/coord-polar.R b/R/coord-polar.R
index 1e30adcd2b..49bc6be202 100644
--- a/R/coord-polar.R
+++ b/R/coord-polar.R
@@ -167,7 +167,7 @@ CoordPolar <- ggproto("CoordPolar", Coord,
     panel_params
   },
 
-  transform = function(self, data, panel_params) {
+  transform_numeric = function(self, data, panel_params) {
     arc  <- self$start + c(0, 2 * pi)
     dir  <- self$direction
     data <- rename_data(self, data)
diff --git a/R/coord-radial.R b/R/coord-radial.R
index 70aa211898..0958562a71 100644
--- a/R/coord-radial.R
+++ b/R/coord-radial.R
@@ -219,7 +219,7 @@ CoordRadial <- ggproto("CoordRadial", Coord,
     panel_params
   },
 
-  transform = function(self, data, panel_params) {
+  transform_numeric = function(self, data, panel_params) {
     data <- rename_data(self, data)
     bbox <- panel_params$bbox %||% list(x = c(0, 1), y = c(0, 1))
     arc  <- panel_params$arc  %||% c(0, 2 * pi)
diff --git a/R/coord-sf.R b/R/coord-sf.R
index 331ca4f1f0..dad1793e6c 100644
--- a/R/coord-sf.R
+++ b/R/coord-sf.R
@@ -78,7 +78,7 @@ CoordSf <- ggproto("CoordSf", CoordCartesian,
     self$params$bbox <- bbox
   },
 
-  transform = function(self, data, panel_params) {
+  transform_numeric = function(self, data, panel_params) {
     # we need to transform all non-sf data into the correct coordinate system
     source_crs <- panel_params$default_crs
     target_crs <- panel_params$crs
diff --git a/R/coord-transform.R b/R/coord-transform.R
index 79d651e8af..c87f29c4f7 100644
--- a/R/coord-transform.R
+++ b/R/coord-transform.R
@@ -127,7 +127,7 @@ CoordTrans <- ggproto("CoordTrans", Coord,
     )
   },
 
-  transform = function(self, data, panel_params) {
+  transform_numeric = function(self, data, panel_params) {
     # trans_x() and trans_y() needs to keep Inf values because this can be called
     # in guide_transform.axis()
     trans_x <- function(data) {
diff --git a/R/geom-.R b/R/geom-.R
index 6d4ed6fc55..ae31171b67 100644
--- a/R/geom-.R
+++ b/R/geom-.R
@@ -145,33 +145,7 @@ Geom <- ggproto("Geom",
 
     # If any after_scale mappings are detected they will be resolved here
     # This order means that they will have access to all default aesthetics
-    if (length(modifiers) != 0) {
-      # Set up evaluation environment
-      env <- child_env(baseenv(), after_scale = after_scale)
-      # Mask stage with stage_scaled so it returns the correct expression
-      stage_mask <- child_env(emptyenv(), stage = stage_scaled)
-      mask <- new_data_mask(as_environment(data, stage_mask), stage_mask)
-      mask$.data <- as_data_pronoun(mask)
-      modified_aes <- lapply(substitute_aes(modifiers),  eval_tidy, mask, env)
-
-      # Check that all output are valid data
-      nondata_modified <- check_nondata_cols(modified_aes)
-      if (length(nondata_modified) > 0) {
-        issues <- paste0("{.code ", nondata_modified, " = ", as_label(modifiers[[nondata_modified]]), "}")
-        names(issues) <- rep("x", length(issues))
-        cli::cli_abort(c(
-          "Aesthetic modifiers returned invalid values",
-          "x" = "The following mappings are invalid",
-          issues,
-          "i" = "Did you map the modifier in the wrong layer?"
-        ))
-      }
-
-      names(modified_aes) <- names(rename_aes(modifiers))
-      modified_aes <- data_frame0(!!!compact(modified_aes))
-
-      data <- cunion(modified_aes, data)
-    }
+    data <- compute_staged_aes(data, modifiers, stage = "after_scale")
 
     # Override mappings with params
     aes_params <- intersect(self$aesthetics(), names(params))
diff --git a/R/geom-boxplot.R b/R/geom-boxplot.R
index 289c10cd97..e9b8af6d63 100644
--- a/R/geom-boxplot.R
+++ b/R/geom-boxplot.R
@@ -277,7 +277,6 @@ GeomBoxplot <- ggproto("GeomBoxplot", Geom,
         shape = outlier.shape %||% data$shape[1],
         size = outlier.size %||% data$size[1],
         stroke = outlier.stroke %||% data$stroke[1],
-        fill = NA,
         alpha = outlier.alpha %||% data$alpha[1],
         .size = length(data$outliers[[1]])
       )
diff --git a/R/layer.R b/R/layer.R
index eb590f8dea..b9d8ee347e 100644
--- a/R/layer.R
+++ b/R/layer.R
@@ -427,6 +427,13 @@ Layer <- ggproto("Layer", NULL,
       return(rep(list(zeroGrob()), n))
     }
 
+    # this is probably not the best way to pass down coord mappings
+    aesthetics <- self$computed_mapping
+    coord_mapping <- aesthetics[is_coord_aes(aesthetics) | is_staged_aes(aesthetics)]
+    for (i in seq_along(layout$panel_params)) {
+      layout$panel_params[[i]]$coord_mapping <- coord_mapping
+    }
+
     data <- self$geom$handle_na(data, self$computed_geom_params)
     self$geom$draw_layer(data, self$computed_geom_params, layout, layout$coord)
   }
diff --git a/R/utilities-unit.R b/R/utilities-unit.R
new file mode 100755
index 0000000000..d55730683d
--- /dev/null
+++ b/R/utilities-unit.R
@@ -0,0 +1,53 @@
+# utilities for grid::unit()
+# mostly these are {vctrs} compatibility functions
+# and should probably go in {vctrs}
+
+
+# proxies -----------------------------------------------------------------
+
+#' @export
+vec_proxy.unit <- function(x, ...) {
+  unclass(x)
+}
+
+#' @export
+vec_restore.unit <- function(x, ...) {
+  # replace NAs (NULL entries) with unit's version of NA
+  is_na <- vapply(x, is.null, logical(1))
+  x[is_na] <- vec_proxy(unit(NA_real_, "native"))
+
+  class(x) <- c("unit", "unit_v2")
+  x
+}
+
+#' @export
+vec_proxy.simpleUnit <- function(x, ...) {
+  # turn a simpleUnit into a unit when proxied, because simpleUnit's format
+  # (a numeric vector with an attribute indicating the type of all entries)
+  # does not work properly with many operations, like binding
+  type <- attr(x, "unit")
+  lapply(unclass(x), function(x_i) list(x_i, NULL, type))
+}
+
+
+# casting -----------------------------------------------------------------
+
+null_unit <- function() {
+  # grid::unit() doesn't allow zero-length vectors,
+  # so we have to do this manually
+  structure(list(), class = c("unit", "unit_v2"))
+}
+
+#' @export
+vec_ptype2.unit.unit <- function(x, y, ...) null_unit()
+#' @export
+vec_ptype2.unit.simpleUnit <- function(x, y, ...) null_unit()
+#' @export
+vec_ptype2.simpleUnit.unit <- function(x, y, ...) null_unit()
+
+#' @export
+vec_cast.unit.unit <- function(x, to, ...) x
+#' @export
+vec_cast.unit.simpleUnit <- function(x, to, ...) vec_restore(vec_proxy(x), null_unit())
+#' @export
+vec_cast.simpleUnit.unit <- function(x, to, ...) vec_restore(vec_proxy(x), null_unit())
diff --git a/man/aes_eval.Rd b/man/aes_eval.Rd
index 827bc6a876..c10c658614 100644
--- a/man/aes_eval.Rd
+++ b/man/aes_eval.Rd
@@ -5,6 +5,7 @@
 \alias{after_stat}
 \alias{stat}
 \alias{after_scale}
+\alias{after_coord}
 \alias{stage}
 \title{Control aesthetic evaluation}
 \usage{
@@ -16,7 +17,9 @@ after_stat(x)
 
 after_scale(x)
 
-stage(start = NULL, after_stat = NULL, after_scale = NULL)
+after_coord(x)
+
+stage(start = NULL, after_stat = NULL, after_scale = NULL, after_coord = NULL)
 }
 \arguments{
 \item{x}{<\code{\link[rlang:topic-data-mask]{data-masking}}> An aesthetic expression
@@ -31,6 +34,9 @@ expression using variables calculated by the stat.}
 
 \item{after_scale}{<\code{\link[rlang:topic-data-mask]{data-masking}}> An aesthetic
 expression using layer aesthetics.}
+
+\item{after_coord}{<\code{\link[rlang:topic-data-mask]{data-masking}}> An aesthetic
+expression using variables calculated by the coord.}
 }
 \description{
 Most \link[=aes]{aesthetics} are mapped from variables found in the data.