本文介绍了找出近似规则网格点的子集的周长的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧! 问题描述 让我们考虑一组2D近似规则的网格。这些网格与相邻网格相邻(相邻网格具有一个或多个相同的顶点)。以下是10个网格的顶点坐标(经度,纬度)如下所示: A lon lat [,1] [,2] [1,] 85.30754 27.91250 [2,] 85.32862 27.95735 [3,] 85.34622 27.89880 [4,] 85.36732 27.94364 [5,] 85.34958 28.00202 [6,] 85.38831 27.98830 [7,] 85.38487 27.88508 [8,] 85.40598 27.92991 [9,] 85.42353 27.87134 [10,] 85.44466 27.91616 [11,] 85.42698 27.97456 [12,] 85.46567 27.96081 [13,] 85.48334 27.90239 [14,] 85.50437 27.94703 [15,] 85.48645 28.00502 [16,] 85.52517 27.99123 [17,] 85.52198 27.88862 [18,] 85.54302 27.93325 [19,] 85.56384 27.97745 以上样本点的散点图(顶点): 网格是con结构如下图所示。 我的问题是如何获得周长(通过所有边界点的红色轮廓)?? 请注意:图中的红色圆圈点(1,3,7,9,10,13,17,18,19,16,15,12,11,6,5,2) 1是边界点。 注:可以观察到网格的边长小于6000米,所有网格的对角线长度 $ b 我在 geosphere distHaversine c $ c $> package函数在R中查找两点之间的距离。 解决方案大纲:所有对点的距离比6000米以网格形式形成图形。构建该图,然后查找所有子图同构于一个正方形(一个大小为4的循环)。由于外部边缘只是一个正方形的一部分(内部边缘由多个正方形共享),所以外部边缘看起来会少于内部边缘。因此,找到所有的内部边缘并放下它们,然后遍历结果图应该是一个简单的循环。 代码: library(igraph);库(geosphere) #主函数 perimeterGrid g = edgeP(makegrid(pts,maxdist = maxdist ,mindist = mindist)) loop = graph.dfs(minimum.spanning.tree(g),1)$ order cbind(V(g)$ x,V(g)$ y)[循环,] } #正弦矩阵距离矩阵 dmat< - 函数(pts){n = nrow(pts) do.call rbind,lapply(1:n,function(i){distHaversine(pts [i,],pts)}))} #使网格单元具有maxdist(和a心理学家停止自我边缘) makegrid dists = dmat(pts)g = graph.adjacency(dists< (g)$ x = pts [,1] V(g)$ y = pts [,2] } #执行网格边缘计数等的巧妙函数 edgeP< - 函数(g){#查找所有简单方块 square = graph.ring(4) subs = graph.get.subis omorphisms.vf2(g,square)#展开所有边 subs = do.call(rbind,lapply(subs,function(s){ rbind(s [1:2] ,s [2:3],s [3:4],s [c(4,1)])}))#绘制所有方格边缘的新图形e = graph.edgelist(subs,directed = FALSE)#添加权重作为边缘计数 E(e)$ weight = count.multiple(e) # (e)$ x = V(g)$ x V(e)$ y = V(g)$ y #删除来源 (e(e)$ weight == 256))e = e-edges(e(e)$ weight == 256)e = simplified(e) #内部边缘现在具有权重256。 b $ b#内部节点如何获得学位0 e = e - 顶点(度(e)== 0)返回(e)} 用法: plot(pts)多边形(perimeterGrid(pts),lwd = 2) 结果: 警告: 这是未经过测试的带网格的带孔洞的网格或者网格单元仅在单角点。也许这不可能发生。另外我不确定算法的效率如何,但看起来相当快。 Let us consider a set of near-regular grids in 2-D. These grids are adjacent (neighbouring grids have one or more same vertices) to the neighbouring grids. Here are the sample of 10 grids with the coordinates of the vertices (longitude,latitude) are as follows A<- lon lat [,1] [,2] [1,] 85.30754 27.91250 [2,] 85.32862 27.95735 [3,] 85.34622 27.89880 [4,] 85.36732 27.94364 [5,] 85.34958 28.00202 [6,] 85.38831 27.98830 [7,] 85.38487 27.88508 [8,] 85.40598 27.92991 [9,] 85.42353 27.87134 [10,] 85.44466 27.91616 [11,] 85.42698 27.97456 [12,] 85.46567 27.96081 [13,] 85.48334 27.90239 [14,] 85.50437 27.94703 [15,] 85.48645 28.00502 [16,] 85.52517 27.99123 [17,] 85.52198 27.88862 [18,] 85.54302 27.93325 [19,] 85.56384 27.97745The scatter-plot of the above sample set of points (vertices):The grids are constructed as in the following picture.My question is how to get the perimeter (red contour passing through all boundary points)??Note that: The red encircled points (1,3,7,9,10,13,17,18,19,16,15,12,11,6,5,2) in figure 1 are the boundary points.Note: It is observed the sides of the grids are less than 6000 metres and length of diagonals of all grids are more than 6000 metres.I am using distHaversine from the geosphere package function in R to find the distance between two points. 解决方案 In outline: all pairs of points closer than 6000m form a graph in the form of grid squares. Construct that graph, and then find all subgraphs isomorphic to a square (a loop of size 4). The external edges will appear less often than internal edges since they are only part of one square (internal edges are shared by multiple squares). Hence find all the internal edges and drop them, then traverse the resulting graph which should be a simple loop.Code:library(igraph); library(geosphere)# main functionperimeterGrid <- function(pts, maxdist=6000, mindist=1){ g = edgeP(makegrid(pts, maxdist=maxdist, mindist=mindist)) loop = graph.dfs(minimum.spanning.tree(g),1)$order cbind(V(g)$x, V(g)$y)[loop,]}# haversine distance matrixdmat <- function(pts){ n=nrow(pts) do.call(rbind,lapply(1:n,function(i){distHaversine(pts[i,],pts)}))}# make the grid cells given a maxdist (and a mindist to stop self-self edges) makegrid <- function(pts, maxdist=6000, mindist=1){ dists = dmat(pts) g = graph.adjacency(dists<maxdist & dists>mindist, mode="undirected") V(g)$x=pts[,1] V(g)$y=pts[,2] g}# clever function that does the grid edge count etcedgeP <- function(g){ # find all the simple squares square=graph.ring(4) subs = graph.get.subisomorphisms.vf2(g,square) # expand all the edges subs = do.call(rbind, lapply(subs, function(s){ rbind(s[1:2], s[2:3], s[3:4], s[c(4,1)]) })) # make new graph of the edges of all the squares e = graph.edgelist(subs,directed=FALSE) # add the weight as the edge count E(e)$weight=count.multiple(e) # copy the coords from the source V(e)$x=V(g)$x V(e)$y=V(g)$y # remove multiple edges e=simplify(e) # internal edges now have weight 256. e = e - edges(which(E(e)$weight==256)) # internal nodes how have degree 0 e = e - vertices(degree(e)==0) return(e)}Usage: plot(pts) polygon(perimeterGrid(pts),lwd=2)Results:Warnings:This is untested on grid fragments with holes or where grid cells are only touching at a single corner points. Maybe that can't happen. Also I'm not sure what the efficiency of the algorithms are, but it seems pretty quick. 这篇关于找出近似规则网格点的子集的周长的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持!
10-29 09:50