Detailed high quality surface-based mouse CAD model suitable for electromagnetic simulations

The computational mouse model prepared specifically for this work was derived from the labeled atlas data of the Digimouse archive [9, 10]. The authors of this atlas collected co-registered positron emission tomography (PET), computed tomography (CT), and cryosection data to fully characterize the internal organs and external surface of a 28 g male mouse. This data resulted in a 380,992,208 cubical voxel matrix; each voxel is cubical with a side length of 0.1 mm. Each voxel in the atlas was given a greyscale intensity value corresponding to one of twenty-one distinct tissues.

This greyscale image stack was imported into ITK-SNAP [16]. After that, each of the twenty-one tissues was segmented from the image stack. The gross features of each tissue were extracted using the semi-autonomous tools of ITK-SNAP; finer features were captured on a voxel-by-voxel basis using a variety of manual segmentation tools. Following segmentation, surface meshes were obtained and saved as stereolithographic (*.STL) files. The resulting mesh files were post-processed using standard algorithms on an as-needed basis. Mesh resurfacing was accomplished using screened Poisson surface reconstruction [17]; mesh size reduction was realized through the use of quadric edge collapse decimation [18]; element smoothing and quality enhancement was achieved using Taubin smoothing [19]. All post-processing was completed using the open-source tool Meshlab [20] and verified as CAD compatible structures (manifold, non-self-intersecting) with the commercial software Ansys SpaceClaim. All 21 surface meshes were 2-manifold and had a high triangle quality (twice the ratio of the incircle radius to the circumcircle radius). However, these meshes were very large in size. For example, the skin surface alone had 2.7 million facets. Plus, the majority of these meshes did intersect in small domains. Figure 1 shows the initial model so constructed.

Zoom In Zoom Out Reset image size

A semi-automated pipeline has been developed for converting the initial composite surface model to a set of strictly non-intersecting meshes while maintaining the closest similarity with the initial meshes with regard to their volume and other geometrical properties. This process includes the following steps:

• 1.  
  All meshes were imported into Ansys SpaceClaim 2023 R1 in STL format.
• 2.  
  A shrink-wrap operation was applied for every mesh with the goal to reduce the faceting while maintaining high mesh quality.
• 3.  
  All meshes were then imported into Ansys Electronics Desktop 2023 R1 and the mesh intersections were resolved semi-automatically using Boolean subtraction operations.

The shrink-wrap operation allows for the specification of a target facet length in which to create a new surface with uniform faceting that is said to be 'wrapped' around the input geometry. The volume of the geometry before and after facet reduction was computed and the facet length applied to each tissue was adjusted to achieve the minimum length needed to successfully retain at least 99% of the original tissue volume. Note that for the 'skin' geometry the large increase in the reduced model was due to encapsulation of folded skin. Using this technique, the total number of facets was able to be reduced by 87% while maintaining high triangle quality and anatomical accuracy. This reduction allows the model to be utilized with commercial EM solvers like Ansys HFSS.

The result of the operations from the previous section was the composite surface mouse model with approximately 670,000 facets in total and with all strictly non-intersecting tissues. However, many of the tissues still had exactly coincident facets. To physically separate the meshes with a small gap, the surface sculpt operation with any and all adjacent meshes using the free triangular surface mesh manipulator Meshmixer v3.5.474 by Autodesk has been performed. Sculpting combines affine transformations to enable mesh bending over a defined area using a controlled strength. Multiple meshes were simultaneously projected so that separation resolution could be visually confirmed prior to checking again using Boolean operations. This separation and validation procedure took approximately 20 manhours to complete.

When constructed in the manner described above, the final model possesses a set of topological characteristics necessary for

• (i)  
  cross-platform compatibility and;
• (ii)  
  computational efficiency.

Each original triangular surface tissue mesh is strictly 2-manifold or thin-shell (no non-manifold faces, no non-manifold vertices, no holes, and no self-intersections). No tissue mesh has any triangular facets in contact with other tissue surfaces. There is always a small gap between the distinct tissue surfaces. This gap physically represents connective tissues or thin membranes separating individual distinct organs from each other.

Numerically, this gap corresponds to an 'average body' container, which encloses organs and tissues, and guarantees compatibility between different CAD formats. Figure 2 shows the model topology using a cross-section view. Note that the 'skin' object is utilized by default as the 'average body' container for this model. The container type can be changed by introducing a sub container.

Zoom In Zoom Out Reset image size

There are many definitions that exist to define the quality of a triangular element as it may relate to deviation from an ideal equilateral triangle [21]. Regardless of the metric used to define it, finite element quality is essential to the success of a given simulation. Studies on the suitability of triangular mesh structure the FEM [22] and maximum internal angle of two- [23] and three-dimensional elements [24], as well as the three-dimensional Joe-Liu metric [25] provide guidelines that drive mesh quality [26, 27] especially as related to the convergence of a given simulation [28].

