9

R パッケージ ラスターを使用して、www.GADM.orgから世界地図データセットをインポートしました。マップのサイズを縮小するために作成したポリゴンにクリップしたいと思います。データを取得し、ポリゴンを問題なく作成できますが、「gIntersection」コマンドを使用すると、あいまいなエラー メッセージが表示されます。

世界地図データセットをクリップする方法について何か提案はありますか?

library(raster)
library(rgeos)

## Download Map of the World ##
WorldMap <- getData('countries')

## Create the clipping polygon
clip.extent <- as(extent(-20, 40, 30, 72), "SpatialPolygons")
proj4string(clip.extent) <- CRS(proj4string(WorldMap))

## Clip the map
EuropeMap <- gIntersection(WorldMap, clip.extent, byid = TRUE)

エラーメッセージ：

RGEOSBinTopoFunc(spgeom1, spgeom2, byid, id, "rgeos_intersection") のエラー: ジオメトリ コレクションには他のジオメトリ コレクションが含まれていない可能性があります さらに: 警告メッセージ: RGEOSBinTopoFunc(spgeom1, spgeom2, byid, id, "rgeos_intersection") 内: spgeom1 および spgeom2異なる proj4 文字列を持つ

4

2 に答える 2

2

How about a little intermediate step? I adopted the following code mainly from R-sig-Geo and I think it should do the trick. You'll need both 'maptools' and 'PBSmapping' packages, so make sure you've got them installed. Here's my code:

# Required packages
library(raster)
library(maptools)
library(PBSmapping)

# Download world map
WorldMap <- getData('countries')
# Convert SpatialPolygons to class 'PolySet'
WorldMap.ps <- SpatialPolygons2PolySet(WorldMap)
# Clip 'PolySet' by given extent
WorldMap.ps.clipped <- clipPolys(WorldMap.ps, xlim = c(-20, 40), ylim = c(30, 72))
# Convert clipped 'PolySet' back to SpatialPolygons
EuropeMap <- PolySet2SpatialPolygons(WorldMap.ps.clipped, close_polys=TRUE)

I just tested it and it worked without any problems. It took some computation time to transform SpatialPolygons to PolySet, though.

Cheers, Florian

于 2013-04-08T16:24:42.700 に答える