Finite-Element Based Image Registration for Endovascular Aortic Aneurysm Repair

In this paper we introduce a new method for the registration between preoperative and intraoperative computerized tomography (CT) images used in endovascular interventions for aortic aneurysm repair. The method relies on a 3D finite-element model (FEM) of the aortic centerline reconstructed from preoperative CT scans. Intraoperative 2D fluoroscopic images are used to deform the 3D FEM and align it onto the current aortic geometry. The method was evaluated on clinical datasets for which a reference CT scan was available to evaluate the registration errors made by our method and to compare them with other registration methods based on rigid transformations. Errors were estimated based on the predicted locations of landmarks positioned at different branch ostia. It appeared that our method always reduced the registration errors of at least 20% compared to gold standard 3D rigid registration and permitted to reach a global precision of 3.8 mm and a renal precision of 2.6 mm, which is a significant improvement compatible with surgical specifications. Finally, the major asset of our method is that it only requires one fluoroscopic intraoperative 2D image to perform the 3D non-rigid registration. This would reduce patient irradiation and cut the costs compared to traditional methods.


Introduction
Cardiovascular diseases are one of the first causes of death in industrialized countries, with 45,000 cases per year in the United States related to abdominal aortic aneurysms (AAA) [1]. Minimally invasive procedures known as Endovascular Aneurysm Repair (EVAR) procedures are particularly well adapted to the treatment of AAA and are conducted daily in modern hospitals. Stent-grafts are inserted through the femoral artery and are deployed using assistance of 2D angiography and fluoroscopy. The intraoperative mortality of EVAR is lower (1.5%) than open surgical repair (4.6%) [2]. However, the aneurysm-related mortality in the long term is similar [3,4]. Sophisticated fenestrated stent grafts have to be used to repair AAA expanding beyond ostia of renal arteries. In that case, the stent-graft should be positioned as precisely as possible, with an error less than 3 mm [5], to facilitate catheterization of the renal arteries. Ostia have a diameter of about 5-7mm and sometimes overlap with secondary renal arteries. Difficulties in positioning the stent graft with fenestrations sitting precisely in front of the ostia may lengthen the surgery, inducing over-irradiation and excessive injection of contrast media. Complications may occur because only 2D imaging is available intraoperatively, whereas a 3D reconstruction would be needed for optimal positioning of stent-grafts.
The key step in image-guided minimally invasive surgery is the registration of pre-and intra-interventional data [6]. There are many 2D/3D registration methods that can be applied to different surgical specialties, such as intensity-based methods using B-splines [7] or convolutional neural network regression [8] for bone registration. A method using a two-step particle filter and ultrasound scanning was proposed by Chen et al. [9] for the registration of flexible interventional catheter. Fusion imaging is commonly used to align the preoperative 3D CT scan onto intraoperative fluoroscopy images [10]. The benefits of fusion imaging are the reduction of contrast agent injection and shortening of procedure duration [11][12][13]. Still, fusion imaging usually requires costly hybrid rooms with fixed imaging systems which automatically adjust the position of the operating table. A hybrid operating room is an advanced surgical space that combines a traditional operating room with a suite of image-guided procedures. This combination enables performing complex and advanced surgical procedures. Conversely, mobile C-arms, which are still very commonplace, do not permit direct fusion imaging. However, fusion imaging in hybrid rooms also increases patient's irradiation as a 3D scan is required at the beginning of the procedure.
The fully automated registration method proposed by Varnavas et al. [14] achieves a satisfactory registration rate thanks to two algorithms based on normalized features and multiple vertebra registrations. However, arteries can be naturally deformed, depending on the patient's body position or artificially modified when the endovascular tools are inserted [15,16]. Therefore, standard 3D/3D fusion imaging based on rigid registration cannot avoid large errors in some cases. Several 2D/3D registration methods based on preliminary computations were proposed. Kaladji et al. [17] have shown the benefit and feasibility of 2D/3D registration under clinical conditions, which reduces the use of contrast agent. They introduced a rigid registration method and proposed to pre-calculate the non-rigid deformations that may occur within the aorta during the procedure to provide a more accurate registration if necessary [18,19]. Gindre et al. [19] showed the possibility of predicting large deformations of vascular geometry with numerical simulations and finite element analysis. Mohammadi et al. [20] proposed precomputing the deformation of the vascular tree during EVAR by simulating the positioning and deployment of stent-grafts. They achieved good accuracy in the positioning of tools and iliac ostia, while avoiding the constraint of computation time. However, pre-calculating model deformations requires perfect prediction of clinical workflows. Any change of tool or any unexpected event that may occur during surgery makes the model obsolete. Precomputing deformations also means not using all the available information from intraoperative imaging. Other studies introduced real-time intraoperative 3D AAA non-rigid registration methods.
Liao et al. [21] presented a method based on two fluoroscopic images and constrained a deformation model with neighboring smoothness terms, length preserving and 2D distance maps. However, this method requires several 2D views to achieve accurate 3D registration. Another method using two fluoroscopic images was proposed by Toth et al. [22], based on skeleton-based deformations which map the 3D pre-operative AAA shape onto intraoperative locations of inserted devices and landmarks. Raheem et al. [23] described a deformable 2D/3D registration method using Thin Plate Spline (TPS) with landmarks. The method also requires at least two fluoroscopic images. This involves stopping the procedure temporarily to take two images, with different angles of incidence, thereby increasing the duration of the procedure. To reduce the duration of interventions, a method based on a single view would be more suitable. Groher et al. [24] proposed a 2D/3D vessel registration method based on the Iterative Closest Point (ICP) approach, TPS, smoothness and length regularization with a single fluoroscopy input. This method shows interesting results with significant improvement of the position error. However, the proposed method is not guaranteed to converge to the right solution in some specific cases, for example when vessels are deformed by tools insertion and does not rely on a mechanically realistic model. It also requires about 5 min to run, which can be too long for surgical purposes. Zheng et al. [25] proposed to graph-match 3D preoperative AAA skeleton with 2D intraoperative skeleton, before performing the registration while respecting skeleton length and smoothness. This method yields interesting results but is sensitive to topology variance. In summary, current 3D/3D fusion imaging systems requiring a hybrid room usually miss non-rigid registration whereas 2D/3D registration methods need several 2D views, a long computation time or are not able to deal with specific cases.
In order to address these issues, we propose a new 2D/3D registration method making a good trade-off between performance, precision and simplicity. The method requires only one x-ray projection to update a finite element model and is able to run on equipment with limited resources. The simplicity of the method's use and the short calculation time enable its integration into an operating workflow without disrupting it and does not require any new hardware. After a comprehensive presentation of the method, we make a proof of concept on clinical datasets for which a reference CT scan was available to evaluate the registration errors made by our method and to compare them with other registration methods based on rigid transformations.

