geom_raster and ggmap

[ted]

Hi,
I am having a problem combining geom_raster and ggmap. They each work independently, but when I try to combine the two I get the error message "Error: geom_raster only works with Cartesian coordinates"

My code follows:
mp <- get_map() 
ggmap(mp)+geom_raster(aes(x=x,y=y,fill=category),data=myRasterData)

Thanks for any help!
Best,
Ted

[David Kahle]

Try inset_raster. (It's analogous to annotation_raster.)  I don't remember if there is an analogue to geom_raster (they take a slightly different argument shape).

Cheers

dk.

--
You received this message because you are subscribed to the ggplot2 mailing list.
Please provide a reproducible example: https://github.com/hadley/devtools/wiki/Reproducibility
 
To post: email ggp...@googlegroups.com
To unsubscribe: email ggplot2+u...@googlegroups.com
More options: http://groups.google.com/group/ggplot2

[ted]

Hi David--
First -- thanks for the great package and your insights.

Playing around more I found out the issue.  ggmap projects the data using the coord_map('mercator') projection.  This is not compatible with geom_raster (hence the error message), and therefore to make  geom_raster work I simply include coord_cartesian() in the ggplot line of code, which forces the resulting image to use the Cartesian coordinate system. 

Best,
ted

Cheers

dk.

To unsubscribe: email ggplot2+unsubscribe@googlegroups.com
More options: http://groups.google.com/group/ggplot2

[David Kahle]

> Hi David--
> First -- thanks for the great package and your insights.

:)

> Playing around more I found out the issue. ggmap projects the data using the coord_map('mercator') projection. This is not compatible with geom_raster (hence the error message), and therefore to make geom_raster work I simply include coord_cartesian() in the ggplot line of code, which forces the resulting image to use the Cartesian coordinate system.

You're right; that's why it doesn't work, and putting that in does indeed make it "work". Unfortunately though, it turns out that's actually not so good of an idea. The map image itself is in a mercator projection already. If you change the coordinates to Cartesian you'll trick ggplot2 into thinking that the map is in the correct projection, but in fact you haven't changed the map, so you still pay the piper in the end. The effect will be that your data overlays will be off relative to the map.

That may be a bit confusing. If more of an explanation would help just let me know and I can explain off-list.

Cheers
david.

[David Kahle]

