Biomech GenIso Mesher
This mesher is used to generate a 3D surface mesh from a 3D segmented image.
The input is a 3D segmented binary image. It contains two values: voxels with a $0$ value are located outside of the object to mesh, voxels with a non-zero value are located inside.
The image must be
.vtk, or event better
.vti, which is more compact. If it is not the right format use Slicer3D, Paraview or my interface
geniso/tools/gui4.py to change it.
First, the geniso tools must be loaded:
import geniso.tools.generalTools as generalTools import geniso.tools.meshingTools as meshingTools import geniso.tools.imagingTools as imagingTools
Then the image must be loaded:
image = generalTools.loadVtkImageXML(imageFileName)
Afterwards, the following set of commands is used to generate a surface mesh of the boundary of the segmented image:
gen = GenerateMesh(generalTools.vtkImageDataToCpp(image)) gen.execute() smooth = SmoothMesh(gen.getSurf()) smooth.smoothWithFctFromPoly(nbr,spx,spy,spz,nmin,interpol,type)
where the parameters are described in the table below.
GenerateMesh(image).execute() creates a triangulation of the boundary of the segmented object in the binary image given. Because of the spatial discretization of the image, the mesh which is generated is not smooth, since the steps of the image are kept. The class
SmoothMesh() is used to eliminate these irregularities while keeping the geometry described in the image. It does so with the construction of an implicit function, described by the parameters
type (for more details: http://hdl.handle.net/2268/136159).
vtkPolydata is retrieved with the command:
poly = generalTools.genisoMeshToPolyData(smooth.getSurf())
TetGen is used to get a tetrahedral volume mesh of the object from the triangular surface mesh of the boundary of the object:
ugrid = meshingTools.callTetgen(poly)
vtkUnstructuredGrid is retrieved. To save it on disk:
generalTools.saveUgridXML('ugrid.vtu',ugrid). Paraview is used to view it. The class
LoadMesh, implemented in the file, is used to include it in a Metafor test case:
groupset.add(Group(1)) mesh = LoadMesh(domain, ugrid, groupset(1) ) mesh.execute()
This option is used to generate a mesh which contains several material areas: two different bones, a brain with ventricles separated from a tumor, a jawbone with a series of teeth…
In the multi-material case, the segment image representing the geometry must contain several labels: 0 for the background (outside of the object), and a series of non-zero positive values (typically 1, 2, 3…) for each material.
To use the GenIso mesher, the implicit function
Fct must first be defined. This function is a continuous equivalency of the discrete segment image, and will allow the mesher to supply a smooth interface between the mesh separating each region.
To define this function, couples of labels of the image are defined, labels between which the mesher must create a triangulation. In addition, these couples are sorted to define the heterogeneous object as a set of closed surfaces on which open surfaces can be attached. In the image below, a tumor (open surface defined by the couple (4, 3)), is tied on the external boundary of the brain (closed surface defined by the union of the couples (3, 0) and (4, 0)).
Consequent, the syntax is as follow:
fct = Fct() surf1 = fct.defineClosedSurface() surf1.mesh(3,0) surf1.mesh(4,0) surf2 = fct.defineOpenSurface() surf2.mesh(4,3) fct.setEpsilon(eps) fct.setNmin(nmin) fct.setType(type) fct.setInterpol(interpol) fct.setVtkImageFromPython(generalTools.vtkImageDataToCpp(image)) fct.build()
Afterwards, the mesh can be generated with the commands:
gen = GenerateMesh(fct) gen.getRegularGrid().setResolution(res) gen.execute() poly = generalTools.genisoMeshToPolyData(gen.getSurf()) // vtkPolyData
From this surface multi-regions mesh, TetGen is used to generate the volume mesh in each of these regions. To do so, TetGen must be given a point, seed situated in each of the regions to be meshed (this point can be found, for example, by loading the
vtkPolydata in the graphical interface
geniso/tools/gui4.py using the option point).
seedbrain= [118.3,160.5,75] seedtumor = [118,196.9,71] regionSeeds = [seedbrain,seedtumor] regionVolmax = [volmax, volmax] ugrid = meshingTools.callTetgenMultipleRegions(poly, regionSeeds, regionVolmax)
Then, the mesh is loaded in Metafor:
labels = [1,2] // the two labels of the regions located in the ugrid (see visu in Paraview) grp = [groupset(1),groupset(2)] mesh = tetgen.LoadMesh(domain, ugrid, grp) mesh.setFillMeshPoints(True) mesh.setGrpLabels(labels) mesh.execute()
|Parameter||Meaning and recommended value|
| || number of iterations of the smoothing function
| || Spacing of the grid used to generate the mesh
recommended = gen.getRegularGrid().getDx(),.getDy(),.getDz()
| || Resolution of the grid used to generate the mesh
recommended: image.GetDimensions() if the grid defined by the initial image is to be kept
| || relative distance allowed between the initial image and the function Fct
recommended: 1./image.GetDimensions(), which corresponds to a distance of a voxel width
| || smoothing degree
recommended: between 15 and 100
| || 0: approximation function, 1: interpolation function
recommended: 0 for a smoother surface, 1 for a geometry closer to the initial image
| || 1: linear, 2: quadratic
recommended: 1 for robustness, 2 for a smoother surface