Data Acquisition
The data available at the beginning of the non-rigid registration method includes the 3D geometry of the aorta extracted from the preoperative scan and the intraoperative 2D images. The geometry extracted from the preoperative scan is an accurate representation of the patient's aorta in the three directions of the patient reference frame. However, the preoperative aortic geometry may be different from the intraoperative aortic geometry. Actually, the aortic configuration can undergo deformations between the day of the preoperative scan and the day of the surgical procedure. This may be simply induced by a different position of the patient or due to the insertion of guidewires. Intraoperative images are a projection of the patient's aorta and thus contain only two-dimensional data, along u and v ( Figure 1). Depth, along the projection axis of cone beam, is unknown. Nevertheless, this information is up to date unlike the preoperative scan. The objective of the non-rigid registration method is therefore to update outdated 3D data with up-to-date 2D data. For this purpose, the 3D geometry is implemented in a mechanical model, which is then solved according to the 2D information. Thus, the model is perfectly up to date according to the projection plane of the 2D images and the mechanical model calculates the resulting deformations according to the depth. In addition, the method must be performed in real time and therefore should rely on a simplified mechanical model with a reduced calculation time.
Modelling 2020, 2, FOR PEER REVIEW 3 variance. In summary, current 3D/3D fusion imaging systems requiring a hybrid room usually miss non-rigid registration whereas 2D/3D registration methods need several 2D views, a long computation time or are not able to deal with specific cases. In order to address these issues, we propose a new 2D/3D registration method making a good trade-off between performance, precision and simplicity. The method requires only one x-ray projection to update a finite element model and is able to run on equipment with limited resources. The simplicity of the method's use and the short calculation time enable its integration into an operating workflow without disrupting it and does not require any new hardware. After a comprehensive presentation of the method, we make a proof of concept on clinical datasets for which a reference CT scan was available to evaluate the registration errors made by our method and to compare them with other registration methods based on rigid transformations.

