Pre- and Post-Processing Tools

In addition to using CARPentry directly for pre- and post-processing a number of specialized standalone tools is also available for specific purposes. The most widely used ones are summarized below.

Mesh generation (mesher)

Section author: Gernot Plank and Aurel Neic

Introduction

mesher is a program for generating simple regular FE meshes. It can also produce element tag regions and fiber definitions during element generation.

Usage

The most important parameters of mesher are

  • size : Set the size of the mesh in each axis direction.
  • resolution : Set the resolution of the mesh in each axis direction.
  • bath : Set the bath size in each axis direction. Negative sizes denote a bath on both sides.
  • fibers : Set fiber rotation.
  • mesh : Set output mesh name.

The unit of the size and resolution parameters is centimeters, while the unit of the resolution is set in micrometers.

Fiber orientation is assigned based on following assumptions:

  • The negative z axis is the transmural direction.
  • The positive y axis is the apico-basal direction.
  • The positive x axis is the circumferential direction.
../_images/Mesher_wedge_layout.png

Convention for the wedge orientation.

Example

The following example creates a myocardial slab of dimensions 0.05 x 0.6 x 0.3 cm immersed in a bath of 0.02 cm on each side:

$ bin/mesher \\
-size[0] 0.05 \\
-size[1] 0.6 \\
-size[2] 0.3 \\
-bath[0] -0.02 \\
-bath[1] -0.02 \\
-bath[2] -0.02 \\
-resolution[0] 100.0 \\
-resolution[1] 100.0 \\
-resolution[2] 100.0 \\
-mesh block

By default, the geometry’s center is located at (0,0,0). As specified, the bounding box of the intracellular domain is

Bbox:
x: (      -250.00 ,       250.00 ), delta 500.00
y: (     -3000.00 ,      3000.00 ), delta 6000.00
z: (     -1500.00 ,      1500.00 ), delta 3000.00

while the bounding box of the overall mesh is

Bbox:
x: (      -450.00 ,       450.00 ), delta 900.00
y: (     -3200.00 ,      3200.00 ), delta 6400.00
z: (     -1700.00 ,      1700.00 ), delta 3400.00

Rule-based fibers can be added with the parameters in the fibers parameter structure. In the previous example, fiber rotation can be added with:

-fibers.rotEndo 60 \\
-fibers.rotEpi -60 \\
-fibers.sheetEndo 90 \\
-fibers.sheetEpi  90 \\

The resulting geometry is:

../_images/Mesher_block.png

A slab mesh produced by mesher.

Tag regions

By default, element tag region 0 is assigned to the bath and region 1 is assigned to the myocardium. Additional regions can be assigned via regdef parameters. A regdef parameter defines a geometric shape. Elements inside this shape are assigned a given element tag.

The supported region definition shapes are:

$ mesher +Help regdef[0].type

regdef[0].type:
    how to define region

    type:   Int
    default:(Int)(0)
    menu: {
            (Int)(2)        cylinder
            (Int)(1)        sphere
            (Int)(0)        block
    }

The following example shows a minimal definition of a block region:

-numRegions      1 \\
-regdeg[0].type  0 \\
-regdef[0].tag   TAG \\
-regdef[0].p0[0] P0X \\
-regdef[0].p0[1] P0Y \\
-regdef[0].p0[2] P0Z \\
-regdef[0].p1[0] P1X \\
-regdef[0].p1[1] P1Y \\
-regdef[0].p1[2] P1Z

Mesh Management (meshtool)

Submesh extraction

Section author: Aurel Neic

With the “extract mesh” mode, a submesh can be extracted from a mesh based on tag regions. The mode’s parameters are as following:

$ meshtool extract mesh

