[R] How to overlay a vector map (polygon) on raster maps?

Ben Tupper btupper @end|ng |rom b|ge|ow@org
Thu Apr 11 18:15:44 CEST 2019


Hi,

It's best to keep all of the replies on the list - you will get better answers and leave a trail for others with similar questions to follow.  If you need more help, I strongly suggest that you start a fresh question on r-sig-geo.

I suppose you could try the panel argument to levelplot().  Using the panel argument will modify each raster - that is each rendering of elements of your raster stack - by performing what ever task you put in the panel function.  It's all a bit mysterious to me how it really works which is why I often gravitate toward the more obvious-to-me layering that the latticeExtra package provides with `levelplot(something) + layer(more stuff)`. 

library(rasterVis)
library(sp)

set.seed(10)
x = runif(2000000, -0.0005, .9875)
y = runif(2000000, -0.0008, .99)
xmat = matrix(x, nrow = 500)
ymat = matrix(x, nrow = 500)
xras = raster(xmat)
yras = raster(ymat)
min_ = min(minValue(xras), minValue(yras))
max_ = max(maxValue(xras), maxValue(yras))
r.range = c(min_, max_)
Poly <- as(raster::extent(xras) + c(.15, -.32, .25, -.10), "SpatialPolygons")
levelplot(stack(xras, yras), 
	col.regions = rev(rainbow(99, start=0, end=1)), 
	colorkey = list(space = "bottom"),
	panel = function(...){
		panel.levelplot(...)
		sp.polygons(Poly, col = 'black', lwd = 3)
		}
	)

Cheers,
Ben