Data Acquisition
The data available at the beginning of the non-rigid registration method includes the 3D geometry of the aorta extracted from the preoperative scan and the intraoperative 2D images. The geometry extracted from the preoperative scan is an accurate representation of the patientʹs aorta in the three directions of the patient reference frame. However, the preoperative aortic geometry may be different from the intraoperative aortic geometry. Actually, the aortic configuration can undergo deformations between the day of the preoperative scan and the day of the surgical procedure. This may be simply induced by a different position of the patient or due to the insertion of guidewires. Intraoperative images are a projection of the patientʹs aorta and thus contain only two-dimensional data, along u and v ( Figure 1). Depth, along the projection axis of cone beam, is unknown. Nevertheless, this information is up to date unlike the preoperative scan. The objective of the nonrigid registration method is therefore to update outdated 3D data with up-to-date 2D data. For this purpose, the 3D geometry is implemented in a mechanical model, which is then solved according to the 2D information. Thus, the model is perfectly up to date according to the projection plane of the 2D images and the mechanical model calculates the resulting deformations according to the depth. In addition, the method must be performed in real time and therefore should rely on a simplified mechanical model with a reduced calculation time. Although their nature and origin are still subject of controversy, natural deformations undergone by the aorta between consecutive CT scans acquired several months apart may be rather large and they constitute a significant case study for the image registration problems considered in this paper. Then we collected retrospectively 3D CT scans of 7 patients who further underwent EVAR (Table 1). These 7 patients all had two preoperative scans available before their EVAR procedure. The two preoperative CT scans were taken at about 7 months interval and were named scan 1 and scan 2. All subjects gave their informed consent for inclusion before they participated in the study. The study was conducted in accordance with the Declaration of Helsinki and the protocol was approved Although their nature and origin are still subject of controversy, natural deformations undergone by the aorta between consecutive CT scans acquired several months apart may be rather large and they constitute a significant case study for the image registration problems considered in this paper. Then we collected retrospectively 3D CT scans of 7 patients who further underwent EVAR (Table 1). These 7 patients all had two preoperative scans available before their EVAR procedure. The two preoperative CT scans were taken at about 7 months interval and were named scan 1 and scan 2. All subjects gave their informed consent for inclusion before they participated in the study. The study was conducted in accordance with the Declaration of Helsinki and the protocol was approved by the Ethics Committee CHU PROM 1408061. All data and images were acquired at University Hospital of Saint-Etienne in clinical settings. The maximum pixel size of all CT scans is 0.9395 × 0.9395 mm 2 and the maximum slice thickness is 2 mm.

Non-Rigid Registration (NR) and Measurements
This section is dedicated to the description of the non-rigid registration method (NR) (Appendix A) and the evaluation of the method (Appendix B).
A detailed flowchart of the proposed framework for the finite element model of the non-rigid registration method is shown in Figure 2. The initial configuration of the model is based on data extracted from preoperative CT scans, the initial centerline composed of a list of 3D points X 0 associated with a connectivity table. The model is constrained with data extracted from intraoperative images. The projection of the initial centerline is registered onto the fluoroscopy. Those 2D constraints are translated in 3D via Lagrange multipliers implemented as complementary equation to the finite element model (FEM) resolution. The aorta is finally reconstructed. We assume that the projection matrix and the patient orientation on the table are known. Thus, errors in the matrix directly affect the performance of the method. Geometric nonlinearities are handled by the co-rotational formulation used for the beam elements representing the arteries. This approach allows to extract the rotations of each element of the model explicitly and to reintroduce them in the calculation of elastic forces. The complete method is detailed in Appendix A included in supplemental materials. First, the geometry of the aorta is extracted from the DICOM stack of the preoperative CT scan with a front collision method implemented in the Vascular Modeling Toolkit (VMTK) [26]. Then, the centerline of the initial geometry is extracted using the Voronoi diagram method also implemented in the VMTK library. The centerline is represented by a matrix of 3D points X0 associated with radii of inscribed spheres and a connectivity table. This centerline is then implemented as the initial geometry of the FEM. The target data consist of a single Digital Subtraction Angiography (DSA) taken

Non-Rigid Registration (NR) Method
First, the geometry of the aorta is extracted from the DICOM stack of the preoperative CT scan with a front collision method implemented in the Vascular Modeling Toolkit (VMTK) [26]. Then, the centerline of the initial geometry is extracted using the Voronoi diagram method also implemented in the VMTK library. The centerline is represented by a matrix of 3D points X 0 associated with radii of inscribed spheres and a connectivity table. This centerline is then implemented as the initial geometry of the FEM. The target data consist of a single Digital Subtraction Angiography (DSA) taken from a single angle of incidence. This DSA is the projected shadow of the target aortic volume. The projection parameters are known. Indeed, the projection parameters are intrinsic to the mobile C-arm used for imaging and were characterized with a dedicated device (calibration target designed for mobile C-arms). The next step is a 2D/2D non-rigid registration of the initial centerline projection onto the target images, which must be translated into 3D deformations. It would have been possible to extract the target centerline directly from the projection of the target volume but the point-to-point association of two different center lines (CLs) projection creates errors and increases the calculation time, due to the lack of information out of the projection plane. The conservation of the initial centerline points throughout the process allows maintaining a perfect match between the 3D centerline points and the 2D centerline points. The points of the 2D initial centerline are back-projected to obtain back-projection lines that will serve as "rails" on which the points of the 3D initial centerline are "guided." Each point of the finite element model is free to move along its corresponding back-projection lines. The deformations needed to reach an equilibrium state are computed by the mechanical model. The centerline is imported into the finite element model as a matrix of 3D points associated with radii of inscribed spheres and a connectivity table. Each point becomes a node of the model and beam elements are created between every consecutive node pair. The elements are linked together with the connectivity table. Each node has 6 degrees of freedom U, 3 translations and 3 rotations. The central line of the aorta is modeled as a series of beam elements within an in-house finite element code written in Matlab ® . This centerline was previously extracted by image segmentation from the pre-operative CT scan. Lagrange multipliers allow new multi-freedom constraints (MFCs) to be prescribed onto the finite element model by adding equations to the stiffness matrix. In our case, Lagrange multipliers are used to constrain nodal displacements to remain along the projection lines, by adding equations to the global stiffness matrix K. Finally, the 3D volume is reconstructed with a VMTK function. 3D spheres associated with each point of the central line are merged to create the mesh corresponding to the new surface, which is consequently defined by a set of 3D points and triangular elements. Surface details of the aorta are lost but the general volumes and most importantly the position of the ostia, are well represented.

