Vector field analysis for surface registration in computer‐assisted ENT surgery

Abstract Background Manual paired‐point registration for navigated ENT‐surgery is prone to human errors; automatic surface registration is often caught in local minima. Methods Anatomical features of the human occiput are integrated into an algorithm for surface registration. A vector force field is defined between the patient and operating room datasets; registration is facilitated through gradient‐based vector field analysis optimization of an energy function. The method is validated exemplarily on patient surface data provided by a mechanically positioned A‐mode ultrasound sensor. Results Successful registrations were achieved within the entire parameter space, as well as from positions of local minima that were found by the Gaussian fields algorithm for surface registration. Sub‐millimetric registration error was measured in clinically relevant anatomical areas on the anterior skull and within the generally accepted margin of 1.5 mm for the entire head. Conclusion The satisfactory behavior of this approach potentially suggests a wider clinical integration.

However, convergence still needs a rough initialization. Improved convergence rates on the human femur and robustness to noise were achieved by a probabilistic variation of ICP, incorporating both positional and orientational information. 35 A combination of principalaxes-based registration and Hausdorff distance minimization 36 was used for automated matching of electro-anatomical and CT data, applicable on closed surfaces only. Manually delineated salient anatomical features 14 lead to lower target registration error (TRE) 37,38 than ICP; however, additional tuning of the anatomical features was required.
Additional information to ICP is usually coded as scalar attributes, or covariance matrices, 11 accounting for the anisotropy of localized surface points. Minimization of anisotropically weighted distances within points' Voronoi regions leads to improved convergence rates and better accuracy than the classical ICP, even in the presence of noise in time-of-flight (ToF) range data; results are highly sensitive to the choice of covariance matrices.
The coherent point drift 39 is a probabilistic approach to surface registration, implying coherent motion of the centroids of Gaussian mixture models. It generalizes well to non-rigid registration and outperforms ICP for brain shift estimation. 40 For rigid registration, however, pre-alignment of the datasets is still necessary, and the basin of convergence is limited to ±70°.
The GoICP algorithm 15  of convergence is further extended (though to a certain limit) through relaxation of the Gaussian aperture, resulting in a higher residual error.
Global convergence of GF cannot be guaranteed and depends on proper detection and weighting of the form attributes in the cost function.

| Brief overview of intraoperative acquisition methods
Various modalities have been used to acquire intraoperative surface data with an optically tracked hand-held probe, such as a mechanical pointer, 16 or ultrasound. [17][18][19][20][21][22][23][24][25] More advanced technologies, like positron imaging 42 or conoscopic holography, 43 achieve improved quality of the data and automatic removal of imaging artifacts. Common drawbacks by these acquisition methods are the requirements for an uninterrupted direct line-of-sight and maintenance of a constant angle of incidence of the scanning beam.
Laser acquisition of the skin surface [26][27][28] is prone to deviations from the pre-operatively generated model, due to skin elasticity. Intraoperative in-situ laser scanning of cartilage surface achieves high precision and accuracy 29 ; however, it is susceptible to stray-light. The localization error is strongly influenced by the angle of incidence and increases at greater depths.
ToF cameras represent a novel method, still under investigation, allowing fast and robust distance measurements on the patient. 30 Its application is still hindered by sensitivity to background light, reflections, and interference between multiple ToF devices. Hybrid methods, combining several modalities, are applied for tracking of inaccessible anatomical areas. 44 Registration of a 3D-model from multi-view stereo reconstructions of the facial relief 45 resulted in clinically relevant accuracy in robotic neurosurgery. The use of curvatures for the detection of geometric features in the generated models is based on a well-developed theory, using differential geometry. The method can be tracked back to early studies, 46 reporting registration of a CT-segmentation to traces of points, intraoperatively acquired with ultrasound.
Curvature-based features in stereo reconstructions from range and ToF images were detected and described through differential geometry and B-spline approximation for an improved accuracy and temporal stability in image-guided radiation therapy. 47 Most of the drawbacks of acquisition methods, eg, changes of patient's anatomy, due to skin-shifts after anaesthetization, can be overcome by A-mode ultrasound, 17-25 enabling intraoperative scanning of the bone surface and thus a rigid body registration, characterized by the highest clinical accuracy. The method is especially suited for navigated surgery of the head, where the relevant anatomy is confined within the skull. The ultrasound beam propagates through soft tissue and thus eliminates the requirement for surgical exposure of the scanned bone surface. As mentioned above, the main challenges by this approach come out of the drawbacks of the optical tracking of a hand-held ultrasound probe.
1.4 | Basis outline of the paper 1. In this contribution a major step towards global convergence of surface registration has been made. A binary energy function is minimized, introducing novelty methods, such as: • Instant center of rotation in the parameterization of the energy function; • Vector field analysis (VFA) for the detection of characteristic points in the optimization.
As to our knowledge, the above-mentioned techniques are unique and unprecedented in the existing literature.
2. The surface registration algorithm was validated on a skull phantom and a test bed, developed for the specific application in navigated ENT-surgery. Intraoperative data acquisition was automated with A-mode ultrasound in the context of registration on the posterior skull.
3. An intuitive tool was developed for visualization of the TRE on a color-coded distance map.  Figure 3 shows the 3D-models, generated from CT-data (left) and from the A-mode ultrasound acquisition (right). The characteristic anatomical relief features (the Lambda fissure and protuberantia occipitalis externa) are manually marked with form attributes.