> On Apr 11, 2019, at 11:50 AM, Kristi Glover <kristi.glover using hotmail.com> wrote:
> 
> Thank you Ben for the link. It has lots of the information. One of the example bellow, here how can we overlay a polygon (let's say river map) on the two raster images? The two raster images are of Maximum and Minimum temperature of a specific area. 
> library(rasterVis)
> 
> set.seed
> (10)
> 
> x 
> = runif(2000000, -0.0005, .9875)
> 
> y 
> = runif(2000000, -0.0008, .99)
> 
> xmat 
> = matrix(x, nrow = 500)
> 
> ymat 
> = matrix(x, nrow = 500)
> 
> xras 
> = raster(xmat)
> 
> yras 
> = raster(ymat)
> 
> min_ 
> = min(minValue(xras), minValue(yras))
> 
> max_ 
> = max(maxValue(xras), maxValue(yras))
> 
> r.range 
> = c(min_, max_)
> 
> 
> levelplot
> (stack(xras, yras), col.regions = rev(rainbow(99, start=0, end=1)), colorkey = list(space = "bottom"))
> thank you so much for your help.
> 
> 
> 
> From: Ben Tupper <btupper using bigelow.org>
> Sent: April 11, 2019 9:43 AM
> To: Kristi Glover
> Cc: r-help using r-project.org
> Subject: Re: [R] How to overlay a vector map (polygon) on raster maps?
>  
> Hi,
> 
> That's great topic to search on RSeek.org
> 
> https://rseek.org/?q=plot+multiple+rasters+with+one+legend
> 
> or to pose a question about on r-sig-geo
> 
> Cheers,
> Ben
> 
>> On Apr 11, 2019, at 11:32 AM, Kristi Glover <kristi.glover using hotmail.com> wrote:
>> 
>> Dear Ben, 
>> Thank you very much for the message. I run it and it produced three separate images with X and y Axis  and a legend for each image. I was thinking to plot all of these three images with a single legend and only X axis value at the bottom's image and y values for each image.  
>> 
>> I added the following code on your code
>> grid.arrange(PP[[1]],PP[[2]],PP[[3]])
>>  
>> But, as I mentioned above, I can get three separate images with its own legend and X and Y axix
>> Thanks,
>> 
>> From: Ben Tupper <btupper using bigelow.org>
>> Sent: April 11, 2019 7:13 AM
>> To: Kristi Glover
>> Cc: r-help using r-project.org
>> Subject: Re: [R] How to overlay a vector map (polygon) on raster maps?
>>  
>> Hi,
>> 
>> I think you want to build a levelplot object with polygon overlaid for each layer.  Like the link below shows but with the added layer per your example.  
>> 
>> https://oscarperpinan.github.io/rastervis/FAQ.html#several_rasters
>> 
>> Also, you will get bucket loads of spatial-centric help using the r-sig-geo mailing list; check it out herehttps://stat.ethz.ch/mailman/listinfo/r-sig-geo.
>> 
>> Cheers,
>> Ben
>> 
>> ### START
>> library(sp)
>> library(raster)
>> library(rasterVis)
>> library(RColorBrewer)
>> 
>> S <- raster::stack(system.file("external/rlogo.grd", package="raster"))
>> 
>> # make a polygon by shrinking the extent and casting object type
>> Poly <- as(raster::extent(S) + c(15, -32, 25, -10), "SpatialPolygons")
>> 
>> # build the layers in a list, adding the polygon to each layer
>> themes <- c("Reds", "Greens", "Blues")
>> PP <- lapply(seq_len(nlayers(S)),
>>         function(i) {
>>                 levelplot(S[[i]], 
>>                         par.settings = rasterVis::rasterTheme(region = RColorBrewer::brewer.pal(9, themes[i])),  
>>                         margin=FALSE) + 
>>                 layer(sp.polygons(Poly, col = "orange", lwd = 2))
>>                 }
>>         )
>>         
>> # print each layer, but specify the location within a layout scheme
>> print(PP[[1]], split=c(1, 1, 1, 3), more = TRUE)
>> print(PP[[2]], split=c(1, 2, 1, 3), more = TRUE)
>> print(PP[[3]], split=c(1, 3, 1, 3), more = FALSE)
>> ### END
>> 
>> 
>> 
>> > On Apr 11, 2019, at 7:28 AM, Kristi Glover <kristi.glover using hotmail.com> wrote:
>> > 
>> > He R users,
>> > I have been struggling to plot a boundary map over the raster maps. I tried using the following example, but the boundary map could not be displayed over the three raster maps.
>> > It works if we plot for a single raster. However when I want to plot the three maps using "levelplot" and add the boundary map it did not work.  I wanted to plot three raster same time because the "levelplot" so that we can compare the maps as they have only one legend.
>> > 
>> > My example code is given below, do you have any suggestions?
>> > 
>> > 
>> > library(gridExtra)
>> > 
>> > library(raster)
>> > 
>> > library(sp)
>> > 
>> > library(rasterVis)
>> > 
>> > library(rgdal)
>> > 
>> > library(maptools)
>> > 
>> > 
>> > boundary<- readShapeSpatial("boundrymap.shp")
>> > 
>> > 
>> > minTemp<-raster("minTemp.tif")
>> > 
>> > maxTemp<-raster("maxTemp.tif")
>> > 
>> > averageTemp<-raster("averageTemp.tif")
>> > 
>> > Temp<-stack(minTemp,maxTemp,averageTemp)
>> > 
>> > 
>> > levelplot(Temp, layout=c(1,3))+ layer(sp.polygons(boundary, col = "yellow"))
>> > 
>> > thanks
>> > 
>> > 
>> >        [[alternative HTML version deleted]]
>> > 
>> > ______________________________________________
>> > R-help using r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> > https://stat.ethz.ch/mailman/listinfo/r-help
>> > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> > and provide commented, minimal, self-contained, reproducible code.
>> 
>> Ben Tupper
>> Bigelow Laboratory for Ocean Sciences
>> 60 Bigelow Drive, P.O. Box 380
>> East Boothbay, Maine 04544
>> http://www.bigelow.org
>> 
>> Ecological Forecasting: https://eco.bigelow.org/
> 
> Ben Tupper
> Bigelow Laboratory for Ocean Sciences
> 60 Bigelow Drive, P.O. Box 380
> East Boothbay, Maine 04544
> http://www.bigelow.org
> 
> Ecological Forecasting: https://eco.bigelow.org/

Ben Tupper
Bigelow Laboratory for Ocean Sciences
60 Bigelow Drive, P.O. Box 380
East Boothbay, Maine 04544
http://www.bigelow.org

Ecological Forecasting: https://eco.bigelow.org/



More information about the R-help mailing list