Evaluation
In this section, we evaluate and compare our 2D/3D non-rigid registration method (NR) to standard rigid registration methods currently used for fusion imaging. More details can be found in Appendix B in supplemental materials. In this study, different set of landmarks were used for different purposes, either to estimate registration errors or to perform the registration. Landmarks were preferred to other image-based methods for the evaluation of registration, although these methods could use data of the entire image. Indeed, the two images are usually acquired with two different scanners and the differences between the two images combines effects of the motion with effects of imaging parameters. It was found that using landmarks was the most straight forward and robust method to handle these issues. All the landmarks were positioned by an expert. A repeatability study was carried out on the position of the markers. The positioning error of the landmarks has been estimated to 0.2mm. The first standard rigid registration method is the 3D/3D registration or spine registration method (R3D). This method aligns two 3D volumes and consequently requires 3D imaging at the beginning of the procedure, which is only possible in hybrid rooms. A rigid transformation T is computed to register the spine lumbar vertebra of the two scanners. The second standard method is the 2D rigid registration method (R2D) performed on the aorta shape, derived from the non-rigid registration process. This method simply consists of rigid registration of the projected initial aorta on the target aorta from the DSA. This registration is representative of a 2D/3D rigid registration method and does not require hybrid rooms. These methods are not able to process non-rigid deformations. By comparing their results with the non-rigid method proposed here, we will determine whether the non-rigid aspect is relevant or not and evaluate the contribution of the mechanical model. Natural deformations undergone by the aorta between consecutive CT scans acquired several months apart constitute a significant case study for the image registration problems considered in this paper. To evaluate the performance of the centerline (CL) deformation method, we used two CT scans for each patient to simulate the clinical data. For each patient, two CT scans and consequently two 3D volumes, were available. The first scan was used as a 3D source image, equivalent to pre-operative image, from which we extracted the geometry implemented in the model. The second scan was used to simulate 2D target images. The images were created by making a 2D projection of the second scan volume, simulating the conical projection of x-rays on the patient during mobile C-arm imaging. More specifically, we simulated the whole process of image generation, including the conical projection associated with the C-arm and the physical laws governing X-rays. The target image is thus similar to clinical fluoroscopy images and equivalent to intra-operative images. Accordingly, the target 3D geometry is known, which allows us to precisely measure the performance of our registration method but also to compare it with a 3D/3D registration method for the validation purpose. This image simulation can be used to make assumptions for possible clinical application-the projection matrix is known and the artery has a perfect circular cross-section. The distance between key points of the aorta is measured to compare the accuracy of the registration methods. This distance characterizes the registration error of the method for the points of interest related to the surgery, that is, the ostia of the abdominal aorta-right and left renal, right and left iliac, superior mesenteric and coeliac arteries. For each CT scan, the 3D position of each ostium is taken 5 times. The mean point is calculated for each ostium and the distance between corresponding ostia of the reference and registered volume is measured. The measurement uncertainty simply results in the difference between the 5 measurements taken for each ostia. Figure 3 shows the global precision of each registration method, in the form of a box plot, representing the median, the first and third quartiles and the maximum and minimum values. Here and for the rest of the document, the term global precision refers to the mean registration errors for all ostia. The mean global precision of the NR method is 3.8 (SD 1.3) mm; R2D method is 5.6 (SD 2.7) mm; R3D method is 5.3 (SD 1.6) mm. Concerning only renal ostia-NR is 2.6 (SD 1.3) mm; R2D is 2.7 (SD 0.9) mm and R3D is 5.3 (SD 3.0) mm. The mean residual RMS on landmarks positioning is 0.8 mm for R3D. Figure 4 shows the precision of the three different registration methods. It represents the distance between the ostia of the registered geometries and the target ostia, for each ostium. Figure 5 shows the reconstructed aortas after the non-rigid registration, compared with the target volumes. all ostia. The mean global precision of the NR method is 3.8 (SD 1.3) mm; R2D method is 5.6 (SD 2.7) mm; R3D method is 5.3 (SD 1.6) mm. Concerning only renal ostia-NR is 2.6 (SD 1.3) mm; R2D is 2.7 (SD 0.9) mm and R3D is 5.3 (SD 3.0) mm. The mean residual RMS on landmarks positioning is 0.8 mm for R3D. Figure 4 shows the precision of the three different registration methods. It represents the distance between the ostia of the registered geometries and the target ostia, for each ostium. Figure 5 shows the reconstructed aortas after the non-rigid registration, compared with the target volumes.   and for the rest of the document, the term global precision refers to the mean registration errors for all ostia. The mean global precision of the NR method is 3.8 (SD 1.3) mm; R2D method is 5.6 (SD 2.7) mm; R3D method is 5.3 (SD 1.6) mm. Concerning only renal ostia-NR is 2.6 (SD 1.3) mm; R2D is 2.7 (SD 0.9) mm and R3D is 5.3 (SD 3.0) mm. The mean residual RMS on landmarks positioning is 0.8 mm for R3D. Figure 4 shows the precision of the three different registration methods. It represents the distance between the ostia of the registered geometries and the target ostia, for each ostium. Figure 5 shows the reconstructed aortas after the non-rigid registration, compared with the target volumes.

