plot raster over google map
On Thu, Feb 28, 2013 at 5:03 PM, Ross Ahmed <rossahmed at googlemail.com> wrote:
I've downloaded the following google map:
library(raster)
myExtent <- matrix(c(-1.87710, -1.772096, 55.61399, 55.682171), ncol=2,
byrow=T)
g <- gmap(myExtent, type='hybrid')
plot(g, interpolate=T)
I'm trying to plot the raster over the google map:
rownames(myExtent) <- c('long', 'lat')
colnames(myExtent) <- c('min', 'max')
rasterExtent <- extent(myExtent)
myRaster <- raster(ext=rasterExtent,crs="+proj=merc +a=6378137 +b=6378137
+lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null
+no_defs")
values(myRaster) <- rnorm(ncell(myRaster))
You've created a raster with an extent in lat-long coordinates, but assigned it a very non lat-long projection CRS. gmap has done a couple of things, it has got the google map for the lat-long extent you requested, but returned a raster in google's CRS. It's not a lat-long grid. > extent(g) class : Extent xmin : -213661.6 xmax : -192526.8 ymin : 7476426 ymax : 7500886 > extent(myExtent) class : Extent xmin : -1.8771 xmax : -1.772096 ymin : 55.61399 ymax : 55.68217
plot(myRaster, add=T, alpha=0.5) However, there's no sign of the raster in the plot. I presume theres a projection issue here. How can I get the raster to plot over the top of the google map?
You need to stop lying to your raster, and tell it it really is a raster in lat-long, then project it to the google map CRS: > projection(myRaster)="+init=epsg:4326" # lat-long WGS84 > pRaster =projectRaster(myRaster, crs=projection(g)) # reproject to the google CRS (no need to type it all in again) > plot(pRaster,add=TRUE) # shazam Barry