r - Extracting the outer-boundary of a set of grid points with some missing grids (holes) -

this question generalization of this question. mentioned question working point set no hole. in present question want perimeter (outer boundary) of subset of near-regular grid of points of grid point in polygon missing (i.e., polygon hole).

the sample data set on grids available here.

i used r-code suggested answer in above mentioned question (with no holes).

the following output of using codes:

now want ignore holes in inside point set , want consider outer boundary of point set required polygon.

any suggestion!! thanks.

this slight variant on previous code works if there's holes finding loops , taking 1 largest x coordinate, must outside loop. unless loops touch... strictly perhaps should take loop largest area... note need use x , y in 1 of functions because of bug (i've reported) in igraph package.

perimetergrid <- function(pts, maxdist=6000, mindist=1){     g = edgep(makegrid(pts, maxdist=maxdist, mindist=mindist))      ## there might holes. find loop 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 grid cells given maxdist (and mindist 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 , y here not x , y due igraph bug     ## these copied output later...     v(g)$x=pts[,1]     v(g)$y=pts[,2]      g }  # clever function grid edge count etc edgep <- function(g){     # find simple squares     square=graph.ring(4)     subs = graph.get.subisomorphisms.vf2(g,square)     # expand 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 edges of squares     e = graph.edgelist(subs,directed=false)     # add weight edge count     e(e)$weight=count.multiple(e)      # copy coords source     v(e)$x=v(g)$x     v(e)$y=v(g)$y      # remove multiple edges     e=simplify(e)      # internal edges 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) }