Discussion
The objective of this study was to introduce and to evaluate our novel non-rigid registration technique, which is based on a mechanical model. For this purpose, the performances of the method were compared to a 3D/3D rigid registration method and a 2D/3D rigid registration method, both

Discussion
The objective of this study was to introduce and to evaluate our novel non-rigid registration technique, which is based on a mechanical model. For this purpose, the performances of the method were compared to a 3D/3D rigid registration method and a 2D/3D rigid registration method, both state-of-the-art in current fusion imaging registration methods. After prior discussion with expert vascular surgeons, the acceptability value has been set below 5mm on ostia. The overall accuracy of the registration changes from low (x > 5mm) for the 2D rigid method to good (3 mm < x < 5 mm) for the non-rigid method according to the scale of Kauffmann et al. [5]. A first qualitative visual comparison between the target aortas and the reconstructed aortas ( Figure 5) gave an overview of the results-the registration in the projection plane (front view) is highly accurate. The registration outside the projection plane is globally less precise but remains acceptable (<5 mm), except for the tortuous geometry of patient 5, where accuracy on the renal arteries is acceptable but should be improved for the mesenteric and celiac trunk arteries, outside the projection plane. When comparing the overall accuracy of the different methods in Figure 3, it appears that the non-rigid method is the most accurate, ahead of 3D rigid registration and 2D rigid registration. The improvement of accuracy between these methods shows the added value of the non-rigid aspect of registration, even for an aorta subjected to low natural deformations. Moreover, the non-rigid method is more robust. The registration errors are contained in a limited range, while the errors approach 10 mm for the 2D rigid registration in some cases. However, for renal arteries, non-rigid registration, although qualified as good according to Kauffmann et al. [5], remains qualitatively below 2D rigid registration. The non-rigid and 2D rigid methods are rigidly registered in the same way on the renal arteries. Consequently, the ostia position of renal arteries is slightly modified during the deformation of the non-rigid model, mainly along the antero-posterior axis, outside the registration plane. This is one of the limitations of the method, which seeks to achieve an overall mechanical balance. Weighting the important points could solve the problem in future developments. The one view 2D/3D non-rigid method achieves a better accuracy than the rigid registration methods presented in this study. However, this simple method, with a calculation time of less than half a second using Matlab ® , could be perfectly adapted to the surgery workflow. Indeed, after discussions with vascular surgeons, the maximum allowed duration for a simulation during a surgical procedure was set to 1 min. Here, calculation time refers to computing the deformation field and to registration, including the finite element analysis. In addition, it requires only one DSA taken from a single incidence, theoretically reducing the patient's irradiation rate. The registration accuracy for both renal arteries (2.5 mm) is qualified as high according to the scale of Kauffmann et al. [5] and is therefore satisfactory for surgical use. The accuracy of the registration obtained is slightly lower than the Zheng method. However, unlike Zheng et al. [25], we do not directly compare the positions of 2D points with 3D points. We compare the position of two sets of 2D points and we perfectly know the correspondence between the projected 2D initial points and 3D initial points. Therefore, we can minimize the registration error. Moreover, our non-rigid method is based on a FEM rather than Zheng et al. [25] and Groher et al. [24]. It integrates the actual mechanical characteristics of the aorta, such as diameter, thickness and Young's modulus. Indeed, the diameter and thickness of the aorta are integrated into the stiffness matrix K via the cross-sectional area of the beam A and the cross section moments of inertia I as described in Appendix A in supplemental materials. Simulation results remain stable for Young's modulus E ranging between 7 and 13 MPa [27,28]. It is therefore possible to implement new functions in the code, such as the navigation of guidewires, in a mechanically consistent model and without significant modifications to the existing code. In summary, it is possible to simulate in real time the deformation of the aorta induced by device insertion during the intervention using the same model.
To evaluate the model with large deformations, the method was tested on an artery deformed by a stent-graft launcher. For the sake of verification, the target artery is deformed using a finite element simulation and does not correspond to actual clinical data but the imposed deformations, reviewed and validated by a vascular surgeon, are similar to a real case. As shown in Figure 6, the method is able to handle large deformations, with a mean accuracy of 6.3 (SD 2.3) mm and an accuracy for the right and left renal arteries of 3.3 and 2.2 mm respectively. Therefore, this recalibration method enables to set up a mechanical model of the aorta, updated according to operating data and suitable for clinical applications. However, more tests with large deformations should be conducted for confirmation. significant modifications to the existing code. In summary, it is possible to simulate in real time the deformation of the aorta induced by device insertion during the intervention using the same model. To evaluate the model with large deformations, the method was tested on an artery deformed by a stent-graft launcher. For the sake of verification, the target artery is deformed using a finite element simulation and does not correspond to actual clinical data but the imposed deformations, reviewed and validated by a vascular surgeon, are similar to a real case. As shown in Figure 6, the method is able to handle large deformations, with a mean accuracy of 6.3 (SD 2.3) mm and an accuracy for the right and left renal arteries of 3.3 and 2.2 mm respectively. Therefore, this recalibration method enables to set up a mechanical model of the aorta, updated according to operating data and suitable for clinical applications. However, more tests with large deformations should be conducted for confirmation. The aorta model still suffers from several limitations that can reduce, in some cases, the accuracy of the method, especially regarding boundary conditions. For example, the anchorage points of the aorta and the surrounding structures are not considered, in particular the spine, which can significantly change the mechanical behavior of the aorta [19]. This should be implemented in a future version of the method. The local stiffness of the aorta is also not considered, although it can vary significantly, for example in the case of different thrombus size or calcification, which are not included in the segmentation of the lumen. As a result, the deformations of the artery would be different, which should be taken into account in a finite element analysis. Therefore, future versions of our methodology should consider these parameters and the underlying local modifications of the stiffness applied to the simplified model of the aorta. However, as our simulations are guided by intraoperative imaging, deformations due to the presence of thrombus and calcification are already partially taken into account, which should minimize errors. The method has other limitations. It depends on the projection matrix and assumes it to be known. Its accuracy depends on the reliability The aorta model still suffers from several limitations that can reduce, in some cases, the accuracy of the method, especially regarding boundary conditions. For example, the anchorage points of the aorta and the surrounding structures are not considered, in particular the spine, which can significantly change the mechanical behavior of the aorta [19]. This should be implemented in a future version of the method. The local stiffness of the aorta is also not considered, although it can vary significantly, for example in the case of different thrombus size or calcification, which are not included in the segmentation of the lumen. As a result, the deformations of the artery would be different, which should be taken into account in a finite element analysis. Therefore, future versions of our methodology should consider these parameters and the underlying local modifications of the stiffness applied to the simplified model of the aorta. However, as our simulations are guided by intraoperative imaging, deformations due to the presence of thrombus and calcification are already partially taken into account, which should minimize errors. The method has other limitations. It depends on the projection matrix and assumes it to be known. Its accuracy depends on the reliability of the matrix. To ensure the stability and robustness of the method, it should be tested with more patients, including cases that are more tortuous. Finally, the evaluation of the methodology used fluoroscopy data taken in similar conditions. We are now considering verification with CT scans acquired at a different time and in different conditions, as our methodology was designed to handle such situations.
As current 3D/3D fusion imaging systems are not optimal for EVAR procedures, we proposed a new 2D/3D registration method making a good trade-off between performances, precision and simplicity. After a comprehensive presentation of the method, we showed its performances on clinical datasets for which a reference CT scan was available to evaluate the registration errors and to compare them with other registration methods based on rigid transformations. The study permitted to evaluate the natural deformations of the aorta observed between two consecutive scan and to differentiate the non-rigid deformations from the rigid displacements. It appeared that a non-rigid registration method is essential to register accurately ostia positions. We also proved that the method was able to handle large deformations induced by tool insertions. Therefore, our method can represent a significant improvement for fusion imaging in EVAR procedures as it has also the potential to reduce patient irradiation and, unlike current 3D/3D methods, it does not require hybrid room facilities. Combined with simulations performed at the planning stage of EVAR [29][30][31][32], this will participate in the uptake of computer simulations for assisting surgical interventions in future medicine. This part is dedicated to the description of the non-rigid registration method (NR). The aim of this method is to find the displacement U to minimize the potential energy π(U) of the aorta. The appendix is an optional section that can contain details and data supplemental to the main text.

