Computerized treatment planning is routinely used in cranio-maxillofacial applications. For these therapeutic applications, it is of critical importance to have a precise model of the structures in question. Improved imaging techniques and advances in software engineering have moved three-dimensional (3D) computer models from the research and development area into routine clinical application. The importance of high-resolution source imaging is well understood by surgeons. The influence of image processing is poorly understood in the surgical community and we hypothesize that this may be a source of significant error. We evaluated the workflow for creating a virtual model using computed tomography data, and the impact that image processing decisions have on final virtual model accuracy. We chose to create a model of the dental surface since it is one of the most complex structures in the area. Individual image processing steps are explained and the magnitudes of their influence on model quality are demonstrated and compared. This study demonstrates that inappropriate image processing can introduce errors of similar magnitude as the use of inadequate source data. Finally, the study shows that errors caused by inappropriate image processing amplify the inaccuracies of low-resolution source imagery and eliminate the benefits of high-resolution source imaging.
Computerized treatment planning is routinely used in dental and in cranio-maxillofacial applications. For these therapeutic applications, it is of critical importance to have a precise model of the structures in question. This is especially important with regard to the facial bones, the occlusal surfaces, and to the inter-occlusal relationship between upper and lower jaw.
Advanced imaging techniques, software, and computerized manufacturing techniques have made three-dimensional (3D) computer models available not only for research and development, but also for routine clinical applications. They could potentially replace classical plaster models, allowing new options for occlusal analysis and treatment workflows. The generation of a virtual model requires imaging and subsequent image processing steps. Image capture is performed by direct scanning of the patient’s occlusion or scanning of a plaster cast using either tomographic imaging modalities or optical surface scanning. The image data is then computed to produce a 3D visualization.
Computed tomography (CT) and cone beam computed tomography (CBCT) imaging capture a consecutive series of two-dimensional (2D) images containing digitized tissue densities of the scanned slice. The data stack is stored as single image slices in DICOM (Digital Imaging and Communication in Medicine) format.
Optical surface scanning captures 3D surface information via the generation of a point cloud rather than the interpolation of tissue densities. These points are extrapolated to a 3D surface mesh to form the shape of an object. In contrast to CT or CBCT, optical scanning only captures surface contours, without any subsurface detail.
CT and CBCT imaging produces volumetric data that subsequently needs to be turned into a surface mesh. In both X-ray based imaging and optical surface scanning, the final result is a 3D surface mesh of the imaged structure. The surface meshes can then be merged together to create a composite skull model. This technique relies on both the accuracy of the models and the precision of the incorporation procedures.
When considering a virtual model, the clinician must be aware of the need for adequate source imaging and the image processing steps required to create the final model. Some of these steps are performed by user interaction, while some are ‘hidden’ procedures within the software. The impact of the image processing steps may subtly influence the accuracy of the final model. The generation of a virtual model from either 2D or 3D data is based on well-established imaging and image processing techniques. Surgical planning requires models of sufficient detail to enable the virtual plan to perform well when put into reality. Inadequate models may result in misleading patient assessments, inappropriate treatment, and unsatisfactory clinical outcomes. Conversely, accurate virtual occlusion modelling may enable more efficient and effective workflows for patient assessment and treatment.
The objective of this study was to show the importance and evaluate the impact of different imaging modalities and image processing steps on the accuracy of the resulting virtual model, and to identify and quantify the most critical steps to ensure accurate model generation. Further goals were to educate clinicians of the implications of these computational procedures and to discuss the clinical impact of the results. We report on the use of high-resolution CBCT imaging and on image processing parameters relevant to creating an accurate virtual model of the occlusal surface.
Materials and methods
A standard dental plaster model was scanned using an experimental high-resolution micro-CBCT scanner (SCANCO Medical AG, Brüttisellen, Switzerland) with an isotropic image resolution of 0.082 mm using standard CBCT image acquisition parameters (60 kVp/40 keV (900 μA)). The DICOM information was transferred to a desktop computer running Windows XP (Microsoft Corporation, Redmond, CA, USA) and post-processed in Amira, a commercial software package for image visualization and data analysis (Visage Imaging GmbH, Berlin, Germany), and Geomagic Studio, a reverse engineering software package (Geomagic Corp., Research Triangle Park, NC, USA).
Image processing and analysis
The workflow for creating a virtual model was evaluated using the CBCT data and divided into six major image processing steps. The first three steps, evaluating the impact of the image resolution (step 1), thresholding procedure (step 2), and smoothing of the labelled threshold data (step 3), are performed on the 2D CT data before surface model generation, as presented in Fig. 1 . Then a 3D triangulated surface model is generated to evaluate the effect of different model generation techniques (step 4), smoothing the surface model (step 5), and reducing the number of triangles of the surfaces (step 6) shown in Fig. 2 .
Step 1—image resolution
A gold standard for step 1 (GS1) CT scan was captured with a 0.082 mm resolution. To simulate lower resolution CT acquisition, the gold standard image data stack was duplicated and down-sampled to a variety of resolutions of interest. A total of eight down-sampling procedures were used to produce the following different image resolutions: four copies were down-sampled in an isotropic manner to represent an image resolution of 0.1–0.4 mm. A further four datasets were down-sampled to represent non-isotropic scanning, all having the same image resolution of 0.4 mm in the x and y axes representing the transverse plane, but varying image resolutions ranging from 0.4 mm to 0.8 mm in the z -axis or axial direction.
Duplication and down-sampling of the models was used to simulate the lower resolution imaging capture, rather than scanning at lower resolution. Simulation eliminated the introduction of error due to malpositioning of the real-world model during repeated scanning at different resolutions. This standardization meant that all deviations from the gold standard could reliably be attributed to image processing steps, rather than errors in acquisition.
From each 2D image stack, an unsmoothed model was generated using the LegoSurfaceGen procedure within Amira. This algorithm produces a surface mesh exactly identical to the voxels, without any automatic smoothing or change involved. For each model an identical threshold window was chosen by an expert and the Amira LegoSurfaceGen applied. Thresholding defines which voxel is in the region of interest on the grey scale. The end result of this procedure was eight virtual models of the same physical object generated from the source image data, of varying resolution. Each model was then compared with the GS1 model using the distance map calculation and visualization procedure in Geomagic. The distance map procedure calculated the differences between the GS1 and the down-sampled models. Distance mapping allowed calculation of maximum and mean differences between the two models ( Fig. 1 A and B).
The down-sampled models were used in the further steps to compare each model of the procedures with the models of its own identical resolution. This enabled the assessment of the impact of the further procedures without the influence of image resolution. The down-sampled models were termed ‘corresponding gold standards’ (CGS).
Every model used in the steps was made from copies from the CGS in order to have them in the exact same position in the virtual space. The imaging procedures were applied to the models in the identical position, thus avoiding the need for model alignment, which could lead to false deviations.
Step 2—threshold segmentation
This procedure defines the borderlines of our models. It converts the continuous greyscale 2D CT data into a binary representation of the pixels of CT data that represent the surface of the model being generated. In most cases it is difficult to decide which voxels belong in the grey values, whether they should be inside or outside the borderline. To mimic this problem, two different threshold values were defined, representing two different mesh border alternatives: a low (inclusive) and a high (exclusive) threshold value. For each of the varying resolution CT scans generated in step 1, two additional models were generated using a high and a low threshold and the LegoSurfaceGen procedure.
For each resolution, the high and low threshold models were separately compared with the medium threshold model of the corresponding resolution generated in step 1 (CGS). This enabled the assessment of the impact of the threshold procedure without the influence of image resolution. The deviations between the models were calculated using the Geomagic distance map calculation and visualization ( Fig. 1 C and D).
Step 3—smoothing 2D CT data
A predefined smoothing procedure was applied to the 2D CT data using the ‘smooth label’ function (parameters: size 5) in the Amira segmentation editor. A LegoSurfaceGen model was created from each of the smoothed 2D CT datasets. The deviations between the models were compared with the models of corresponding resolution obtained in step 1, using the Geomagic distance map calculation and visualization ( Fig. 1 E and F).
Step 4—3D model generation using the marching cubes algorithm
For each of the resolutions established in step 1, a model was created using the Amira standard SurfaceGen module, rather than LegoSurfaceGen. The SurfaceGen procedure creates meshed surface models according to the marching cubes algorithm, smoothing the voxel structure into a triangulated meshed surface topology. The deviations between the SurfaceGen model and the identical resolution LegoSurfaceGen models obtained in step 1 were compared using the Geomagic distance map calculation and visualization ( Fig. 2 A and B).
Step 5—smoothing of generated 3D models
Two different smoothing settings were chosen in Geomagic for each model obtained in step 4: setting 1, corresponding to moderate smoothing (‘smoothness level’ to a low level ‘3’ and ‘curvature priority’ to a high level ‘8’), and setting 2, corresponding to extensive smoothing (‘smoothness level’ to a high level ‘8’ and ‘curvature priority’ to a low level ‘3’). For each resolution, the models with moderate and extensive smoothing were separately compared with the original unsmoothed model of identical resolution generated in step 4. The differences between models were analyzed using the Geomagic distance map calculation and visualization ( Fig. 2 C and D).
Step 6—reducing the number of mesh triangles
In the eight surface models generated in step 4, the number of mesh triangles was calculated in Geomagic. The number of triangles for each model was reduced by 50% to represent a common complexity reduction procedure performed in model generation workflows. For each resolution, the resulting simplified model was compared with the identical resolution model obtained in step 4. The differences between models were analyzed using the Geomagic distance map calculation and visualization ( Fig. 2 E and F).
The impact of CT scan resolution and the five subsequent workflow steps was assessed ( Fig. 3 ).