For the purposes of this work, mesh quality is defined as [29]

$$\begin{eqnarray}&&\begin{array}{c}q\equiv \frac{2{r}_{{in}}}{{r}_{{out}}}=\frac{\left(b+c-a\right)\left(c+a-b\right)\left(a+b-c\right)}{{abc}}\\ =2\left(\cos \alpha +\cos \beta +\cos \gamma -1\right)\end{array}\,\end{eqnarray} \tag{ 1 }$$

where ${r}_{in}$ is the radius of the largest inscribed circle within the element, ${r}_{out}$ is the radius of the smallest circumscribed circle that encapsulates the element, and a, b, and c are the corresponding triangle line segments. These radii and line segments are shown for a right-angled isosceles triangle in the lowest section of figure 2.

Under this definition, an equilateral triangle would have a quality factor of $q=1$ and degenerate triangles would deviate from this value to zero. While the value of $q$ for a triangular element suitable for numerical analysis is subjective, a threshold of $q\gt 0.4$ shall be adopted. The quality of the meshed representation of an entire structure (e.g., a bone, muscle, etc) shall be assessed based on the worst individual value of $q.$

Figure 2 shows the resulting composite surface model with 0.67 million facets in total. Table 1 gives the tissue mesh roster along with the minimum triangle quality for every mesh. Figure 3 illustrates the model topology along with the brain composition. The brain is divided into several sub-regions including the cerebellum striatum, external cerebrum, and the 'rest of the brain'.

Zoom In Zoom Out Reset image size

Table 2 provides the quality value for each mesh as determined by equation (1) as well as the minimum mesh edge length.

For the present study, we employ frequency ranges as defined by the International Telecommunication Union and assert that the meshes created herein are suitable for applications from Very Low Frequency (VLF) (3–30 kHz) to Very High Frequency (VHF) (30–300 MHz). The flexibility to be able to apply a single set of meshes across a broad range of frequencies (thereby enabling the simulation of many relevant applications) is due to the ease with which triangular surface meshes may be manipulated. Consider as an example the requirements of a mesh as frequency increases from VLF to VHF and the corresponding wavelength decreases. If we stipulate that twenty points are necessary to capture the behavior of an electromagnetic wave, the maximum allowable edge length of any triangular facet in a mesh may be given as:

$$\begin{eqnarray}&&{maximum}\,{edge}\,{length}=\frac{\lambda }{20}\end{eqnarray} \tag{ 2 }$$

where $\lambda$ is the wavelength in the medium under consideration (free space, tissue, bone, etc). If further refinement is required, barycentric subdivision, depicted in figure 4, or a similar refinement methodology may be employed to ensure the threshold defined by equation (2) is met. This refinement process may indeed be iteratively employed to enable simulation of applications at frequencies beyond those cited above. If lower resolution can be tolerated, faster simulation execution times can be realized through triangular surface mesh decimation via quadric edge collapse [18] or another standard decimation process.

Zoom In Zoom Out Reset image size

The model has also been augmented with tissue properties mimicking the respective human body properties and adapted from the IT'IS database [30], which tabulates values of electric permittivity and permeability from 10 Hz to 100 GHz and has become the industry standard from which countless studies have referenced. For the low-frequency simulations, we use the data at 10 kHz. For high-frequency simulations, we use the data at 600 MHz.

To model TMS exposure of the mouse brain, we applied the recently developed charge-based boundary element fast multipole method or BEM-FMM [31, 32]. This method is more accurate than the first-order FEM [33]. The entire project was constructed and disseminated using the MATLAB^® platform. The model executed in 30–50 seconds on a 2.6 GHz machine for one coil position using an accurate iterative solution.

A Cool 40 Rat coil of MagVenture, Denmark with dI/dt = 4.7e7 A/m targeting the mouse brain in an anterior position is shown in figures 5(a). (b), (c) show the computed total inducted electric field just inside the brain interfaces (cerebellum + external cerebrum + the whole brain).

Zoom In Zoom Out Reset image size