| Surface registration through vector field analysis
The anatomical structures of the posterior skull ( Figure 3) were utilized for surface registration. Unique form attributes were assigned to the 3D points, belonging to each of the sutures of the Lambda fissure and to the protuberantia occipitalis externa, respectively. A form attribute of zero was assigned to the rest of the points. Extending the Euclidean coordinates with a form attribute coordinate, two sets of four-dimensional points were defined: q n x n ; y n ; z n ; a n ð and q n is formulated as The shortest distance (3)  . Assuming a spring constant k = 1, the total potential energy of the system in the initial position is By proper scaling of the coordinates, through tuning the Gaussian aperture σ, the range of convexity and differentiability of Equation 5 can be extended. 41 Expressing d n 2 from Equation 5, the squared closest distance becomes The real positive constant σ 2 in front of the brackets in Equation 6 can be omitted without loss of generality and influence on the convergence properties. Then, a Gaussian energy function, expressing the total potential energy of the system, is formulated as Our registration method is based on Boughorbel et al, 41 Respectively, the total potential energy of the system is Applying linear filtering on the potential field, 51 defined in Equation 8, its gradient is computed as Equation 10 defines a vector field at the points in the moving dataset.
The negative gradients are interpreted as forces of attraction, tending to fit the moving dataset onto the fixed dataset, expressed through the vector function: The curl 52 of the vector force field is computed as The symbol "×" in Equation 12 indicates the cross-product of two vectors. Analysis of the vector force field between the preoperatively and the intraoperatively generated 3D-models allows adaptation of the algorithm to the specifics of the clinical application through the depiction of the optimal center of rotation in the optimization process. A detected point with zero curl (a vortex) would gain no (rotational) velocity by the rigid motion, initiated through the vector force field.
Thus, it would act as an instant center of rotation. 50 Under the influence of force vectors, the moving dataset would tend to rotate around its instant center of rotation. The latter was used in the parameterization of the transformation for minimization of the energy function (7).

| Implementation details
The energy function (7)  The quasi-Newton method 33,41 was preferred to other gradientbased optimization methods, due to its efficient function evaluation at various positions in the parameter space, achieved through an approximation of the Hessian matrix (second-order partial derivatives). 9,33 Further, a "backtracking" strategy 33 was applied in determining the correct step-length in the gradient (Newtonian) direction, thus enabling a "stepping-out" of local minima.   Figure 5, pointed by the arrow); however, they did not influence the overall convergence rate.

Surface registration was implemented as a plug-in
Even starting from a local minimum leads to global convergence due to the "backtracking" 33 technique of the optimizer.

| Registration
Registrations from arbitrary initial positions (including the local minimum of E σ , Figure 5 The initial positions, listed in Table 1, converged to local minima of the GF algorithm (original formulation). Starting from the same initial positions, VFA achieved correct registrations. Table 1  Registrations, using VFA, were performed with datasets of a CT segmentation and an A-mode ultrasound scan of the skull phantom, with varying numbers of points with assigned form attributes. Table 2 contains the number of iterations and the times, needed for reaching a visually correct alignment (coarse registration) and for reaching the global minimum of the energy function (fine registration). All registrations were started from initial position 4 in Table 1.
Minimization of the energy function (7) with the quasi-Newton method for datasets with 411 and 315 spatial points is shown in

| Registration accuracy
For qualitative evaluation of the clinical accuracy, the registration with the proposed method is visualized in Figure 8. Two 3D models from CT segmentation of the skull phantom are in initial position 4 from     parameter space is also used by the GoICP algorithm, 15 where global convergence is achieved through a branching and nesting strategy.
The correct bounding (computed for pure rotation) is decisive for convergence.
GoICP, 15 like the classical ICP, 5 minimizes a cost function in Euclidean space. We performed minimization in scale space, where the spatial distances were robustified through the addition of form attributes. The convexity of the binary energy function in scale space facilitates optimization, while in the neighborhood of the registered position its values approximate Euclidean distance.
As the registration datasets are acquired on the posterior skull, the registration error would tend to increase in the direction of surgical areas on the anterior skull, due to a lever effect. 38 Further validations, using different imaging modalities (eg, 3D-reconstruction from a video stream) on other anatomical areas, such as patient's face, are foreseen.
As already mentioned, the purpose of this contribution is to show the feasibility of the method, allowing a unique identification of the global minimum for intraoperative registration of patient's surface data to their preoperative radiological datasets.

| CONCLUSION
An innovative approach for surface registration was successfully validated on an experimental test bed with mechanically positioned Amode ultrasound. It proved to be suitable for clinical application in navigated ENT surgery and to be generalized over other surgical and imaging domains. The automatic and reliable patient registration equipped with intuitive guiding means is assistive to the surgeon and facilitates treatment quality.