sc_neighbours fn for TRI

[mpadge]

It'll be necessary to implement a function to extract neighbours from a TRI object, or maybe even include them as part of it?

Also some work to do in TRI-land:

> rbenchmark::benchmark (
                       TRI(inlandwaters),
                       nbs <- sc_coord (inlandwaters) %>%
                           dplyr::distinct (x_, y_) %>%
                           dplyr::rename (x = x_, y = y_) %>%
                           tripack::tri.mesh () %>%
                           tripack::neighbours (),
                       replications = 5)$relative
# [1] 4.5  1.0

or something like that.

[mdsumner]

tripack doesn't respect input edges, and since decido is path based it must work per feature, but I'll see if it can be faster :)

Decido is now on CRAN!

Some neighbour stuff here needs update from old PRIMITIVE model

https://rpubs.com/cyclemumner/neighbours02

And thanks, good idea!

[mdsumner]

Actually, do you mean object neighbours or individual triangle neighbours? What's your immediate use-case?

The first is a summary of object_link_triangle, the second requires conversion to SC and then summary of object_link_edge (possibly also while carrying the triangle id). I think the first is straightforward, but the second starts to imply the need to "activate" a different table(in tidygraph sense), with object being an obvious default.

I'm not sure that activate makes sense for these extractor verbs, usually it would imply an operation on the whole graph IIUC.

I think we do need some kind of union-class of SC and TRI, but it's not obvious to me yet. It might be the true super class SC, the simplicial complex - and the current SC is more of an EDGE (that SC extends along with TRI), and then still a two-table UNIVERSAL or something ... A problem is that we can't triangulate a set of edges without RTriangle via anglr::DEL.

[mpadge]

great that decido is on CRAN! Nice work! I'll have a ponder of the other stuff tomorrow and get back to you

[mdsumner]

Whoa, tripack does allow adding constraints. I guess I didn't understand the topic well enough first time around. Ah but license is restrictive also...

[mpadge]

I just mean a list of all points in all triangles which include a given point. That would translate in TRI$vertex to finding all vertices which occur in same row as a given vertex. For the blue point in the image, that's all green points. The return from tripack::neighbours is, imho, supremely useful and simple. Something like that would be great. (And the fig is coz I'm building custom spatially-constrained spanning trees based on maximal covariance between points.)

[mdsumner]

oh, tripack "constrained delaunay" just means a boundary around it to replace the convex hull, it's not a constraint on the mesh ...

I'll take on their neighbours approach

[mdsumner]

Is this best done as SC::TRI? (that was recently added in mike-dev)

library(silicate)
library(dplyr)
x <- SC(TRI(minimal_mesh))
neigh0 <- do.call(rbind, purrr::map(purrr::transpose(x$edge), ~cbind(match(.x$.vertex0, x$vertex$vertex_), 
                                        match(.x$.vertex1, x$vertex$vertex_))))
nb <- tibble(vertex = c(t(as.matrix(x$edge[c(".vertex0", ".vertex1")]))), neighbour = c(x$vertex$vertex_[neigh0])) %>% 
  dplyr::distinct() %>% dplyr::filter(!vertex == neighbour)


i <- 10
plot(x)
points(x$vertex[i, c("x_", "y_")], pch = 19, cex = 2, col = "dodgerblue")
nb %>% dplyr::filter(vertex == x$vertex$vertex_[i]) %>% 
  dplyr::inner_join(x$vertex, c("neighbour" = "vertex_")) %>% 
  dplyr::select(x_, y_) %>% points()

I consider that dodgr-level to map out all the pairings efficiently rather than with transpose(), but that seems to fit the bill. Note that your figure shows the result with the internal triangles that fill holes - in the convex hull points-only triangulation - those are available from anglr::DEL (but gives the wrong result if I replace TRI/DEL above, for some reason ...).

Good example!

[mdsumner]

With https://github.com/hypertidy/anglr/commit/3c29470d3cb1c9e01f028b90c0014e765c6d87e2 this code now works for either model. With a DEL triangulation we get an analogous result as above https://github.com/hypertidy/silicate/issues/62#issuecomment-372853096

(Note it's not the same as tripack doesn't have shape-preserving edges, and the lower quad part of the hole triangulates differently).

library(silicate)
library(dplyr)
##d <- anglr::DEL(minimal_mesh)
d <- silicate::TRI(minimal_mesh)
x <- SC(d)
neigh0 <- do.call(rbind, purrr::map(purrr::transpose(x$edge), ~cbind(match(.x$.vertex0, x$vertex$vertex_), 
                                                                     match(.x$.vertex1, x$vertex$vertex_))))
nb <- tibble(vertex = c(t(as.matrix(x$edge[c(".vertex0", ".vertex1")]))), 
             neighbour = c(x$vertex$vertex_[t(neigh0[,2:1])])) %>% 
  dplyr::distinct() %>% dplyr::filter(!vertex == neighbour)


i <- 10
plot(x)
points(x$vertex[i, c("x_", "y_")], pch = 19, cex = 2, col = "dodgerblue")
nb %>% dplyr::filter(vertex == x$vertex$vertex_[i]) %>% 
  dplyr::inner_join(x$vertex, c("neighbour" = "vertex_")) %>% 
  dplyr::select(x_, y_) %>% points()

PATH based triangulation

SC based triangulation

As an added bonus is we can nuance the mesh with RTriangle (and showing me that my edges aren't coming out right because done per-path):

## max_area maps to RTriangle::triangulate(, a = max_area)
d <- anglr::DEL(minimal_mesh, max_area = 0.01)

[mpadge]

Now that is impressive! I'll at least replace tripack with decido for the mo, and implement the other stuff down the line. Thanks for the superb examples!

[mpadge]

Further thoughts on this that actually reflect a separate issue, but am writing here to keep within this excellent flow of discussion: Why not offer a point-based (Delaunay) triangulation method for silicate::TRI, beyond the edge-based methods currently implemented? There are likely to be many use cases where this would be useful, and not having it may restrict the generality of silicate. My spatialcluster package is an example - if silicate::TRI could be point-based, then I'd use that straight off.

The generality I refer to would apply to all simple point-based analyses. My suspicion is that we've both been so focussed on complex and code-breaking use-cases of MULTI-BLAH and above that we've neglected the most simple case of points alone. The way I see it, it can't make sense to implement edge-based triangulation in those cases, because the edges convey no information - they just need to be constructed as ignorantly as possible, and that's a Delaunay-type scheme. Right?

[mdsumner]

I hadn't even thought of this and it's an excellent idea, and I've conveniently made a comparison of available options:

http://rpubs.com/cyclemumner/369623

My pick would be geometry, which I think is ok as qhull is GPL-3 and the only other deps are abind and magic. No other choice can work given dependencies or speed. That has the added benefit of delaunayn which is useful for wrapping XYZ points in a closed mesh-shell ;)

Would you require a method to be point-based, or provide an escape hatch argument to ignore input edges? It might mean another model. It might be a wrapper package or function on decido, with a function that takes xy, and calls earcut if indicated that it's a path. It does bring up the bogey of the degenerate PATH model for points, which probably should be fixed.

I'm not at all bound to any of this (the core defs are likely to change a bit anyway, so anything is possible).

[mdsumner]

Check out https://github.com/davidcsterratt/RTriangle/pull/13 - would have to be implemented in anglr with the RTriangle dep