[I've added this back to the list because I think it will be helpful for others.]

It's not quite that, actually.  In fact, ggplot knows about 1 and 2.

What ggplot2 doesn't know about is the map.  You can think of the map itself as a rectangle (just set each pixel in the image to gray or something); it could have been created by geom_rect.  As a rectangle, it lives in the plot and has to abide by the full conventions of the grammar in the sense that it's dimensions live in the x-y coordinate system of the graph.  However, the map, unlike the rectangle, has a meaningful fill value which cannot be interpreted with a kind of scale.  

The map is not really part of the plot.  It is more like an annotation, like the title of the plot or the legend.  But in reality it is an amphibian – both plot element and annotation.  Its plot component is where it is situated in the coordinate system.  Its annotation component is the map part – the fill colors on the rectangle which have meaning but not according to a scale other than identity and even then not because of that scale.  However, the x-y aesthetics of the fill colors are still meaningful and have a scale which is set by the source of the image (Google, OSM, etc.) and not the ggplot2 user, even when we do use coord_cartesian().  Rather, you will have a picture whose image is in a projected coordinate system forced into a coordinate system which is not projected.  There will appear to be no mistake, because the map image has no visible coordinate system.  But it's there.

Perhaps the best way to see it is through example :

library(ggmap)

gc <- geocode('the white house')

# the bad way  

qmap('washington, dc', zoom = 10) + # ok

  geom_point(colour = 'red', data = gc) +

  coord_cartesian()

  

qmap('washington, dc', zoom = 7) + # slightly north?

  geom_point(colour = 'red', data = gc) +

  coord_cartesian()

  

qmap('washington, dc', zoom = 4) + # in northern pennsylvania

  geom_point(colour = 'red', data = gc) +

  coord_cartesian()  

  

qmap('washington, dc', zoom = 1) + # potus = santa claus

  geom_point(colour = 'red', data = gc) +

  coord_cartesian() 

  

  

  

  

# the good way

qmap('washington, dc', zoom = 10) + # ok

  geom_point(colour = 'red', data = gc)

  

qmap('washington, dc', zoom = 7) + # still ok

  geom_point(colour = 'red', data = gc)  

  

qmap('washington, dc', zoom = 4) + # still ok

  geom_point(colour = 'red', data = gc)    

  

qmap('washington, dc', zoom = 3) + # still ok

  geom_point(colour = 'red', data = gc)    

  

qmap('washington, dc', zoom = 2) + # there are latitude restrictions on projection

  geom_point(colour = 'red', data = gc)      

Does that help?

david.

On Jun 29, 2012, at 1:48 PM, Ted Rosenbaum wrote:

HI David--
Thanks again for your reply.

I guess a bit of clarification on your point about the projections would be helpful.  As far as I understand it, there are two issues here:
1. That coordinate system that the data is stored with internally.
2. How the data is displayed on the screen -- (mapping the (x,y) coordinates from the data onto (x,y) coordinates on the screen.)

GGplot presumably knows nothing about 1 (the type of coordinate system) and just stores (x,y) data. 
What (I think) GGmap does is takes the latitude and longitude coordinates of the google map and maps it using a mercator projection (but just stores the latitude and longitude).  However, I could also map those (lat,long) coordinates using any other projection -- including cartesian.  This would impact the way it looks, but not much else. 

Therefore, if I overlay data using geom_raster, where the underlying coordinates are in latitudes and longitudes the pictures should overlay correctly, and the only issue the the way in which the (correctly overlayed map) appears on the screen. 

So therefore, based on this logic I think what I did is fine.  But -- I'm probably missing something somewhere, and I look forward to hearing what it is :-).

Thanks!
Best,
Ted

[Ted Rosenbaum]

Yes that makes sense. 
It also makes it seem like a workaround is hard. I guess the only "workaround" is for geom_raster to support projected coordinate systems (which I realize is somewhat complicated). 
Thanks for the help.

[David Kahle]

Well, if you reshaped your data you could use inset_raster no problem.  Otherwise, if you really must use geom_raster, run the following and then it should work. (With any luck.)

library(proto)

library(grid)

geom_raster <- function (mapping = NULL, data = NULL, stat = "identity", position = "identity", hjust = 0.5, vjust = 0.5, interpolate = FALSE, ...) { 

  stopifnot(is.numeric(hjust), length(hjust) == 1)

  stopifnot(is.numeric(vjust), length(vjust) == 1)

  

  ggplot2:::GeomRaster$new(mapping = mapping, data = data, stat = stat, position = position, hjust = hjust, vjust = vjust, interpolate = interpolate, ...)

}

GeomRaster <- proto(ggplot2:::Geom, {

  objname <- "raster"

  

  reparameterise <- function(., df, params) {

    hjust <- params$hjust %||% 0.5

    vjust <- params$vjust %||% 0.5

    

    w <- resolution(df$x, FALSE)

    h <- resolution(df$y, FALSE)

    

    df$xmin <- df$x - w * (1 - hjust)

    df$xmax <- df$x + w * hjust

    df$ymin <- df$y - h * (1 - vjust)

    df$ymax <- df$y + h * vjust

    df

  }

  

  # This is a dummy function to make sure that vjust and hjust are recongised

  # as parameters and are accessible to reparameterise.

  draw <- function(vjust = 0.5, hjust = 0.5) {}

  

  draw_groups <- function(., data, scales, coordinates, interpolate = FALSE, ...) {

    data <- ggplot2:::remove_missing(data, TRUE, c("x", "y", "fill"), 

      name = "geom_raster")

    data <- ggplot2:::coord_transform(coordinates, data, scales)

    # Convert vector of data to raster

    x_pos <- as.integer((data$x - min(data$x)) / resolution(data$x, FALSE))

    y_pos <- as.integer((data$y - min(data$y)) / resolution(data$y, FALSE))

    

    nrow <- max(y_pos) + 1

    ncol <- max(x_pos) + 1

    raster <- matrix(NA_character_, nrow = nrow, ncol = ncol)

    raster[cbind(nrow - y_pos, x_pos + 1)] <- alpha(data$fill, data$alpha)

    

    # Figure out dimensions of raster on plot

    x_rng <- c(min(data$xmin, na.rm = TRUE), max(data$xmax, na.rm = TRUE))

    y_rng <- c(min(data$ymin, na.rm = TRUE), max(data$ymax, na.rm = TRUE))

    

    rasterGrob(raster, x = mean(x_rng), y = mean(y_rng), 

      width = diff(x_rng), height = diff(y_rng), 

      default.units = "native", interpolate = interpolate)

  }

  default_stat <- function(.) ggplot2:::StatIdentity

  default_aes <- function(.) aes(fill = "grey20", alpha = NA)

  required_aes <- c("x", "y")

  guide_geom <- function(.) "polygon"

})

I might have missed some packages you need, so let me know how it goes.

david.