Mesh extract error: Insufficient parameters provided.
extract mesh: a submesh is extracted from a given mesh based on given element tags
parameters:
-msh=<path>         ... (input) path to basename of the mesh to extract from
-tags=tag1,tag2     ... (input) ","-seperated list of tags
-ifmt=<format>      ... (optional) mesh input format. may be: carp_txt, carp_bin, vtk, vtk_bin, mmg, stellar
-ofmt=<format>      ... (optional) mesh output format. may be: carp_txt, carp_bin, vtk, vtk_bin, vtk_polydata, mmg, stellar
-submsh=<path>      ... (output) path to basename of the submesh to extract to

Therefore a typical call would be:

meshtool extract mesh -msh=/home/aneic/meshes/TBunnyC/mesh/TBunnyC -submsh=TBC_125 -tags=125

Both input and output mesh formats are the default CARP-text formats. In addition to the submesh, meshtool also writes two *.nod and *.eidx files. These hold the node and element mappings between the submesh and the original mesh and are needed by modes that relate submeshes with their associated meshes (e.g. the “insert” and “map” modes).

Surface extraction

Meshtool offers a potent interface for mesh surface extraction. The mode’s help description is:

$ meshtool extract surface

Surface extract error: Insufficient parameters provided.
extract surface: extract a sequence of surfaces defined by set operations on element tags
parameters:
-msh=<path>           ... (input) path to basename of the mesh
-surf=<path>          ... (output) list of names associated to the given operations
-op=operations        ... (optional) list of operations to perform. By default, the surface of the full mesh is computed.
-ifmt=<format>        ... (optional) mesh input format. may be: carp_txt, carp_bin, vtk, vtk_bin, mmg, purk, stellar
-ofmt=<format>        ... (optional) mesh output format. may be: carp_txt, carp_bin, vtk, vtk_bin, vtk_polydata, mmg, stellar
-edge=<deg. angle>    ... (optional) surface elements connected to sharp edges will be removed. A sharp edge is defined
                          by the nodes which connect elements with normal vectors at angles above the given threshold.
-coord=<xyz>:<xyz>:.. ... (optional) restrict surfaces to those elements reachable by surface edge-traversal
                          from the surface vertices closest to the given coordinates.
                          If -edge= is also provided, sharp edges will block traversal, thus limit what is reachable.
-size=<float>         ... (optional) surface edge-traversal is limited to the given radius from the initial index.

The format of the operations is:
tagA1,tagA2,..[+-:]tagB1,tagB2../tagA1,..[+-:]tagB1../..
Tag regions separated by "," will be unified into submeshes. If two submeshes are separated by "-", the rhs mesh surface
will be removed from the lhs surface (set difference). Similarly, using "+" will compute the surface union.
If the submeshes are separated by ":", the set intersection of the two submesh surfaces will be computed.
Individual operations are separated by "/".
The number of names provided with "-surf=" must match the number of operations. If no operations are provided,
the surface of the whole geometry will be extracted. Individual names are separated by ",".
Further restrictions can be added to the surface extraction with the -edge= , -coord= , -size= options.

The surface of a whole mesh can therefore be extracted with:

$ meshtool extract surface -msh=TBC_i -surf=TBC_i
Reading mesh: TBC_i.*
Reading elements in text CARP format:                       [====================]
Reading points in text CARP format:                         [====================]
Reading fibers (2 directions) in text CARP format:          [====================]
Done in 7.47963 sec
Computing surface 1/1
Done in 0.870943 sec
Writing surface TBC_i.surf
Writing vertices TBC_i.surf.vtx
Writing neubc TBC_i.neubc
Done in 0.103977 sec

Here, a surface TBC_i.surf is generated for the mesh TBC_i.elem, TBC_i.pts. Further, *.surf.vtx file (holding the unique surface vertices) and a *.neubc file (containing data needed for mechanics Neumann boundary conditions) are written. If “-ofmt” is specified, the surface is also outputted as a standalone mesh in the desired format.

More complex surfaces can be extracted with set operations on tag sets. For this, the “-op” parameter needs to be specified. In the following example, the set intersection of the surface of tag region 125 and the surface of tag region 150 is computed:

meshtool extract surface -msh=TBC_i -surf=TBC_i_125:150 -op=125:150

The result looks like:

../_images/Meshtool_TBC_i.125:150.surf.png

The surface intersection between tags 125 and 150.

Finally, surfaces can be computed by traversing the surface nodes connected to a given node list. Sharp edges can be set to block surface traversal. The “-edge” parameter defines the sharp edge angle threshold between surface elements. For example the base of the TBunnyC myocardium can be computed as following:

$ meshtool extract surface -msh=TBC_i -surf=TBC_i.top -idx=506478 -edge=10
Reading mesh: TBC_i.*
Reading elements in text CARP format:                       [====================]
Reading points in text CARP format:                         [====================]
Reading fibers (2 directions) in text CARP format:          [====================]
Done in 7.49793 sec
Setting up n2e / n2n graphs ..
Computing connectivity graphs:                              [====================]
Computing surface 1/1
Done in 0.859898 sec
Computing connectivity graphs:                              [====================]
Nodes reached:                                              [====================]
Writing surface TBC_i.top.surf
Writing vertices TBC_i.top.surf.vtx
Writing neubc TBC_i.top.neubc
Done in 0.460613 sec

The result looks like:

../_images/Meshtool_TBC_i.top.png

The surface of the TBunnyC myocardium base.

Meshtool query function

Query list of tags

meshtool query tags -msh=MESHNAMEWITHOUTEXTENSION -ifmt=%

where % is either carp_txt, carp_bin, vtk, vtk_bin depending on the mesh to be checked. In case you are checking a vtk file make sure that the name of the tag variable in the vtk file is “elemTag”. This command prints a list of all tags comprised in the mesh. An example output is shown below:

Reading mesh: /home/plank/meshes/TBunnyC2/mesh/TBunnyC2.*
Reading elements in text CARP format:                       [====================]
Reading points in text CARP format:                         [====================]
Reading fibers (2 directions) in text CARP format:          [====================]
Done in 0.638392 sec
All tags:
0 1 2 100 101 102 103 104 105 106 150 151 152 153 154 155 156 200 201 202 203 204 205 206
Extracting myocardium ..
Done in 0.220971 sec
Myocardium tags:
100 101 102 103 104 105 106 150 151 152 153 154 155 156 200 201 202 203 204 205 206

Retrieve resolution statistics of a mesh

Mean and variability of edge lengths in a FEM mesh is an important metric as spatial resolution influences to some degree numerical computation. An edge lengths statistics is retrieved by

meshtool query edges -msh=MESHNAMEWITHOUTEXTENSION -ifmt=%

This plots a text-based edge length histograms to the console:

Reading mesh: /home/plank/meshes/TBunnyC/mesh/TBunnyC.*
Reading elements in binary CARP format:                     [====================]
Reading points in binary CARP format:                       [====================]
Done in 2.92863 sec
Setting up n2e / n2n graphs ..
Computing connectivity graphs:                              [====================]
Done in 3.47954 sec

Number of connections

Average: 13.8348, (min: 4, max: 31)

 == Histogramm ==

         --------------------------------------------------------------------------------------------
   +0.26 |                                   *                                                      |
         |                               *****                                                      |
   +0.22 |                             **     *                                                     |
         |                           **       *                                                     |
         |                       ****         **                                                    |
         |                     ***             *                                                    |
   +0.15 |                     *               *                                                    |
         |                    **                *                                                   |
         |                    *                 *                                                   |
         |                   **                 **                                                  |
   +0.07 |                   *                   *                                                  |
         |                  **                   *                                                  |
         |                  *                    ***                                                |
         |               ****                      *****                                            |
   +0.00 |****************                             ****************************************     |
         --------------------------------------------------------------------------------------------
          4.0             8.9            13.7            18.6            23.4            28.3

Edge lengths

