【发布时间】：2014-04-15 18:27:39
【问题描述】：
我正在尝试在底图上绘制属性的彩色网格。这里也提出了类似的问题，包括我最近问自己的一个问题，但早期的解决方案通常涉及ggplot2。我希望在不使用ggplot2 的情况下找到理想的解决方案，因为我最终的实际地图会非常复杂，并且到目前为止我使用plot 取得了实质性进展。
以下是绘制科罗拉多州并创建虚假属性数据网格的代码。我可以使用SpatialPoints 将网格转换为空间点并绘制这些点。但是，我怀疑我需要转换为空间多边形，可能使用SpatialPolygons。当我尝试使用SpatialPolygons 时，我得到如下所示的错误。
一旦我将网格单元图层添加到地图中，我想裁剪网格以仅在科罗拉多州内显示（或者，例如，如果我在地图上添加另一个州并使用更大的属性，则在科罗拉多州和怀俄明州内显示网格）。不过，这可能是一个后续问题。
这是一个较早的问题，其答案使用plot，但要求发布者提供更多信息：R Plot Filled Longitude-Latitude Grid Cells on Map
这是我之前提出的一个类似的问题，只有一个 ggplot2 解决方案：color grid cells in the United States and Canada
我下面的方法是基于我之前的另一个问题的解决方案：R: creating a map of selected Canadian provinces and U.S. states
感谢您的帮助。
library(rgdal)
library(maptools)
library(ggplot2)
library(plyr)
library(RColorBrewer)
library(classInt)
library(raster)
# NW(long,lat) SE(long,lat)
mapExtent <- rbind(c(-115, 43), c( -97, 35))
# assign projection
newProj <- CRS("+proj=longlat +datum=NAD83 +no_defs +ellps=GRS80 +towgs84=0,0,0")
## Project map extent
mapExtentPr <- spTransform(SpatialPoints(mapExtent,
proj4string=CRS("+proj=longlat")), newProj)
us1 <- getData('GADM', country="USA", level=1)
colorado <- us1[(us1$NAME_1 %in% c('Colorado')),]
## Project Colorada layer
coloradoPr <- spTransform( colorado, newProj)
# Create grid cells containing fake attribute data
# using approach found here:
# http://www.numbertheory.nl/2011/11/08/drawing-polar-centered-spatial-maps-using-ggplot2/
set.seed(1234)
xlim = c(-113, -99)
ylim = c( 35, 45)
dat.grid = expand.grid(x = xlim[1]:xlim[2], y = ylim[1]:ylim[2])
dat.grid$z = runif(nrow(dat.grid))
## Project grid layer
dat.gridPr <- spTransform(SpatialPoints( dat.grid, proj4string=CRS("+proj=longlat")), newProj)
dat.gridPr2 <- spTransform(SpatialPolygons(dat.grid, proj4string=CRS("+proj=longlat")), newProj)
#Error in spTransform(SpatialPolygons(dat.grid, proj4string = CRS("+proj=longlat")), :
# error in evaluating the argument 'x' in selecting a method for function 'spTransform': Error in SpatialPolygons(dat.grid, proj4string = CRS("+proj=longlat")) :
# cannot get a slot ("area") from an object of type "integer"
## Plot each projected layer, beginning with the projected extent
plot(mapExtentPr, pch=NA)
plot(coloradoPr , border="white", col="lightgrey", add=TRUE)
plot(dat.gridPr , border="white", col="lightgrey", add=TRUE)
plot(dat.gridPr , fill=dat.grid$z, add=TRUE)
plot(dat.gridPr , fill=dat.gridPr$z, add=TRUE)
编辑
这篇文章描述了如何一次添加一个网格单元，如果需要我可以这样做：
【问题讨论】：
标签： r map plot geospatial shapefile