Patents
(u,v) are the 2D coordinates of the key points extracted from intraoperative images. X 0 are the initial 3D coordinates of the key points from preoperative CT scan. U stands for the 3D displacements solutions of Equation (A1). P is the projection matrix. With: w is the strain energy density function; V the volume of the elastic body; T i the ith component of the surface traction; u i the ith component of the deformation; F i the ith component of a body force; S T the portion of the body surface over which tractions are prescribed. This minimization problem is solved with a traditional finite element model. A detailed flowchart of the proposed framework for the finite element model of the non-rigid registration method is shown in Figure 2. The initial configuration of the model is based on data extracted from preoperative CT scans, the initial centerline composed of a list of 3D points X 0 associated with a connectivity table. The model is constrained with data extracted from intraoperative images. The projection of the initial centerline is registered onto the fluoroscopy. Those 2D constraints are translated in 3D via Lagrange multipliers and the related equations M are implemented as complementary equation to the finite element model (FEM) resolution (Equation (A12)). The model properties include the stiffness matrix K, involving 2 constants (Young's modulus E, Poisson's ratio ν) and 2 dimensional parameters (the diameter of the aorta d and the wall thickness e). The aorta is finally reconstructed. We assume that the projection matrix and the patient orientation on the table are known. Thus, errors in the matrix directly affect the performance of the method. Geometric nonlinearities are neglected in the finite element model. First, the geometry of the aorta is extracted from the DICOM stack of the preoperative CT scan with a front collision method implemented in VMTK [26]. This method consists in defining several seeding points on the image. Fronts propagate from the seeds with their speed proportional to the voxel intensity. The region where the two fronts collide is then the initial model. The segmented aorta is extracted as a 3D mesh composed of triangular elements. Then, the centerline of the initial geometry is extracted using the Voronoi diagram method implemented in the VMTK library. Centerlines are determined as weighted shortest paths traced between two extremal points. In order to ensure that the final lines do not lie anywhere in space, they are bound to run on the Voronoi diagram of the vessel model. For every point belonging to the Voronoi diagram, there is a sphere centered in that point that is a maximal inscribed sphere (the information relative to the radius is therefore defined everywhere on the Voronoi diagram). This initial centerline is then implemented in the FEM. The centerline is represented by a matrix of 3D points X 0 associated with radius of inscribed sphere and a connectivity table. This centerline is then implemented as the initial geometry of the FEM.