Average: 326.493, (min: 37.0675, max: 1840.16)

 == Histogramm ==

         --------------------------------------------------------------------------------------------
   +0.53 |                                                                                          |
         |        **                                                                                |
   +0.46 |        ***                                                                               |
         |       ** *                                                                               |
         |       *  **                                                                              |
         |       *   **                                                                             |
   +0.30 |      **    *                                                                             |
         |      *     **                                                                            |
         |      *      **                                                                           |
         |     **       *                                                                           |
   +0.15 |     *         *                                                                          |
         |     *         **                                                                         |
         |    *           **                                                                        |
         |  ***            ****                                                                     |
   +0.00 |***                 *****************************************************************     |
         --------------------------------------------------------------------------------------------
         37.1           361.2           685.4          1009.5          1333.7          1657.8

Retrieve mesh quality statistics

A statistics on mesh quality can also be retrieved using the query command. So far, quality calculation is only implemented for tetrahedral elements. In tissue scale simulations it is a frequent source of errors that meshes are generated without properly checking element quality.

meshtool query quality -msh=MESHNAMEWITHOUTEXTENSION -ifmt=%

Shown below is a mesh quality histogram as output to the console. The quality metric outputs the quality of an element on a scale ranging from 0 to 1 where 0 indicates best quality and 1 is poorest quality.

Reading mesh: /home/plank/meshes/TBunnyC/mesh/TBunnyC.*
Reading elements in binary CARP format:                     [====================]
Reading points in binary CARP format:                       [====================]
Done in 0.491191 sec
Computing mesh quality ..
Done in 0.202803 sec

-------- mesh statistics ---------
Number of elements: 5082272
Number of nodes: 862515

-------- element quality statistics ---------
Used method: tet_qmetric_weightedShapeDist with 0 == best, 1 == worst
Min: 3.98307e-05 Max: 0.873355 Avrg: 0.125359

Histogramm of element quality:

         --------------------------------------------------------------------------------------------
   +0.12 |                                                                                          |
   +0.12 |      ****                                                                                |
         |     **  **                                                                               |
         |    **    **                                                                              |
         |    *      *                                                                              |
   +0.09 |   **      **                                                                             |
         |   *        *                                                                             |
         |  **        **                                                                            |
         |  *          *                                                                            |
   +0.07 | **          **                                                                           |
         | *            *                                                                           |
         | *            **                                                                          |
         | *             *                                                                          |
   +0.05 |*              **                                                                         |
         |*               **                                                                        |
         |*                *                                                                        |
         |*                **                                                                       |
   +0.02 |                   **                                                                     |
         |                    **                                                                    |
         |                     ***                                                                  |
         |                       *****                                                              |
   +0.00 |                           ************************************************************** |
         --------------------------------------------------------------------------------------------
          0.0             0.2             0.4             0.5             0.7             0.9

Mesh quality improvement

In the “clean quality” mode, the quality of a mesh can be improved by vertex displacement and by volumetric smoothing. The effectiveness of either methods depends on the reason for the bad element quality (e.g. mesh has been clipped, mesh has been smoothed).

Vertex displacement restores quality for elements that have been “flattened-out”, either explicitly by mesh smoothing, or implicitly by the mesh generation algorithm. It’s drawback is that it destroys surface smoothness.

Vertex smoothing redistributes vertices in equal distance of each other. This improves mesh quality if the mesh has been modified by clipping.

The two methods can be used alone or combined. Depending on the underlying mesh problem, they can complement or act against each other. For example, if the mesh generation algorithm was configured to produce very smooth surfaces, it has probably also produced flat elements at the interfaces of the surfaces. This flat elements can only by improved by vertex displacement which will reduce surface smoothness. Vertex smoothing is not useful for this case, since it would only counter the vertex displacement by restoring surface smoothness.

In the following example, the mesh quality is restored for a LV mesh where an artery has been clipped. Here, mesh smoothing will be generate the main quality improvement. Therefore, the smoothing coefficient (smth=0.3 means that a vertex will be moved 30 percent towards its center-point) and number of iterations (iter=150) are set rather high. Surface files are passed via the “surf” parameter, in order to protect them and their edges (edge detection is turned on via edge=30).

