当我使颜色透明时，SpatialPolygonsDataFrame 无法正确绘制

Question

这几天我一直在努力解决这个问题。我有一个包含邮政编码多边形的 shapefile。我正在尝试用透明颜色 (alpha=0.7) 绘制其中的某个部分，但绘图边缘的一些多边形最终没有颜色。如果alpha=1.0，所有多边形都可以。 shapefile 存储为大SpatialPolygonsDataFrame。

# Code that you can use (Thanks jbaums)

library(rgdal);
download.file('http://www.naturalearthdata.com/http//www.naturalearthdata.com/do‌​wnload/110m/cultural/ne_110m_admin_0_countries.zip', {f <- tempfile()}); 
unzip(f, exdir=tempdir()); 
shp <- readOGR(tempdir(), 'ne_110m_admin_0_countries'); 
plot(shp, col='#FF000070', xlim=c(-100, 25))

作为参考，我的代码如下所示：

#~~~ Read zip shapefile
  in_zip.shp  = "./Data/USZIP11_wSTCntyLLFIPS.shp"
  zip_poly = readShapePoly(in_zip.shp)

#~~~ Set-up the mapping parameters
  xlim = c(-82.0, -80.7)
  ylim = c(25.5, 26.8)

#~~~ Works well
  output_figure <- "./Good.png"
  png(output_figure, width=5, height=4, units="in", res=1200)
  par(mar=c(1.5, 1.8, 2.0, 0.7))
  plot(xlim, ylim, type="n", xlim=xlim, ylim=ylim, axes=FALSE, xlab="", ylab="")

  plot(zip_poly, col=rgb(1, 0, 0, alpha=1), xlim=xlim, ylim=ylim, lty=1.3, add=T)
  box()
  title(main="Good")
  dev.off()

#~~~ Doesn't work well
  output_figure <- "./Bad.png"
  png(output_figure, width=5, height=4, units="in", res=1200)
  par(mar=c(1.5, 1.8, 2.0, 0.7))
  plot(xlim, ylim, type="n", xlim=xlim, ylim=ylim, axes=FALSE, xlab="", ylab="")

  plot(zip_poly,col=rgb(1,0,0,alpha=0.7),xlim=xlim,ylim=ylim,lty=1.3,add=T)
  box()
  title(main="Bad")
  dev.off()

任何提示将不胜感激。我附上了数字，这样你就可以看到我的代码有什么问题。谢谢！

Answer 1

首先，我不知道问题的原因。它适用于常规颜色但不适用于具有 alpha 通道的颜色，这对我来说似乎很奇怪。

其次，这里是使用rgeos 包中的gIntersection 将多边形裁剪为感兴趣的大小（即xlim 和ylim）的变通方法示例。绘制裁剪的地图似乎适用于 alpha 颜色。

library(rgdal)
library(rgeos)  # for gIntersection
library(raster) # for extent()

# download world map
download.file('http://www.naturalearthdata.com/http//www.naturalearthdata.com/download/110m/cultural/ne_110m_admin_0_countries.zip', {f <- tempfile()}); 
unzip(f, exdir=tempdir()); 
shp <- readOGR(tempdir(), 'ne_110m_admin_0_countries'); 

# define map limits
xlim = c(0, 90)
ylim = c(40, 60)

# red with alpha
mycol = rgb(1,0,0,0.5)

# plot world map
plot(shp)

# create rectangle from xlim and ylim
CP <- as(extent(c(xlim,ylim)), "SpatialPolygons")
proj4string(CP) <- proj4string(shp)

# add to plot
plot(CP, add=T, col=rgb(1,0,0,0.5))

# set up window for 2 plots, comparing 2 methods
par(mfrow=c(2,1))

# plot map with xlim and ylim: does not work the way we want
plot(shp, xlim=xlim, ylim=ylim, col=mycol, axes=T, las=1)

# clip map to xlim and ylim, then plot: works
out <- gIntersection(shp, CP, byid=TRUE)
plot(out, col=mycol, axes=T, las=1)