Model Constrains
The target data consist of a single Digital Subtraction Angiography (DSA) taken from a single angle of incidence. This DSA is the projected shadow of the target aortic volume. The projection parameters are supposed to be known. A projection matrix is composed of the intrinsic parameter matrix C which contains the properties of the camera.
α α and α v . express the focal length, in number of pixels (respectively in horizontal and vertical directions); u 0 and v 0 are the coordinates of the principal point, expressed in the pixel reference frame. The position and the orientation of the camera are expressed by a 3D rotation matrix R and a 3D translation matrix t. The 3 × 4 projection matrix P is: It would have been possible to extract the target centerline directly from the projection of the target volume but the point-to-point association of two different CLs projection creates errors and increases the calculation time, due to the lack of information out of the projection plane. The conservation of the initial centerline points throughout the process allows maintaining a perfect match between the 3D centerline points and the 2D centerline points. Thus, the first step is a 2D/2D non-rigid registration of the projection of the initial centerline on the target images.
The geometry implemented in the FEM model, that is the centerline of preoperative aorta, is then projected according to the projection parameters. First, a 3D/2D rigid registration is performed between the projection of the initial centerline and the intraoperative shadow of the target aorta from the DSA, to initiate the non-rigid registration. This was done manually for this study but many studies presented 3D/2D rigid registration method, for example Groher et al. [24]. Then the centerline is projected again, to obtain a new set of projected points S CL (x cl ,y cl ), to perform the non-rigid registration. The rigidly registered centerline is then interpolated to be aligned onto the shape of the target, in 2D. We do not use the CL of the 2D target, which is not extracted. A set of landmarks is positioned on the projection of the preoperative aorta, named S1 2D (x1 2D ,y1 2D ). Another set of landmarks is positioned on the target image at the same physiological points, named S2 2D (x2 2D ,y2 2D ). S1 2D and S2 2D are used to obtain two displacement vectors: The data are interpolated with the griddata Matlab ® function, using biharmonic spline interpolation, based on Green functions. It may be written: v q = griddata x, y, v, x q , y q where (x,y) are the former coordinates related to the sample values v; v q are the updated values interpolated on the new coordinates (x q ,y q ). So here: 2D , v y , x cl1 , y cl1 (A7) v xi and v yi are the interpolated displacement vectors of the new centerline 2D points. Green functions of the biharmonic operator, in one and two dimensions, are used for minimum curvature interpolation of irregularly spaced data points. The interpolating curve is a linear combination of Green functions centered at each data point [33]. The result is a 2D/2D non-rigid registration of the projection of the initial centerline onto the target, which must be translated in 3D deformations.
Appendix A.1.4. Back Projection Lines Then, the 2D registration must be transferred to the 3D model. The points of the 2D initial centerline are back-projected to obtain back-projection lines that will serve as "rails" on which the points of the 3D initial centerline are "guided." Each point of the finite element model is free to move along his corresponding back-projection lines. The deformations necessary to reach an equilibrium state are computed by the mechanical model. Displacements are assigned to the 3D model via Lagrange multipliers in addition to the stiffness matrix as shown previously. The equations of the projection lines result from basic Cartesian geometry. The X-ray source is the origin of the global coordinate system. Therefore, all projection lines have a common point O (0,0,0). The projection parameters are known, along with the position of the target image in 3D space, representing the layout of the X-ray source and the flat panel detector. Thus, each target point P(x,y,z) of the image is associated with a back-projection line running through this point and through the origin. The normalized vector of the line is therefore: The points of the model are forced to follow these lines, prescribe as constraints using Lagrange multipliers.
matrix Ke calculated in local coordinates that relates angular and spatial position of each end of a beam element to the forces and torques applied to them [34]: where E is the Young's modulus and ν the Poisson's ratio; Iy and Iz are cross-section moments of inertia; A is the cross-sectional area of the beam, l its length; ϕy and ϕz represent shear deformation parameters and are defined as ϕy = 12EIz GAsyl 2 and ϕz = 12EIz GAsyl 2 with A sy and A sz the shear area in the y and z directions. Elementary stiffness matrices Ke of each beam elements are then assembled together into a global stiffness matrix K. Since K e is initially calculated in local coordinates, a transformation matrix is applied to express Ke in the global coordinate system. As each element DoFs are in global coordinates and as the topology is a simple sequence of segments, the resulting global stiffness matrix K is block-tri-diagonal. This type of system can be solved very efficiently [31]. This centerline is the simplified representation of the aorta and will be deformed during the next steps according to target data. The boundary conditions are limited to zeroing the translations of a single reference point positioned between the renal ostia.  Finally, the 3D volume is reconstructed with a VMTK function. 3D spheres associated with each point of the central line are merged to create the mesh corresponding to the new surface, which is consequently defined by a set of 3D points and triangular elements. Surface details of the aorta are lost but the general volumes and most importantly the position of the ostia, are well represented.

