Level Sets Framework Refactoring Arnaud Gelas Kishore Mosaliganti

Level Sets Framework Refactoring Arnaud Gelas, Kishore Mosaliganti, Nicolas Rannou, Lydie Souhait, Sean Megason Boston 02/03/2011

Outline • Goals • Principles • Status o Traits o Fast Marching § Image / Mesh § Stopping Criterion § Constrained Topology § Shortest Path Computation § Isotropic / Anisotropic o Level Sets / GPU • Plan • Requirements

Goal: Generic Level Set Framework • Multi-level set support • simultaneous evolution of level sets • Multi-channel support • Vector or tensor data segmentation • Level set representation • Mesh-based (unstructured), image, or parametric • Terms used in the PDE • Add/delete terms in the update equation • Topological constraints • Stopping criterion • RMS, Iterations, target points

Discrete Level Set STATUS

Traits • TInput. Domain • TNode • TOutput. Domain • TSuperclass

Traits - Base Class template< class TInput. Domain, class TNode, class TOutput. Domain, class TSuperclass > class Level. Set. Traits { public: typedef [. . . ] class Node. Pair : public std: : pair< Node. Type, Output. Pixel. Type > [. . . ] };

Traits - Image Specialization template<unsigned int VDimension, class TInput. Pixel, typename TOutput. Pixel > class Image. Level. Set. Traits : public Level. Set. Traits< Image< TInput. Pixel, VDimension >, Index< VDimension >, Image< TOutput. Pixel, VDimension >, Image. To. Image. Filter< Image< TInput. Pixel, VDimension >, Image< TOutput. Pixel, VDimension > > >

Traits - Mesh Specialization template< unsigned int VDimension, typename TInput. Pixel, class TInput. Mesh. Traits, typename TOutput. Pixel, class TOutput. Mesh. Traits > class Mesh. Level. Set. Traits : public Level. Set. Traits< Mesh< TInput. Pixel, VDimension, TInput. Mesh. Traits >, typename TInput. Mesh. Traits: : Point. Identifier, Mesh< TOutput. Pixel, VDimension, TOutput. Mesh. Traits >, Mesh. To. Mesh. Filter< Mesh< TInput. Pixel, VDimension, TInput. Mesh. Traits >, Mesh< TOutput. Pixel, VDimension, TOutput. Mesh. Traits > >

Fast Marching • Code available at the following address: https: //github. com/arnaudgelas/itk. Fast. Marching • Numbers: o 36 tests o Tested § Fedora 13, 14 (64 bits) § Ubuntu 10. 10 (64 bits) § Mac OS-X 10. 5, 10. 6 o Coverage: 80. 49%

Stopping Criterion - Base Class class Stopping. Criterion. Base : public Object { public: virtual bool Is. Statisfied() const = 0; virtual const std: : string Get. Description() const = 0; };

Stopping Criterion - Examples • Threshold on the current value o Equivalent to the current implementation of itk: : Fast. Marching. Image. Filter<> • Reached Target Nodes (One, Some, All), with possible overshoot offset o Equivalent to the current implementation of itk: : Fast. Marching. Up. Wind. Gradient. Image. Filter<>

Constrained Topology [Tustisson'10 - Insight Journal 778] Escher's Ants as Metaphor: Topological Marching for the Well-Composed, Genus Zero Crowd

Minimal Path Extraction [Mueller'08 - Insight Journal 213] Fast Marching Minimal Path Extraction in ITK

Isotropic / Anisotropic Isotropic • Can be solved using current implementation Anisotropic • Several possible schemes • which one the best? • make it easy to implement any of these methods

Requirements Possible performance improvement Update. Neighbors() calls 2 * Image. Dimension Update. Value() Thread Pool ?

Fast Marching - Process ? • Integration Process? o Should we struggle for its integration (backward compatibility) ? o Should we struggle a second time when integrating new level sets framework?

Fast Marching - Process ? • update software guide? o When ? o How ? o Any constraint?

FUTURE WORK

Plan • Git repository • Discrete Representations o o o Domain Traits Iterators Dense § Term container § Propagation § Advection § Curvature § Chan & Vese energy § Multithread § Reinitialization § Stopping Criterion

Plan • Discrete Representations • Sparse – Constrained Topology – Multithread • Real time algorithm [Shi] • Parametric Representations • Splines • RBF

Discrete Level Sets - simplified view (a) while( ! m_Stopping. Criterion->Is. Satisfied() ) (b) { (c) for each level set ls_i in the level set container (d) { (e) for each nodes n_j in the domain of ls_i (f) { (g) for each term t_k in the term container (h) { (i) Compute Term Value t_k( n_j, ls_i ) (j) Compute Term Contribution for time step computation (k) } (l) Evaluate the updated level set function ( delta( ls_i( n_j) ) ) (m) } (n) } (o) Compute time step from CFL Condition (p) for each level set ls_i in the level set container (q) { (r) Update the level set function ls_i (s) Reinitialize to signed distance function (if requested by user) (t) } (u) }

GPU Involvement - 1 (a) while( !m_Stopping. Criterion->Is. Satisfied() ) (b) { (c) for each level set ls_i in the level set container (d) { (e) for each nodes in the domain of ls_i (f) { (g) for each term t_k in the term container (h) { (i) Compute Term Value t_k( n_j, ls_i ) GPU implementation during pixel updates at (i): • Pixel neighborhood in image and level set is copied to GPU memory • Terms are evaluated in the GPU function • Each term will have a CPU and GPU implementation • A term factory will call the GPU implementation (i) Advantages: • Minimal changes in the current proposed design • Drawbacks: • very bad according to performance

GPU Involvement -2 (a) while ( !m_Stopping. Criterion->Is. Satisfied() ) (b) { (c) for each level set ls_i in the level set container (d) { (e) for each nodes in the domain of ls_i (f) { (g) for each term t_k in the term container (h) { • Compute Term Value t_k( n_j, ls_i ) • Entire while loop iteration (a) in GPU • Everything is copied inside the GPU memory • Advantages: • Fastest solution in terms of performance • Downside: • memory limitation of the GPU (<2 Gb) • Code duplication: CPU and GPU • Note: Copy b/w memory 4 Gb/s

GPU Involvement -3 In the last scenario, the code nesting is different: (a) while( !m_Stopping. Criterion->Is. Satisfied() ) (b) { (c) for each level set ls_i in the level set container (d) { (e) for each term in the term container (f) { (g) for each nodes in the domain of ls_i (h) { (i) Evaluate the updated level set function [. . . ] In this one the GPU Implementation will occur for the most nested for loop (g) • Keep copying the level set and image in each iteration in the GPU • Second most optimal implementation for GPU • No code duplication

Questions and Comments ?