提取一些缺少網格（孔）的網格點的外邊界

Question

這個問題是關於this問題的泛化。上面提到的問題對於沒有漏洞的點來說效果很好。在目前的問題中，我希望獲得缺少多邊形中某些網格點的點（即帶有孔的多邊形）的近似規則網格子集的邊界（外邊界）。

網格上的樣本數據集可用here。

我在上面提到的問題中使用了R代碼作爲答案（沒有漏洞）。

以下是使用這些代碼的輸出：

現在我想它忽略設置點內孔要考慮設置爲所需的多邊形點的外邊界。

任何建議！謝謝。

Answer 1

如果通過查找所有循環並僅取得具有最大X座標的一個循環（這必須是外部循環），那麼以前代碼的這種輕微變體纔有效。除非循環觸及...嚴格地說，它可能需要佔用最大面積的循環......注意，由於igraph軟件包中存在一個錯誤（我已經報道），所以需要在其中一個函數中使用X和Y.

perimeterGrid <- function(pts, maxdist=6000, mindist=1){ 
    g = edgeP(makegrid(pts, maxdist=maxdist, mindist=mindist)) 

    ## there might be holes. Find the loop with the largest X coordinate. 
    parts = components(g) 
    outer = which.max(tapply(V(g)$x, parts$membership, function(x){max(x)})) 

    g = induced_subgraph(g, which(parts$membership==outer)) 

    loop = graph.dfs(minimum.spanning.tree(g),1)$order 
    cbind(V(g)$x, V(g)$y)[loop,] 
} 

# haversine distance matrix 
dmat <- 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") 
    ## use X and Y here not x and y due to igraph bug 
    ## these get copied to the output later... 
    V(g)$X=pts[,1] 
    V(g)$Y=pts[,2] 

    g 
} 

# clever function that does the grid edge count etc 
edgeP <- 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 back 
    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) 
}