[R-sig-Geo] plot overlay two spatialpolygons.

Michael Sumner mdsumner at gmail.com
Wed Aug 4 12:35:56 CEST 2010


Try add = TRUE with the second plot:

library(sp)
## plot of SpatialPolygonsDataFrame, using grey shades
library(maptools)
nc1 <- readShapePoly(system.file("shapes/sids.shp", package="maptools")
[1], proj4string=CRS("+proj=longlat +datum=NAD27"))
names(nc1)
rrt <- nc1$SID74/nc1$BIR74
brks <- quantile(rrt, seq(0,1,1/7))
cols <- grey((length(brks):2)/length(brks))
dens <- (2:length(brks))*3
plot(nc1, col=cols[findInterval(rrt, brks, all.inside=TRUE)])


## plot of SpatialPolygonsDataFrame, using line densities

nc <- readShapePoly(system.file("shapes/sids.shp", package="maptools")
[1], proj4string=CRS("+proj=longlat +datum=NAD27"))
names(nc)
rrt <- nc$SID74/nc$BIR74
brks <- quantile(rrt, seq(0,1,1/7))
cols <- grey((length(brks):2)/length(brks))
dens <- (2:length(brks))*3
plot(nc, density=dens[findInterval(rrt, brks, all.inside=TRUE)], add = TRUE)


On Wed, Aug 4, 2010 at 8:27 PM, Nikhil Kaza <nikhil.list at gmail.com> wrote:
> How can I overlay two spatialpolygons and plot them? The following
> produces two maps since a new plot is called. I understand that if it
> is any other type, we can simply call lines or points and overlay them
> on polygons.
>
> e.g.
> library(sp)
>
> ## plot of SpatialPolygonsDataFrame, using grey shades
> library(maptools)
> nc1 <- readShapePoly(system.file("shapes/sids.shp", package="maptools")
> [1], proj4string=CRS("+proj=longlat +datum=NAD27"))
> names(nc1)
> rrt <- nc1$SID74/nc1$BIR74
> brks <- quantile(rrt, seq(0,1,1/7))
> cols <- grey((length(brks):2)/length(brks))
> dens <- (2:length(brks))*3
> plot(nc1, col=cols[findInterval(rrt, brks, all.inside=TRUE)])
>
> library(sp)
>
> ## plot of SpatialPolygonsDataFrame, using line densities
> library(maptools)
> nc <- readShapePoly(system.file("shapes/sids.shp", package="maptools")
> [1], proj4string=CRS("+proj=longlat +datum=NAD27"))
> names(nc)
> rrt <- nc$SID74/nc$BIR74
> brks <- quantile(rrt, seq(0,1,1/7))
> cols <- grey((length(brks):2)/length(brks))
> dens <- (2:length(brks))*3
> plot(nc, density=dens[findInterval(rrt, brks, all.inside=TRUE)])
>
> Nikhil Kaza
> Asst. Professor,
> City and Regional Planning
> University of North Carolina
>
> nikhil.list at gmail.com
>
>
>        [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>



-- 
Michael Sumner
Institute for Marine and Antarctic Studies, University of Tasmania
Hobart, Australia
e-mail: mdsumner at gmail.com



More information about the R-sig-Geo mailing list