[ITK] [ITK-users] watershed on surface meshs

Matt McCormick matt.mccormick at kitware.com
Sat Dec 2 10:11:30 EST 2017


Hi Roman,

Another possible idea is to use ITK's mesh fast marching capabilities:

  https://itk.org/ITKExamples/src/Filtering/FastMarching/ComputeGeodesicDistanceOnMesh/Documentation.html

with curvature as a speed function in same optimization approach used
to compute minimal paths on images:

  http://www.insight-journal.org/browse/publication/213

HTH,
Matt

On Fri, Nov 24, 2017 at 10:16 PM, Richard Beare <richard.beare at gmail.com> wrote:
> 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.
>
> 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.
>
> library(liqueueR)
>
> ## Mostly a copy of PriorityQueue from liqueueR
> StablePriorityQueue <- setRefClass("StablePriorityQueue",
>                                    contains = "Queue",
>                                    fields = list(
>                                      count = "numeric",
>                                      entries = "numeric",
>                                      priorities = "numeric",
>                                      prioritise = "function"
>                                    ),
>                                    methods = list(
>                                      sort_ = function() {
>                                        order = order(priorities, entries,
> decreasing = TRUE, partial = size():1)
>                                        #
>                                        data <<- data[order]
>                                        priorities <<- priorities[order]
>                                        entries <<- entries[order]
>                                      },
>                                      push = function(item, priority = NULL)
> {
>                                        'Inserts element into the queue,
> reordering according to priority.'
>                                        callSuper(item)
>                                        #
>                                        if (is.null(priority)) priority =
> prioritise(item)
>                                        #
>                                        priorities <<- c(priorities,
> priority)
>                                        entries <<- c(entries, count)
>                                        count <<- count - 1
>                                        #
>                                        sort_()
>                                      },
>                                      pop = function(N = 1) {
>                                        # 'Removes and returns head of queue
> (or raises error if queue is empty).'
>                                        if (N > size()) stop("insufficient
> items in queue!")
>                                        priorities <<-
> priorities[-seq_len(N)]
>                                        entries <<- entries[-seq_len(N)]
>                                        callSuper(N)
>                                      },
>                                      initialize = function(prioritise =
> NULL, ...) {
>                                        'Creates a PriorityQueue object.'
>                                        callSuper(...)
>                                        #
>                                        ## to enforce FIFO order
>                                        count <<- 0
>                                        if (is.null(prioritise))
>                                          .self$prioritise = function(x) 0
>                                        else
>                                          .self$prioritise = prioritise
>                                        #
>                                        .self
>                                      }
>                                    )
> )
>
>
>
> igraph.watershed <- function(Gr, labelfield, unlabelled, vertexid, alltimes)
> {
>   Gres <- Gr
>   ## Watershed, without marking boundary (Beucher)
>   ## 1. find all marker nodes that have a background neighbour
>   lablist <- which(vertex_attr(Gr, labelfield) != unlabelled)
>   nlist <- ego(Gr, 1, lablist)
>   uu <- which(map_lgl(nlist, ~any(vertex_attr(Gr, labelfield,
> .x)==unlabelled)))
>   ## uu indexes the labelled vertexes
>   ## need indexes into all vertexes
>   uu <- lablist[uu]
>   ## create priority queue
>   qq <- StablePriorityQueue$new()
>   ## Insert the boundary markers
>   kk <- lapply(uu, qq$push, priority=0)
>   all.labels <- get.vertex.attribute(Gr, labelfield)
>   all.ids <- get.vertex.attribute(Gr, vertexid)
>   alltimes <- subset(alltimes, from %in% all.ids, select=c("from",
> "Hospital", "minutes"))
>   dd <- duplicated(alltimes[, c("from", "Hospital")])
>   alltimes <- alltimes[!dd,]
>   alltimes.wide <- spread_(alltimes, key="Hospital", value="minutes")
>   ## Make the order the same as all.ids - so now we'll be able to index by
> number
>   rownames(alltimes.wide) <- alltimes.wide$from
>   alltimes.wide <- alltimes.wide[all.ids, ]
>   cc <- 1:ncol(alltimes.wide)
>   names(cc) <- colnames(alltimes.wide)
>   while (qq$size() > 0) {
>     vid <- qq$pop()
>     ## Get the neighbours
>     nb <- neighborhood(Gr, 1, vid, mode="all", mindist=1)[[1]]
>     nlabs <- all.labels[nb]
>     ## Are any neighbours unknown?
>     ul <- nlabs==unlabelled
>     if (any(ul)) {
>       this.label <- all.labels[vid]
>       nbb <- nb[ul]
>       all.labels[nbb] <- this.label
>       this.label.idx <- cc[this.label]
>       priorities <- alltimes.wide[[this.label.idx]][nbb]
>       kk <- map2(nbb, priorities, ~qq$push(.x, priority = .y* -1))
>     }
>   }
>   return(data.frame(PlaceID=all.ids, Hospital=all.labels, stringsAsFactors =
> FALSE))
> }
>
>
>
>
> On Fri, Nov 24, 2017 at 10:26 PM, Grothausmann, Roman Dr.
> <grothausmann.roman at mh-hannover.de> wrote:
>>
>> Dear mailing list members,
>>
>>
>> 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:
>>
>> Partitioning 3D surface meshes using watershed:
>>
>> http://teacher.en.rmutt.ac.th/ktw/Resources/Full%20paper%20PDF/Partitioning%203D%20surface%20meshes%20using%20watershed%20segmentation.pdf
>>
>> Surface Segmentation Using Morphological Watersheds:
>>
>> 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
>>
>> 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?
>>
>> Thanks for any help or hints.
>> Roman
>>
>>
>> --
>> Dr. Roman Grothausmann
>>
>> Tomographie und Digitale Bildverarbeitung
>> Tomography and Digital Image Analysis
>>
>> Medizinische Hochschule Hannover
>> Institut für Funktionelle und Angewandte Anatomie
>> OE 4120, Carl-Neuberg-Str. 1, 30625 Hannover, Deutschland
>>
>> Tel. +49 511 532-2900
>> grothausmann.roman at mh-hannover.de
>> http://www.mh-hannover.de/anatomie.html
>>
>> The ITK community is transitioning from this mailing list to
>> discourse.itk.org. Please join us there!
>> ________________________________
>> Powered by www.kitware.com
>>
>> Visit other Kitware open-source projects at
>> http://www.kitware.com/opensource/opensource.html
>>
>> Kitware offers ITK Training Courses, for more information visit:
>> http://www.kitware.com/products/protraining.php
>>
>> Please keep messages on-topic and check the ITK FAQ at:
>> http://www.itk.org/Wiki/ITK_FAQ
>>
>> Follow this link to subscribe/unsubscribe:
>> http://public.kitware.com/mailman/listinfo/insight-users
>>
>
>
> The ITK community is transitioning from this mailing list to
> discourse.itk.org. Please join us there!
> ________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.html
>
> Kitware offers ITK Training Courses, for more information visit:
> http://www.kitware.com/products/protraining.php
>
> Please keep messages on-topic and check the ITK FAQ at:
> http://www.itk.org/Wiki/ITK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://public.kitware.com/mailman/listinfo/insight-users
>
The ITK community is transitioning from this mailing list to discourse.itk.org. Please join us there!
________________________________
Powered by www.kitware.com

Visit other Kitware open-source projects at
http://www.kitware.com/opensource/opensource.html

Kitware offers ITK Training Courses, for more information visit:
http://www.kitware.com/products/protraining.php

Please keep messages on-topic and check the ITK FAQ at:
http://www.itk.org/Wiki/ITK_FAQ

Follow this link to subscribe/unsubscribe:
http://public.kitware.com/mailman/listinfo/insight-users


More information about the Community mailing list