$ meshtool clean quality -msh=OPBG010.cut -thr=0.9 -surf=OPBG010.cut,OPBG010.cut.blood -smth=0.3 -iter=150 -edge=30 -outmsh=OPBG010.cut.cln -ofmt=vtk_bin
Using OpenMP parallelization with 6 threads.
Reading mesh: OPBG010.cut.*
Reading elements in text CARP format:                       [====================]
Reading points in text CARP format:                         [====================]
Reading fibers (1 direction) in text CARP format:           [====================]
Mesh consists of 2191988 elements, 383188 nodes.
Correcting inside-out tets ..
Pass 1: Corrected 0 inside-out elements.
Done in 3.48826 sec
Reading surface: OPBG010.cut.surf
Reading elements in text CARP format:                       [====================]
Reading surface: OPBG010.cut.blood.surf
Reading elements in text CARP format:                       [====================]
Setting up n2e / n2n graphs for mesh surface ..
Computing connectivity graphs:                              [====================]
Setting up line interfaces for mesh surface ..
Computing connectivity graphs:                              [====================]
Computing connectivity graphs:                              [====================]
Done in 5.89264 sec

Checking mesh quality ..
Mesh quality statistics:
min: 6.69769e-05, max: 2.78184, avg: nan
Number of elements above threshold: 1381
Done in 0.15244 sec

Improving mesh ..
Setting up n2e / n2n graphs for mesh ..
Computing connectivity graphs:                              [====================]
Starting with optimization ..
Iter 0: Number of elements left above threshold: 1381
Smoothing progress:                                         [====================]
Smoothing progress:                                         [====================]
Smoothing progress:                                         [====================]
Smoothing progress:                                         [====================]
Iter 1: Number of elements left above threshold: 51
Smoothing progress:                                         [====================]
Smoothing progress:                                         [====================]
Smoothing progress:                                         [====================]
Smoothing progress:                                         [====================]
Iter 2: Number of elements left above threshold: 64
Mesh quality statistics:
min: 6.69769e-05, max: 0.997674, avg: 0.105546
Done in 12.1368 sec
Correcting inside-out tets ..
Pass 1: Corrected 20 inside-out elements.
Pass 2: Corrected 0 inside-out elements.
Writing mesh: OPBG010.cut.cln.*
Writing vtk file in binary format:                          [====================]
Done in 1.15729 sec

The original mesh had element of such a poor quality that the element quality couldn’t even be computed, while the resulting mesh is of useable quality.

The result looks like:

../_images/Meshtool_cutcomp.png

Comparison between original mesh and imporved mesh for a clipped geometry.

igbutils

For all tools, help is output with the -h option.

igbhead

Operations to perform on the headers of IGB files. Most commonly, it it used to print the header information but can also be usd to modify it. This utility can also be used to change the storage type of the file, or put an IGB header on to ASCII and binary data file.

igbapd

A simple utility to compute action potential durations.

igbops

This utility can perform basic math operations on one or two IGB files. There are preset functions as well as a parser for arbitrary functions. For example, one can difference two IGB files, find maxima over time, or apply a box filter.

igbextract

This is a tool to extract a hyperslab of data from an IGB file. The format of the output can be chosen between ASCII, binary or IGB.

igbdft

This is a tool for simple frequency domain techniques. It contains a limited number of operations which include Butterworth filtering, DFTs, and computing phase.

FEM tools

Tools recognize the -h to display detailed help.

GlElemCenters

Compute the element centers of a nodally defined mesh

GlFilament

Compute phase singularities and filaments.

GlGradient

Compute gradients of grid data and conduction velocities if activation times are used.

GlInterpolate

Interpolate data restitution-protocol-definition, from nodes to element centers or vice versa.

GlRuleFibers

Compute the fibres for a biventricular mesh.

GlTransfer

Transfer data between spatially overlapping grids.