从嵌套的 SpatialPolygonsDataFrame 中提取父多边形或从父多边形中提取 'Dissolving' 个孔
Extracting the parent polygon from a nested SpatialPolygonsDataFrame or 'Dissolving' holes from a parent ploygon
编辑 经过更多研究但仍然没有解决方案,我正在添加大量编辑以及 link 到 .shp 文件。
The shape file is included here
我有一个包含 9 个多边形的 SpatialPolygonsDataFrame,每个多边形还包含多个嵌套多边形 - 'holes'。数据摘要在这里。
> summary(data)
Object of class SpatialPolygonsDataFrame
Coordinates:
min max
x 483298.9 643204.4
y 4782172.1 4997248.3
Is projected: TRUE
proj4string :
[+proj=utm +zone=12 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0]
Data attributes:
Id IndID
Min. :0 BHS_011_A:1
1st Qu.:0 BHS_015_A:1
Median :0 BHS_083_A:1
Mean :0 BHS_089_A:1
3rd Qu.:0 BHS_091_A:1
Max. :0 BHS_129_A:1
(Other) :3
数据 structure
的示例如下。
Formal class 'SpatialPolygonsDataFrame' [package "sp"] with 5 slots
..@ data :'data.frame': 9 obs. of 2 variables:
.. ..$ Id : int [1:9] 0 0 0 0 0 0 0 0 0
.. ..$ IndID: Factor w/ 9 levels "BHS_011_A","BHS_015_A",..: 1 2 3 4 5 6 7 8 9
..@ polygons :List of 9
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 5
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 513497 4986246
.. .. .. .. .. .. ..@ area : num 76614017
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:287, 1:2] 509244 507384 507214 507010 506899 ...
.. .. .. .. .. .. .. ..- attr(*, "dimnames")=List of 2
.. .. .. .. .. .. .. .. ..$ : NULL
.. .. .. .. .. .. .. .. ..$ : chr [1:2] "x" "y"
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 509678 4979511
.. .. .. .. .. .. ..@ area : num 1462398
.. .. .. .. .. .. ..@ hole : logi TRUE
.. .. .. .. .. .. ..@ ringDir: int -1
.. .. .. .. .. .. ..@ coords : num [1:7, 1:2] 509301 509269 509194 509007 509412 ...
.. .. .. .. .. .. .. ..- attr(*, "dimnames")=List of 2
.. .. .. .. .. .. .. .. ..$ : NULL
.. .. .. .. .. .. .. .. ..$ : chr [1:2] "x" "y"
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 515572 4988493
.. .. .. .. .. .. ..@ area : num 1579348
.. .. .. .. .. .. ..@ hole : logi TRUE
.. .. .. .. .. .. ..@ ringDir: int -1
.. .. .. .. .. .. ..@ coords : num [1:10, 1:2] 514520 514570 514684 516501 515996 ...
.. .. .. .. .. .. .. ..- attr(*, "dimnames")=List of 2
.. .. .. .. .. .. .. .. ..$ : NULL
.. .. .. .. .. .. .. .. ..$ : chr [1:2] "x" "y"
如下例(九个之一)所示,父多边形有多个孔。
data <- readOGR(".", "IndLineBuff")
plot(data[data$IndID == "MTG_005_A",])
这是我第一次涉足 sp()
、rgdal()
、rgeos()
和其他空间包,我发现了很多关于使用运算符提取的有用帖子区域等,但问题仍然存在。虽然 this post 提供了接近的解决方案,但我似乎无法调整代码以满足此处描述的需求。
我想获得一个 SpatialPolygonsDataFrame,它只包含来自每组子列表(即 data@polygon)的父(最大)多边形。看来我应该只能提取父多边形或 'dissolve' 孔。
最终结果是 9 个多边形,每个都是 9 个列表的父项,我可以将其导出为 ESRI shapefile。
如有任何建议,我们将不胜感激。
以下 hack 可能有效:
data@polygons = lapply(data@polygons, function(x) { x@Polygons = x@Polygons[1]; x})
plot(data[data$IndID == "MTG_005_A",])
在你的情况下,你似乎已经缓冲了线性特征,并希望消除由此产生的多边形中的孔洞。 9 个多边形中的每一个似乎都由主多边形和后跟一组整体组成。该函数选择此列表中的第一个,返回所需的列表列表。要理解 class 结构,请仔细研究 ?"SpatialPolygons-class"
、?"Polygons-class"
和 ?"Polygon-class"
。
当我第一次学习 sp 和相关包时,我发现 Barry Rowlingson 的网站非常有用。
它包括:
作弊sheet
非常有用的参考 sheet,其中包含有关如何做事的提示
使用 R 包\
http://www.maths.lancs.ac.uk/~rowlings/Teaching/UseR2012/cheatsheet.html
Spatial 简介
spatialDataframe 结构的解释
包 sp。还解释了 CRS(坐标参考系统)和
投影(proj4string 和 proj4.Project)。最后还有
涵盖光栅图像的使用和操作,以及提取
来自他们的数据使用 R\
http://www.maths.lancs.ac.uk/~rowlings/Teaching/UseR2012/static/talk1.pdf
回答你的问题,这很混乱,因为你正在处理一连串的事情,但这是可以做到的。
以下代码似乎有效(我已经在我已有的 shapefile 上测试了它(位),而不是下载你的)。如果在读入后将 NSWACTA
替换为 shapefile 的名称,它 应该 适合您的情况。虽然没有承诺 :-)。
第一部分是介绍性的,展示了 类 的所有列表的情况。希望这可以帮助您弄清楚发生了什么。它是 dropHole
函数和它的方法,它实际上做你想做的事。
class(NSWACTA)
# NSWACTA has class sp::SpatialPolygonsDataFrame
slotNames(NSWACTA)
# The SpatialPolygonsDataFrame has a slot called @polygons
class(NSWACTA@polygons)
# The @polgyons slot contains a list
class(NSWACTA@polygons[[1]])
# Elements of that list are of class sp::Polygons
slotNames(NSWACTA@polygons[[1]])
# An sp::Polygons class has several slots, one of which is @Polygons
class(NSWACTA@polygons[[1]]@Polygons)
# The @Polygons slot contains a list
class(NSWACTA@polygons[[1]]@Polygons[[1]])
# The entries in the @Polygons list are of class sp::Polygon
slotNames(NSWACTA@polygons[[1]]@Polygons[[1]])
# An sp::Polygon has several slots.
# @coords holds the coordinates which define the boundary.
# @hole is a logical field indicating whether or not the sp::Polygon is a hole
setGeneric('dropHole',def=function(poly, ...){
standardGeneric('dropHole')
})
# Drop a single sp::Polygon if it is holey
setMethod('dropHole',signature=signature('Polygon'),
def=function(poly) {
#return only Polygons which are not holes
if (poly@hole) NULL else poly
}
)
# Drop holey sp::Polygon entries in the @Polygons list of an sp::Polygons class
setMethod('dropHole', signature = signature('Polygons'),
def = function(poly) {
noHoles <- lapply(poly@Polygons, dropHole)
#Remove the NULL entries from the list
noHoles <- Filter(Negate(is.null), noHoles)
# Turn back into a (single) Polygons
# The generator function (sp::Polygons) fills in the other slots!
# return the new sp:Polygons object
sp::Polygons(noHoles, ID = poly@ID)
}
)
# Drop holey parts of sp::Polygons in the @polygons list
# of an sp::SpatialPolygonsDataFrame
setMethod('dropHole', signature = signature('SpatialPolygonsDataFrame'),
def = function(poly) {
noHoles <- lapply(poly@polygons, dropHole)
# Put the un holey Polygons list back into the @polygons slot
poly@polygons <- noHoles
#return the modified SpatialPolygonsDataFrame
poly
}
)
newmap <- dropHole(NSWACTA)
这是第一部分的输出。 注意注意名字大小写、单复数的区别! (即多边形、多边形和多边形)
d> class(NSWACTA)
[1] "SpatialPolygonsDataFrame"
attr(,"package")
[1] "sp"
d> # NSWACTA has class sp::SpatialPolygonsDataFrame
d> slotNames(NSWACTA)
[1] "data" "polygons" "plotOrder" "bbox" "proj4string"
d> # The SpatialPolygonsDataFrame has a slot called @polygons
d> class(NSWACTA@polygons)
[1] "list"
d> # The @polgyons slot contains a list
d> class(NSWACTA@polygons[[1]])
[1] "Polygons"
attr(,"package")
[1] "sp"
d> # Elements of that list are of class sp::Polygons
d> slotNames(NSWACTA@polygons[[1]])
[1] "Polygons" "plotOrder" "labpt" "ID" "area"
d> # An sp::Polygons class has several slots, one of which is @Polygons
d> class(NSWACTA@polygons[[1]]@Polygons)
[1] "list"
d> # The @Polygons slot contains a list
d> class(NSWACTA@polygons[[1]]@Polygons[[1]])
[1] "Polygon"
attr(,"package")
[1] "sp"
d> # The entries in the @Polygons list are of class sp::Polygon
d> slotNames(NSWACTA@polygons[[1]]@Polygons[[1]])
[1] "labpt" "area" "hole" "ringDir" "coords"
d> # An sp::Polygon has several slots.
d> # @coords holds the coordinates which define the boundary.
d> # @hole is a logical field indicating whether or not the sp::Polygon is a hole
编辑 经过更多研究但仍然没有解决方案,我正在添加大量编辑以及 link 到 .shp 文件。
The shape file is included here
我有一个包含 9 个多边形的 SpatialPolygonsDataFrame,每个多边形还包含多个嵌套多边形 - 'holes'。数据摘要在这里。
> summary(data)
Object of class SpatialPolygonsDataFrame
Coordinates:
min max
x 483298.9 643204.4
y 4782172.1 4997248.3
Is projected: TRUE
proj4string :
[+proj=utm +zone=12 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0]
Data attributes:
Id IndID
Min. :0 BHS_011_A:1
1st Qu.:0 BHS_015_A:1
Median :0 BHS_083_A:1
Mean :0 BHS_089_A:1
3rd Qu.:0 BHS_091_A:1
Max. :0 BHS_129_A:1
(Other) :3
数据 structure
的示例如下。
Formal class 'SpatialPolygonsDataFrame' [package "sp"] with 5 slots
..@ data :'data.frame': 9 obs. of 2 variables:
.. ..$ Id : int [1:9] 0 0 0 0 0 0 0 0 0
.. ..$ IndID: Factor w/ 9 levels "BHS_011_A","BHS_015_A",..: 1 2 3 4 5 6 7 8 9
..@ polygons :List of 9
.. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. .. .. ..@ Polygons :List of 5
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 513497 4986246
.. .. .. .. .. .. ..@ area : num 76614017
.. .. .. .. .. .. ..@ hole : logi FALSE
.. .. .. .. .. .. ..@ ringDir: int 1
.. .. .. .. .. .. ..@ coords : num [1:287, 1:2] 509244 507384 507214 507010 506899 ...
.. .. .. .. .. .. .. ..- attr(*, "dimnames")=List of 2
.. .. .. .. .. .. .. .. ..$ : NULL
.. .. .. .. .. .. .. .. ..$ : chr [1:2] "x" "y"
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 509678 4979511
.. .. .. .. .. .. ..@ area : num 1462398
.. .. .. .. .. .. ..@ hole : logi TRUE
.. .. .. .. .. .. ..@ ringDir: int -1
.. .. .. .. .. .. ..@ coords : num [1:7, 1:2] 509301 509269 509194 509007 509412 ...
.. .. .. .. .. .. .. ..- attr(*, "dimnames")=List of 2
.. .. .. .. .. .. .. .. ..$ : NULL
.. .. .. .. .. .. .. .. ..$ : chr [1:2] "x" "y"
.. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. .. .. ..@ labpt : num [1:2] 515572 4988493
.. .. .. .. .. .. ..@ area : num 1579348
.. .. .. .. .. .. ..@ hole : logi TRUE
.. .. .. .. .. .. ..@ ringDir: int -1
.. .. .. .. .. .. ..@ coords : num [1:10, 1:2] 514520 514570 514684 516501 515996 ...
.. .. .. .. .. .. .. ..- attr(*, "dimnames")=List of 2
.. .. .. .. .. .. .. .. ..$ : NULL
.. .. .. .. .. .. .. .. ..$ : chr [1:2] "x" "y"
如下例(九个之一)所示,父多边形有多个孔。
data <- readOGR(".", "IndLineBuff")
plot(data[data$IndID == "MTG_005_A",])
这是我第一次涉足 sp()
、rgdal()
、rgeos()
和其他空间包,我发现了很多关于使用运算符提取的有用帖子区域等,但问题仍然存在。虽然 this post 提供了接近的解决方案,但我似乎无法调整代码以满足此处描述的需求。
我想获得一个 SpatialPolygonsDataFrame,它只包含来自每组子列表(即 data@polygon)的父(最大)多边形。看来我应该只能提取父多边形或 'dissolve' 孔。
最终结果是 9 个多边形,每个都是 9 个列表的父项,我可以将其导出为 ESRI shapefile。
如有任何建议,我们将不胜感激。
以下 hack 可能有效:
data@polygons = lapply(data@polygons, function(x) { x@Polygons = x@Polygons[1]; x})
plot(data[data$IndID == "MTG_005_A",])
在你的情况下,你似乎已经缓冲了线性特征,并希望消除由此产生的多边形中的孔洞。 9 个多边形中的每一个似乎都由主多边形和后跟一组整体组成。该函数选择此列表中的第一个,返回所需的列表列表。要理解 class 结构,请仔细研究 ?"SpatialPolygons-class"
、?"Polygons-class"
和 ?"Polygon-class"
。
当我第一次学习 sp 和相关包时,我发现 Barry Rowlingson 的网站非常有用。
它包括:
作弊sheet
非常有用的参考 sheet,其中包含有关如何做事的提示 使用 R 包\ http://www.maths.lancs.ac.uk/~rowlings/Teaching/UseR2012/cheatsheet.htmlSpatial 简介
spatialDataframe 结构的解释 包 sp。还解释了 CRS(坐标参考系统)和 投影(proj4string 和 proj4.Project)。最后还有 涵盖光栅图像的使用和操作,以及提取 来自他们的数据使用 R\ http://www.maths.lancs.ac.uk/~rowlings/Teaching/UseR2012/static/talk1.pdf
回答你的问题,这很混乱,因为你正在处理一连串的事情,但这是可以做到的。
以下代码似乎有效(我已经在我已有的 shapefile 上测试了它(位),而不是下载你的)。如果在读入后将 NSWACTA
替换为 shapefile 的名称,它 应该 适合您的情况。虽然没有承诺 :-)。
第一部分是介绍性的,展示了 类 的所有列表的情况。希望这可以帮助您弄清楚发生了什么。它是 dropHole
函数和它的方法,它实际上做你想做的事。
class(NSWACTA)
# NSWACTA has class sp::SpatialPolygonsDataFrame
slotNames(NSWACTA)
# The SpatialPolygonsDataFrame has a slot called @polygons
class(NSWACTA@polygons)
# The @polgyons slot contains a list
class(NSWACTA@polygons[[1]])
# Elements of that list are of class sp::Polygons
slotNames(NSWACTA@polygons[[1]])
# An sp::Polygons class has several slots, one of which is @Polygons
class(NSWACTA@polygons[[1]]@Polygons)
# The @Polygons slot contains a list
class(NSWACTA@polygons[[1]]@Polygons[[1]])
# The entries in the @Polygons list are of class sp::Polygon
slotNames(NSWACTA@polygons[[1]]@Polygons[[1]])
# An sp::Polygon has several slots.
# @coords holds the coordinates which define the boundary.
# @hole is a logical field indicating whether or not the sp::Polygon is a hole
setGeneric('dropHole',def=function(poly, ...){
standardGeneric('dropHole')
})
# Drop a single sp::Polygon if it is holey
setMethod('dropHole',signature=signature('Polygon'),
def=function(poly) {
#return only Polygons which are not holes
if (poly@hole) NULL else poly
}
)
# Drop holey sp::Polygon entries in the @Polygons list of an sp::Polygons class
setMethod('dropHole', signature = signature('Polygons'),
def = function(poly) {
noHoles <- lapply(poly@Polygons, dropHole)
#Remove the NULL entries from the list
noHoles <- Filter(Negate(is.null), noHoles)
# Turn back into a (single) Polygons
# The generator function (sp::Polygons) fills in the other slots!
# return the new sp:Polygons object
sp::Polygons(noHoles, ID = poly@ID)
}
)
# Drop holey parts of sp::Polygons in the @polygons list
# of an sp::SpatialPolygonsDataFrame
setMethod('dropHole', signature = signature('SpatialPolygonsDataFrame'),
def = function(poly) {
noHoles <- lapply(poly@polygons, dropHole)
# Put the un holey Polygons list back into the @polygons slot
poly@polygons <- noHoles
#return the modified SpatialPolygonsDataFrame
poly
}
)
newmap <- dropHole(NSWACTA)
这是第一部分的输出。 注意注意名字大小写、单复数的区别! (即多边形、多边形和多边形)
d> class(NSWACTA)
[1] "SpatialPolygonsDataFrame"
attr(,"package")
[1] "sp"
d> # NSWACTA has class sp::SpatialPolygonsDataFrame
d> slotNames(NSWACTA)
[1] "data" "polygons" "plotOrder" "bbox" "proj4string"
d> # The SpatialPolygonsDataFrame has a slot called @polygons
d> class(NSWACTA@polygons)
[1] "list"
d> # The @polgyons slot contains a list
d> class(NSWACTA@polygons[[1]])
[1] "Polygons"
attr(,"package")
[1] "sp"
d> # Elements of that list are of class sp::Polygons
d> slotNames(NSWACTA@polygons[[1]])
[1] "Polygons" "plotOrder" "labpt" "ID" "area"
d> # An sp::Polygons class has several slots, one of which is @Polygons
d> class(NSWACTA@polygons[[1]]@Polygons)
[1] "list"
d> # The @Polygons slot contains a list
d> class(NSWACTA@polygons[[1]]@Polygons[[1]])
[1] "Polygon"
attr(,"package")
[1] "sp"
d> # The entries in the @Polygons list are of class sp::Polygon
d> slotNames(NSWACTA@polygons[[1]]@Polygons[[1]])
[1] "labpt" "area" "hole" "ringDir" "coords"
d> # An sp::Polygon has several slots.
d> # @coords holds the coordinates which define the boundary.
d> # @hole is a logical field indicating whether or not the sp::Polygon is a hole