<div dir="ltr">Just remembered something that may be useful for testing. I wrote the following for doing watersheds on geospatial graphs using R/igraph. Certainly isn't efficient, but may be useful for testing if you can turn your mesh+curvature into something like a text file.<div><br></div><div>This was set up for allocating points on maps to destinations based on travel time. It performed a kind of smoothing to smooth non-sensical travel time estimates, but the basic algorithm is there. It will take some fiddling. You could probably do something equivalent with python, probably doing better than this priorty queue, but using much the same approach with igraph.</div><div><br></div><div><div>library(liqueueR)</div></div><div><br></div><div><div>## Mostly a copy of PriorityQueue from liqueueR</div><div>StablePriorityQueue <- setRefClass("StablePriorityQueue",</div><div> contains = "Queue",</div><div> fields = list(</div><div> count = "numeric",</div><div> entries = "numeric",</div><div> priorities = "numeric",</div><div> prioritise = "function"</div><div> ),</div><div> methods = list(</div><div> sort_ = function() {</div><div> order = order(priorities, entries, decreasing = TRUE, partial = size():1)</div><div> #</div><div> data <<- data[order]</div><div> priorities <<- priorities[order]</div><div> entries <<- entries[order]</div><div> },</div><div> push = function(item, priority = NULL) {</div><div> 'Inserts element into the queue, reordering according to priority.'</div><div> callSuper(item)</div><div> #</div><div> if (is.null(priority)) priority = prioritise(item)</div><div> #</div><div> priorities <<- c(priorities, priority)</div><div> entries <<- c(entries, count)</div><div> count <<- count - 1</div><div> #</div><div> sort_()</div><div> },</div><div> pop = function(N = 1) {</div><div> # 'Removes and returns head of queue (or raises error if queue is empty).'</div><div> if (N > size()) stop("insufficient items in queue!")</div><div> priorities <<- priorities[-seq_len(N)]</div><div> entries <<- entries[-seq_len(N)]</div><div> callSuper(N)</div><div> },</div><div> initialize = function(prioritise = NULL, ...) {</div><div> 'Creates a PriorityQueue object.'</div><div> callSuper(...)</div><div> #</div><div> ## to enforce FIFO order</div><div> count <<- 0</div><div> if (is.null(prioritise))</div><div> .self$prioritise = function(x) 0</div><div> else</div><div> .self$prioritise = prioritise</div><div> #</div><div> .self</div><div> }</div><div> )</div><div>)</div><div><br></div><div><br></div><div><br></div><div>igraph.watershed <- function(Gr, labelfield, unlabelled, vertexid, alltimes)</div><div>{</div><div> Gres <- Gr</div><div> ## Watershed, without marking boundary (Beucher)</div><div> ## 1. find all marker nodes that have a background neighbour</div><div> lablist <- which(vertex_attr(Gr, labelfield) != unlabelled)</div><div> nlist <- ego(Gr, 1, lablist)</div><div> uu <- which(map_lgl(nlist, ~any(vertex_attr(Gr, labelfield, .x)==unlabelled)))</div><div> ## uu indexes the labelled vertexes</div><div> ## need indexes into all vertexes</div><div> uu <- lablist[uu]</div><div> ## create priority queue</div><div> qq <- StablePriorityQueue$new()</div><div> ## Insert the boundary markers</div><div> kk <- lapply(uu, qq$push, priority=0)</div><div> all.labels <- get.vertex.attribute(Gr, labelfield)</div><div> all.ids <- get.vertex.attribute(Gr, vertexid)</div><div> alltimes <- subset(alltimes, from %in% all.ids, select=c("from", "Hospital", "minutes"))</div><div> dd <- duplicated(alltimes[, c("from", "Hospital")])</div><div> alltimes <- alltimes[!dd,]</div><div> alltimes.wide <- spread_(alltimes, key="Hospital", value="minutes")</div><div> ## Make the order the same as all.ids - so now we'll be able to index by number</div><div> rownames(alltimes.wide) <- alltimes.wide$from</div><div> alltimes.wide <- alltimes.wide[all.ids, ]</div><div> cc <- 1:ncol(alltimes.wide)</div><div> names(cc) <- colnames(alltimes.wide)</div><div> while (qq$size() > 0) {</div><div> vid <- qq$pop()</div><div> ## Get the neighbours</div><div> nb <- neighborhood(Gr, 1, vid, mode="all", mindist=1)[[1]]</div><div> nlabs <- all.labels[nb]</div><div> ## Are any neighbours unknown?</div><div> ul <- nlabs==unlabelled</div><div> if (any(ul)) {</div><div> this.label <- all.labels[vid]</div><div> nbb <- nb[ul]</div><div> all.labels[nbb] <- this.label</div><div> this.label.idx <- cc[this.label]</div><div> priorities <- alltimes.wide[[this.label.idx]][nbb]</div><div> kk <- map2(nbb, priorities, ~qq$push(.x, priority = .y* -1))</div><div> }</div><div> }</div><div> return(data.frame(PlaceID=all.ids, Hospital=all.labels, stringsAsFactors = FALSE))</div><div>}</div></div><div><br></div><div><br></div><div><br></div></div><div class="gmail_extra"><br><div class="gmail_quote">On Fri, Nov 24, 2017 at 10:26 PM, Grothausmann, Roman Dr. <span dir="ltr"><<a href="mailto:grothausmann.roman@mh-hannover.de" target="_blank">grothausmann.roman@mh-hannover.de</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Dear mailing list members,<br>
<br>
<br>
I need to separate a mesh at "curved corners" (see attached PNG, using the colored labels from a facet analysis do not suffice but go in the right direction). So my current thought is to run vtkCurvature to get a Gaussian curvature value per point/vertex and then try to separate regions of positive values around local maxima. Just thresholding the result of vtkCurvature does not fully separate each local max from neighboring ones, but to my understanding a surface watershed would. I found two publications by Mangan and Whitaker on this:<br>
<br>
Partitioning 3D surface meshes using watershed:<br>
<a href="http://teacher.en.rmutt.ac.th/ktw/Resources/Full%20paper%20PDF/Partitioning%203D%20surface%20meshes%20using%20watershed%20segmentation.pdf" rel="noreferrer" target="_blank">http://teacher.en.rmutt.ac.th/<wbr>ktw/Resources/Full%20paper%20P<wbr>DF/Partitioning%203D%20surface<wbr>%20meshes%20using%20watershed%<wbr>20segmentation.pdf</a><br>
<br>
Surface Segmentation Using Morphological Watersheds:<br>
<a href="https://www.google.de/url?sa=t&rct=j&q=&esrc=s&source=web&cd=4&cad=rja&uact=8&ved=0ahUKEwjD0by1lafWAhVUGsAKHZ2MAbUQFgg_MAM&url=http%3A%2F%2Fciteseerx.ist.psu.edu%2Fviewdoc%2Fdownload%3Fdoi%3D10.1.1.464.2788%26rep%3Drep1%26type%3Dpdf&usg=AFQjCNGX-p9-ElQFcpsUyBRO0pCjBKCmNg" rel="noreferrer" target="_blank">https://www.google.de/url?sa=t<wbr>&rct=j&q=&esrc=s&source=web&cd<wbr>=4&cad=rja&uact=8&ved=0ahUKEwj<wbr>D0by1lafWAhVUGsAKHZ2MAbUQFgg_<wbr>MAM&url=http%3A%2F%2Fciteseerx<wbr>.ist.psu.edu%2Fviewdoc%<wbr>2Fdownload%3Fdoi%3D10.1.1.464.<wbr>2788%26rep%3Drep1%26type%<wbr>3Dpdf&usg=AFQjCNGX-p9-ElQFcpsU<wbr>yBRO0pCjBKCmNg</a><br>
<br>
Does anybody know about an implementation for this in VTK/ITK or another open-source library? If not, would it be possible to transfer the ITK watershed implementation for voxel data to mesh data, e.g. to crate a VTKmorphWatershedFilter?<br>
<br>
Thanks for any help or hints.<br>
Roman<br>
<br>
<br>
-- <br>
Dr. Roman Grothausmann<br>
<br>
Tomographie und Digitale Bildverarbeitung<br>
Tomography and Digital Image Analysis<br>
<br>
Medizinische Hochschule Hannover<br>
Institut für Funktionelle und Angewandte Anatomie<br>
OE 4120, Carl-Neuberg-Str. 1, 30625 Hannover, Deutschland<br>
<br>
Tel. <a href="tel:%2B49%20511%20532-2900" value="+495115322900" target="_blank">+49 511 532-2900</a><br>
<a href="mailto:grothausmann.roman@mh-hannover.de" target="_blank">grothausmann.roman@mh-hannover<wbr>.de</a><br>
<a href="http://www.mh-hannover.de/anatomie.html" rel="noreferrer" target="_blank">http://www.mh-hannover.de/anat<wbr>omie.html</a><br>
<br>The ITK community is transitioning from this mailing list to <a href="http://discourse.itk.org" rel="noreferrer" target="_blank">discourse.itk.org</a>. Please join us there!<br>
______________________________<wbr>__<br>
Powered by <a href="http://www.kitware.com" rel="noreferrer" target="_blank">www.kitware.com</a><br>
<br>
Visit other Kitware open-source projects at<br>
<a href="http://www.kitware.com/opensource/opensource.html" rel="noreferrer" target="_blank">http://www.kitware.com/<wbr>opensource/opensource.html</a><br>
<br>
Kitware offers ITK Training Courses, for more information visit:<br>
<a href="http://www.kitware.com/products/protraining.php" rel="noreferrer" target="_blank">http://www.kitware.com/<wbr>products/protraining.php</a><br>
<br>
Please keep messages on-topic and check the ITK FAQ at:<br>
<a href="http://www.itk.org/Wiki/ITK_FAQ" rel="noreferrer" target="_blank">http://www.itk.org/Wiki/ITK_<wbr>FAQ</a><br>
<br>
Follow this link to subscribe/unsubscribe:<br>
<a href="http://public.kitware.com/mailman/listinfo/insight-users" rel="noreferrer" target="_blank">http://public.kitware.com/<wbr>mailman/listinfo/insight-users</a><br>
<br></blockquote></div><br></div>