As another computational example, figure 6 shows computational results for the total electric field distribution just inside mouse bone (figure 6(a)) and in two principal planes (figures 6(b),(c) given a primary induced coil field of $E=-\partial A/\partial t=100{\rm{V}}/{\rm{m}}$ in the anterior-posterior direction. A custom coil design is located posterior to the body center. The coil field is applied predominantly to the lower extremities of a sleeping animal. Only the average body shell/volume, skeleton, and heart were included in these computations. It follows from figure 6(a) that the field just inside bone is on the order of 25%–50% of the primary coil field, which is a useful initial estimate. Figures 6(b),(c) also show that the electric field in soft tissues (and heart) is relatively small which indicates one potential advantage of magnetic stimulation compared to electrical/ultrasound stimulation.

Zoom In Zoom Out Reset image size

A typical high frequency application of the computational mouse model is modeling and design of a radio frequency (RF) coil for MRI studies (cf [34, 35]). Here, we first consider a simple ring-shaped capacitively loaded surface RF coil resonating at 600 MHz shown along with the mouse model in figure 7(a). The simulation engine is Ansys Electronics Desktop (HFSS) 2023 R1 which was used to compute the fields induced in and around the mouse by the RF surface coil as well as tissue loading effects on the coil circuit network parameters (port S-parameters).

Zoom In Zoom Out Reset image size

An electromagnetic and circuit co-simulation methodology implemented as a linked finite element method (FEM) model and linear network (S-Parameter) circuit model is used to tune and match the RF coil to operate at the target Larmor frequency of 600 MHz.

Figure 7(b) shows the equivalent circuit and a co-simulation linear network circuit analysis setup for coil tuning and matching. It includes four tuning capacitors and one matching capacitor. Figure 7(c) is the result of the port matching of the surface coil placed above the mouse's head and loaded with the head tissues. Figure 8 shows the induced electric a) and magnetic b) field magnitude plot through a cross-section of the model at 600 MHz when 2.5 mW of incident power is applied to the input of the coil. Since meshes that were constructed for this work are entirely surface based shells, consisting solely of contiguous, high-quality triangles that define the mesh surface, Ansys Electronics Desktop HFSS, the software package used to solve Maxwell's equations in 3D throughout the solution domain, uses these surfaces as a basis from which to construct tetrahedral volume meshes. Ansys Electronics Desktop HFSS then performs an evaluation of the gradient of the fields following solution to check if there are higher than desired changes across elements. If this is the case, Ansys HFSS uses its inherent tetrahedral mesh refinement to subdivide tetrahedral elements in areas of high gradients until a user defined threshold is obtained. This results in a much smaller minimum (tetrahedral) edge length employed by the FEM based solver, as shown in the rightmost column of table 2. After 6 iterations of tetrahedral mesh refinement, the final mesh utilized 4.3 million tetrahedra, required 44.5GB RAM, and had a run-time of 7 h and 6 min using 6x cores of an Intel i7–11850H CPU. Simulation convergence was evaluated by ensuring subsequent adaptively refined meshes of increasing density did not change any component of the solved S-parameter matrix by more than 0.01. The final convergence was 0.0038 after 6 solved solution passes.

Zoom In Zoom Out Reset image size

Another example shown in figure 9 is a design and tuning a 2-port 16 rung high pass birdcage rodent coil at 600 MHz. The coil has the diameter of 40 mm (the shield diameter is 44 mm) and the coil length is 30 mm. Figure 9(a) shows the mouse model with the MRI birdcage rodent coil in Ansys Electronics Desktop HFSS and a co-simulation linear network circuit analysis model for coil tuning and matching. Figure 9(b) demonstrates port matching of the birdcage coil with the mouse model inside.

Zoom In Zoom Out Reset image size

Whenever medical images are used to construct models, regardless of the simulation methodology or application, one must always consider the anatomical accuracy of the resulting meshes, as this gives the modeler a feel for how realistic the results are. Unfortunately, there is no particularly concise process or metric to accomplish this comparison. However, given how good the human eye is at locating differences between two sets of images, we can superimpose the meshes onto the original medical source data to identify any discrepancies and to ascertain if these discrepancies should be addressed.

Figure 10 provides two examples where the meshes constructed using the processes described herein are overlayed on top of the original Digimouse cryosection images. On the left of figure 10, cryosection slice 54 is used as a reference; on the right of figure 10, cryosection 101, which is higher in the axial (transverse) plane of the mouse, is the reference. We see an overall agreement between the data and the meshes with some minor variations likely due to two factors. First, the segmentation masks employed were constructed from a combination of registered computed tomography (CT), positron emission tomography (PET), and cryosection data [9]. The comparisons shown in figure 10 are lacking this other source data, but still provide a strong indication that the mesh accuracy is high. Second, some misalignments due to both cryosection slice preparation and the fact that the animal may not have been fully frozen during CT and PET imaging were noted by the original research team [9]. While techniques including piecewise volume registration were employed to properly align the various data types, it is not possible to entirely rule out some minor misalignments. However, the end results show an excellent overall agreement.

Zoom In Zoom Out Reset image size