Dear 3D Slicer Team,
I have generated a volumetric heart mesh and now need to extract certain geometry information, that is I want to export indices of nodes on the epicardium, indices of nodes on the left ventricular surface and indices of nodes on the right ventricular surface and save them in three separate files. Is it possible to perform such extraction in the framework of 3D Slicer?
Thank you very much for your attention!
This should be all doable, but how exactly depends on what would you like to do with it. What is your overall goal? Biomechanical analysis (FEM), electrophysiology mapping, volumetric analysis,…?
By “nodes” do you mean point indices? Have you separated the epicardium, left ventricle, and right ventricle in the segmentation (into different segments)?
Dear Andras,
Thank you for a quick reply. My overall goal is to run a model of myocardial ischemia with heart-torso mesh. Those indices are required for computational reasons (used by mesh constructor and for some algorithm implementation).
Yes, indeed. By “nodes” I mean “points”, so “node index” equals “point index”.
I have not separated epi, lv and rv in the segmentation yet, but will positively do!
What kind of model? Biomechanical, electrophysiology, geometrical, visualization? Do you need surface or volumetric mesh? If surface mesh then open or closed?
It is an electrophysiological model. I need a volumetric mesh (Tetgen format).
You can create the mesh you need by segmenting each structure as a separate segment (see for example this tutorial) and then generate volumetric mesh using Segment Mesher extension. I would recommend Cleaver over Tetgen as it works directly from the labelmap image and you don’t need to go through surface mesh representation. Cleaver can output tetgen file format directly.
What software do you use for the EP modeling?
Sure, that’s the way I usually do, and I use Cleaver over Tetgen as well.
For the EP modelling I use Chaste. It uses Tetgen format, so after mesh generation (using Segment Mesher) I convert resulting |.vtk| file into |.node|, |.ele|, |.face| files using script provided with Chaste.
By the way, I have separated epi, lv and rv into different segments in the segmentation and generated a volumetric mesh using Segment Mesher. What should be the next step to obtain the respective indices?
I know that in VTK file output of segment mesher there is a cell array that specifies which cell belong to which input segment. Probably something like that exists in the tetgen output, too.
Yes, indeed. There is an array that specifies which cell belongs to which input segment but it is not of any use since any array looks like this
…
2 2 2 2 2 2 1 1 1
1 3 3 3 3 3 3 3 1
3 3 3 3 3 3 3 1 1
…
Furthermore, I need to know which point belongs to which input segment but, unfortunately, VTK file output does not contain that information.
If I reformulate my basic question, it will sound like this. For instance, a biventricular model consists of 300k points (and of some number of cells). I need to know that, shall we say, epicardium includes points [1; 70k], left ventricular surface includes points [70001; 85k], right ventricular surface includes points [85001; 100k] and everything else includes points [100001; 200k]. Is there a way to extract such data?
In general, points are shared between cells, so segment identifier in a volumetric mesh can only be cell data.
Yes, of course, you can do this.
The straightforward implementation would be to get all the cell IDs corresponding to a segment, for example using numpy functions (e.g., np.where) then get point IDs of those cells using VTK methods.
If you have tens of thousands of cells then iterating through all of them in Python would be probably too slow. So, a smarter method would be this:
Thank you, I’ll try that.
In case I got completely stuck, would there be a simpler method with using masks, for example?
The second method that I described is probably 2-3 magnitudes faster then the naive implementation (using for loops in Python, taking fraction of a second instead of minutes) and probably less than 20 lines of code in total (and 10 lines is already available in the example that I gave above). There is nothing Slicer specific, it is just connecting a few VTK filters from Python. Give it a try and let us know if you get stuck.
Dear Andras,
I also found another solution. It is rather similar to what you’ve suggested but can be performed directly in ParaView using Data Analysis. To do so, one should import a volumetric mesh in VTK format into ParaView and then use the relative filters. It might be done through scripting interface and takes less than 20 lines of code as well.
Thank you for a very helpful conversation!
Both ParaView and Slicer uses the exact same underlying library (VTK) for manipulating meshes, so you can use either software for processing or visualizing meshes. ParaView has more mesh manipulation and visualization features than Slicer, but it is very limited in what it can do with images, transforms, and markups.