Ns

gbufferを使用してRの(地理)空間ポイントをバッファリングする



Buffer Spatial Points R With Gbuffer



解決:

もう少し掘り下げてみると、「投影された」座標参照系の使用は次のように簡単であることがわかります。

#Statscan CRSを取得するには、ここを参照してください:#http://spatialreference.org/ref/epsg/3347/ pc<- spTransform( sampledf, CRS( '+init=epsg:3347' ) )  

STATSCAN(カナダの住所に適しています)で使用されるEPSG3347は、ランベルト正角円錐図法を使用します。 NAD83は不適切であることに注意してください。これは、「投影された」CRSではなく「地理的な」CRSです。ポイントをバッファリングするには



pc100km<- gBuffer( pc, width=100*distm, byid=TRUE ) # Add data, and write to shapefile pc100km <- SpatialPolygonsDataFrame( pc100km, [email protected] ) writeOGR( pc100km, 'pc100km', 'pc100km', driver='ESRI Shapefile' )   

@MichaelChiricoが指摘したように、使用するデータを投影するrgeos :: gBuffer()は注意して適用する必要があります。私は測地学の専門家ではありませんが、このESRIの記事(測地線バッファリングについて)から理解した限りでは、投影してから適用しますgBufferは実際に生産することを意味します ユークリッド とは対照的にバッファ 測地線 もの。ユークリッドバッファは、投影された座標系によって導入された歪みの影響を受けます。これらの歪みは、分析に広いバッファーが含まれている場合、特に広い範囲の緯度の範囲が広い場合に心配になる可能性があります(カナダが適切な候補であると思います)。

しばらく前に同じ問題に遭遇し、gis.stackexchange(Rのユークリッドおよび測地線バッファリング)に質問を向けました。そのときに提案したRコードと、与えられた回答は、ここでもこの質問に関連していると思います。



主なアイデアはを利用することですgeosphere :: destPoint()。詳細とより高速な代替手段については、上記のgis.stackexchangeリンクを参照してください。これがあなたの2つのポイントに適用された私の古い試みです:

library(geosphere)library(sp)pts lon lat#> 1 -53.20198 47.05564#> 2 -52.81218 47.31741 make_GeodesicBuffer<- function(pts, width) { # A) Construct buffers as points at given distance and bearing --------------- dg <- seq(from = 0, to = 360, by = 5) # Construct equidistant points defining circle shapes (the 'buffer points') buff.XY <- geosphere::destPoint(p = pts, b = rep(dg, each = length(pts)), d = width) # B) Make SpatialPolygons ------------------------------------------------- # Group (split) 'buffer points' by id buff.XY <- as.data.frame(buff.XY) id <- rep(1:dim(pts)[1], times = length(dg)) lst <- split(buff.XY, id) # Make SpatialPolygons out of the list of coordinates poly <- lapply(lst, sp::Polygon, hole = FALSE) polys <- lapply(list(poly), sp::Polygons, ID = NA) spolys <- sp::SpatialPolygons(Srl = polys, proj4string = CRS('+proj=longlat +ellps=WGS84 +datum=WGS84')) # Disaggregate (split in unique polygons) spolys <- sp::disaggregate(spolys) return(spolys) } pts_buf_100km <- make_GeodesicBuffer(as.matrix(pts), width = 100*10^3) # Make a kml file and check the results on Google Earth library(plotKML) #>plotKMLバージョン0.5-9(2019-01-04)#> URL:http://plotkml.r-forge.r-project.org/ kml(pts_buf_100km、file.name = 'pts_buf_100km.kml')#> KMLファイル書き込み用に開かれました...#> KMLへの書き込み...#>終了pts_buf_100km.kml

reprexパッケージ(v0.2.1)によって2019-02-11に作成されました

ここに画像の説明を入力してください




そして、おもちゃで関数をパッケージにラップしました-geobuffer

次に例を示します。

#install.packages( 'devtools')#devtoolsがない場合は、インストールしますdevtools :: install_github( 'valentinitnelav / geobuffer')library(geobuffer)pts<- data.frame(lon = c(-53.20198, -52.81218), lat = c(47.05564, 47.31741)) pts_buf_100km <- geobuffer_pts(xy = pts, dist_m = 100*10^3)  

reprexパッケージ(v0.2.1)によって2019-02-11に作成されました

他の人はより良い解決策を思い付くかもしれませんが、今のところ、これは私の問題に対してうまく機能し、うまくいけば他の問題も解決できるでしょう。