Appendix B.1. Evaluation
The previous part presents the non-rigid registration method based on a FEM of the centerline of the aorta. Three registration methods are now compared in this section to evaluate the performance of the non-rigid registration method.
The first method is the main 2D/3D non-rigid registration method (NR) itself. This method is then compared to rigid registration methods currently used for fusion imaging. The surgeons we are collaborating with currently use these methods; therefore, we use them as a standard to evaluate the performance of our approach. The first rigid registration method is the 3D/3D registration already described below as the spine registration method (R3D). This method aligns two 3D volumes and consequently requires 3D imaging at the beginning of the procedure, which is only possible in hybrid rooms. The second one is the 2D rigid registration method (R2D) performed on the aorta shape, derived from the non-rigid registration process. This registration is representative of a 2D/3D rigid registration method and does not require hybrid rooms. These methods are not able to process non-rigid deformations. By comparing their results with the non-rigid method proposed here, we will determine whether the non-rigid aspect is relevant or not and evaluate the contribution of the mechanical model. Natural deformations undergone by the aorta between consecutive CT scans acquired several months apart constitute a significant case study for the image registration problems considered in this paper. To evaluate the performance of the centerline (CL) deformation method, we used two CT scans for each patient to simulate the clinical data. The first CT scan represents the preoperative CT scan used to extract the 3D geometry implemented in the model. The second CT scan was used to simulate intra-operative 2D AAA target fluoroscopies. Thus, the target 3D geometry is known, which allows us to precisely measure the performance of our registration method but also to compare it with a 3D/3D registration method. This image simulation can be used to make assumptions for possible clinical application: the projection matrix is known and the artery has a perfect circular cross-section.