Image distortion correction and robust phantom detection

Abstract

The invention relates to a method for detecting a phantom, comprising the steps of: arranging a phantom with respect to an object, acquiring at least one image of said object by means of an x-ray apparatus, such that the image contains projections of the object and projections of at least three first calibration fiducials of the phantom, detecting the projections of the at least three first calibration fiducials in said at least one image, and establishing a correspondence between the 2D image coordinates of said projections of the at least three first calibration fiducials and the 3D coordinates of said at least three first calibration fiducials in a local coordinate system of the phantom for computing the projection matrix at least up to a scale factor.

Claims

1. A method for detecting a phantom, comprising the steps of: arranging a phantom with respect to an object, the phantom comprising a plurality of first calibration fiducials in a first plane, acquiring at least one image of said object by means of an x-ray apparatus, such that the image contains projections of the object as well as projections (IM.sub.i.sup.1, IM.sub.i.sup.2, IM.sub.i.sup.3) of at least three first calibration fiducials, detecting the projections of the at least three first calibration fiducials in said at least one image, and establishing a correspondence between 2D image coordinates (I.sub.x, I.sub.y) of said projections of the at least three first calibration fiducials and 3D coordinates (x, y, z) of said at least three first calibration fiducials in a local coordinate system of the phantom for computing a projection matrix (P) at least up to a scale factor (α), which projection matrix (P) relates said 2D image coordinates (I.sub.x, I.sub.y) to said corresponding 3D coordinates (x, y, z) in a local coordinate system of the phantom according to α[I.sub.x, I.sub.y, 1].sup.T=P[x,y,z,1].sup.T, wherein the at least three first calibration fiducials (M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3) corresponding to said detected projections are successively arranged on an associated line of said first plane, wherein the ratio r.sub.i=|M.sub.i.sup.1M.sub.i.sup.2|/|M.sub.i.sup.2M.sub.i.sup.3| is different from 1.

2. The method according to claim 1, wherein said correspondence is established using a pinhole camera model for a projection of the x-ray apparatus with a reference position corresponding to extrinsic parameters R.sub.0 and T.sub.0 according to
α[I.sub.x, I.sub.y, 1].sup.T=K(R.sub.0(R.sup.xR.sup.yR.sup.z[x,y,z].sup.T+T)+T.sub.0), where K is the an intrinsic matrix, R.sup.x, R.sup.y, R.sup.z are the rotation matrices and T a translation vector from an arbitrary acquisition position to the reference position (R.sub.0,T.sub.0), respectively, expressed in the local coordinates of the phantom (x,y,z), wherein from the 3D coordinates of the at least three first calibration fiducials (M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3) the intrinsic matrix K, R.sub.0 and T.sub.0 are determined once for the x-ray apparatus (3) being positioned in the reference position (R.sub.0,T.sub.0).

3. The method according to claim 1, wherein the projections (IM.sub.i.sup.1, IM.sub.i.sup.2, IM.sub.i.sup.3) of the at least three first calibration fiducials (M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3) are transformed and normalized such that the 2D image coordinates (I.sub.x, I.sub.y) of the projections (IM.sub.i.sup.1, IM.sub.i.sup.2, IM.sub.i.sup.3) of the at least three first calibration fiducials (M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3) merely depend on rotations (R.sup.x,R.sup.y) about an x-axis and a y-axis, wherein said transformation transforms a 2D coordinate system having an origin 0 located at (IM.sub.i.sup.1+IM.sub.i.sup.2)/2 and an x-axis defined along the direction origin 0.fwdarw.IM.sub.i.sup.3 to a coordinate system with its origin at (0, 0) and an x-axis along direction (1, 0), and wherein said transformation normalizes the length of the vector IM.sub.i.sup.1.fwdarw.IM.sub.i.sup.3 to 1.

4. The method according to claim 3, wherein a lookup table is generated, wherein said lookup table contains all the 2D Image coordinates (I.sub.x, I.sub.y) of the detected, transformed and normalized projections of the at least three first calibration fiducials for different combinations of rotations (R.sup.x,R.sup.y) about said x- and y-axis.

5. The method according to claim 4, wherein, for each of the detected, transformed and normalized 2D image coordinates (I.sub.x, I.sub.y) of the projections of the at least three first calibration fiducials the corresponding 3D coordinates (x, y, z) are taken from the lookup table.

6. The method according to claim 1, wherein projections of at least three further calibration fiducials are detected, wherein at least one of said further calibration fiducials is a second or third calibration fiducial being arranged in a second or third plane of the phantom (20) running parallel to the first plane, wherein the second or third fiducial(s) comprise a smaller volume than the first calibration fiducials.

7. The method according to claim 6, wherein the projection matrix P is computed using the 2D image coordinates (I.sub.x, I.sub.y) of projections of the at least three first calibration fiducials and the at least three further calibration fiducials as well as the 3D coordinates (x, y, z) of the at least three first calibration fiducials and the 3D coordinates (x, y, z) of the at least three further calibration fiducials.

8. The method, according to claim 1, wherein the method further comprises the following steps, which are performed before acquiring the at least one image of said object: fixing a fiducial device comprising fiducials with respect to an image detector of an x-ray apparatus, said x-ray apparatus comprising an x-ray generating means-opposing said image detector along a viewing direction of the x-ray apparatus, and wherein the fiducial device comprises a plate carrying said fiducials, arranging the viewing direction -in different orientations and positions with respect to the earth's magnetic field, acquiring a model building image ({I.sup.n}) for each orientation, said model building image containing projections of all fiducials, detecting positions ({P.sup.i.sub.n}) of the projections of the fiducials in the respective model building image ({I.sup.n}), computing a displacement vector D.sub.n=[D.sup.1.sub.n,x,D.sup.1.sub.n,y, . . . , D.sup.N′.sub.n,x,D.sup.N′.sub.n,y].sup.T for each acquired model building image (I.sup.n), wherein the vectors D.sup.i.sub.n=[D.sup.i.sub.n,x,D.sup.i.sub.n,y].sup.T extend from the nominal position (P.sup.i) of the respective fiducial (i) to the detected position (P.sup.i.sub.n) of the projection of said fiducial (i), and wherein a mean displacement map (M.sub.D) as the average of all displacement vectors (D.sub.n) is computed and subtracted from each displacement vector (D.sub.n) and a principal component analysis is applied to the resulting variation of the displacement vectors (D.sub.n) in order to compute eigen displacement maps ({P.sup.i.sub.D}) and eigen values ({λ.sup.i}) of said variation.

9. The method according to claim 8, wherein a part of said plate comprising some of the fiducials is removed from said fiducial device so that, when an image is acquired by means of the x-ray apparatus, the image comprises a corresponding region being free from projections of the fiducials of said part as well as a further region containing projections of remaining boundary fiducials of the fiducial device, which are arranged around said region.

10. The method according to claim 8, wherein at least one image of said object is acquired by means of the x-ray apparatus such that the at least one image contains projections of the object as well as projections of boundary fiducials.

11. The method according to claim 10, wherein the method further comprises the steps of: detecting positions of the projections of the boundary fiducials in the at least one image, computing a corresponding boundary displacement vector B.sub.n=[B.sup.1.sub.n,x,B.sup.1.sub.n,y, . . . , B.sup.N′.sub.n,x,B.sup.N′.sub.n,y].sup.T, wherein the vectors B.sup.i.sub.n=[B.sup.i.sub.n,x,B.sup.i.sub.n,y].sup.T extend from the nominal position of the respective boundary fiducial (i) to the detected position of the projection of the respective boundary fiducial (i), and fitting a linear combination M.sub.D+Σ.sub.jb.sub.jP.sup.j.sub.D, where M.sub.D is a mean displacement map and {P.sup.i.sub.D} are eigen displacement maps, to the boundary displacement vector (B.sub.n) for computing an estimate for the displacement vector (D.sub.n) of all fiducials.

12. The method according to claim 11, wherein the method further comprises the step of undistorting the at least one image using the estimated displacement vector (D.sub.n) of the at least one image.

13. A Mobile Phantom for an x-ray device, comprising: a plurality of first calibration fiducials (M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3), a carrier for carrying said first calibration fiducials (M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3), wherein the phantom comprises three first calibration fiducials (M.sub.1.sup.1, M.sub.1.sup.2, M.sub.1.sup.3) being arranged along a first line, wherein the phantom comprises three further first calibration fiducials (M.sub.3.sup.1, M.sub.3.sup.2, M.sub.3.sup.3) being arranged on a third line running parallel to the first line, and wherein the phantom comprises a further first calibration fiducial (M.sub.2.sup.2) arranged on a second line running from a first calibration fiducial (M.sub.1.sup.2) that is arranged between two first calibration fiducials (M.sub.1.sup.1, M.sub.1.sup.3) on the first line to a first calibration fiducial (M.sub.3.sup.2) that is arranged between two first calibration fiducials (M.sub.3.sup.1, M.sub.3.sup.3) on the third line, wherein the second line runs perpendicular to the first and the third line, such that also three first calibration fiducials (M.sub.2.sup.1, M.sub.2.sup.2, M.sub.2.sup.3) are arranged on the second line, wherein the seven first fiducials are all arranged in a first plane, and wherein the three ratios r.sub.i=|M.sub.i.sup.1M.sub.i.sup.2|/|M.sub.i.sup.2M.sub.i.sup.3|, i=1, 2, 3, are different from each other.

14. A Fiducial device being designed to be fixed to an image detector of an x-ray apparatus, wherein the fiducial device comprises a plate carrying fiducials, wherein said fiducials are arranged in a single extension plane along which the plate extends, wherein said fiducials are arranged on lattice points of a rectangular lattice, wherein the plate comprises a part containing some of the fiducials, which part is designed to be released from the plate leaving a recess in the plate.

15. The method according to claim 10, wherein said at least one image further contains said projections of the at least three first calibration fiducials.

16. The mobile phantom according to claim 13, wherein the three ratios r.sub.i=|M.sub.i.sup.1M.sub.i.sup.2|/|M.sub.i.sup.2M.sub.i.sup.3|, i=1, 2, 3, are different from 1, and wherein the phantom comprises second calibration fiducials arranged in a second plane running parallel to the first plane, and wherein the phantom comprises third calibration fiducials being arranged in a third plane running parallel to the first and the second plane, wherein the first calibration fiducials (M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3) have a larger diameter than the second or third calibration fiducials.

17. A method for building a statistical model for image distortion correction, wherein the method comprises the following steps: fixing a fiducial device comprising fiducials with respect to an image detector of an x-ray apparatus comprising an x-ray generating means opposing said image detector along a viewing direction of the x-ray apparatus, wherein the fiducial device comprises a plate carrying said fiducials, arranging the viewing direction in different orientations and positions with respect to the earth's magnetic field, acquiring a model building image ({I.sup.n}) for each orientation, said model building image containing projections of all fiducials, detecting positions ({P.sup.i.sub.n}) of projections of the fiducials in the respective model building image ({I.sup.n}), computing a displacement vector D.sub.n=[D.sup.1.sub.n,x,D.sup.1.sub.n,y, . . . , D.sup.N′.sub.n,x,D.sup.N′.sub.n,y].sup.T for each acquired model building image (I.sup.n), wherein the vectors D.sup.i.sub.n=[D.sup.i.sub.n,x,D.sup.i.sub.n,y].sup.T extend from the nominal position (P.sup.i) of the respective fiducial (i) to the detected position (P.sup.i.sub.n) of the projection of said fiducial (i), and wherein a mean displacement map (M.sub.D) as the average of all displacement vectors (D.sub.n) is computed and subtracted from each displacement vector (D.sub.n) and a principal component analysis is applied to the resulting variation of the displacement vectors (D.sub.n) in order to compute eigen displacement maps ({P.sup.i.sub.D}) and eigen values ({λ.sup.i}) of said variation.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) Further features and advantages of the invention shall be described by means of detailed descriptions of embodiments with reference to the Figures, wherein

(2) FIG. 1 shows the fiducial device comprising a plate from which a part that is releasably fastened to the plate can be removed;

(3) FIG. 2 shows the fiducial device according to FIG. 1 without said part;

(4) FIG. 3 shows a schematical side view of an x-ray apparatus (C-arm), wherein for building the statistical model of the image distortion pattern (before image calibration), an image for each considered orientation with respect to the axes A.sub.grav and A.sub.pitch is taken from all fiducials of the fiducial device attached to the x-ray generating means;

(5) FIG. 4 shows a statistical model of an image distortion pattern of an x-ray device, i.e., the mean displacement map (top left), and the first three eigen displacement maps (top right: 1st; bottom left: 2nd, bottom right: 3rd eigen displacement map) computed by means of PCA. During image acquisition, such a statistical model of the distortion pattern will be used together with the acquired images of the plate with recess (hole) to first estimate the image distortion and then to correct the image distortion by an image dewarping process;

(6) FIG. 5 shows an example for a statistical model based image distortion correction pipeline: detection of boundary fiducial projection from the acquired image of the plate with recess (hole) (top left), estimation the statistical model of the image distortion pattern (top right), image distortion correction (bottom left) and fiducial projection inpainting (bottom right);

(7) FIG. 6 shows a perspective view of a mobile phantom with calibration fiducials;

(8) FIG. 7 shows a schematic view of the look-up table based method for 2D-3D correspondence establishment;

(9) FIG. 8 shows two examples of phantom detection;

(10) FIG. 9 a possible flow chart of the method for building a statistical model for image distortion correction; and

(11) FIG. 10 a possible flow chart of the method for phantom detection (with prior image distortion correction).

DETAILED DESCRIPTION

(12) The present invention particularly relates to an x-ray apparatus (e.g. C-arm) calibration approach that is based on a fiducial device 10, also denoted as distortion correction plate, having fiducials 100 arranged in a single extension plane (cf. FIGS. 1 and 2) as well as a mobile phantom 20 (cf. FIG. 6) that can be placed to anywhere close to the patient anatomy being imaged.

(13) Particularly, the invention relates to calibration and tracking of digital x-ray apparatuses 3 (cf. FIG. 3) such as those used in conventional computer tomographs or cone-beam computer tomographs as well as those equipped with C-arm and O-arm type stands. According to FIG. 3 such an x-ray apparatus 3 comprises an x-ray generating means 30 for generating x-rays as well as an image detector 31 for detecting the x-rays, which opposes the x-ray generating means 30 along a viewing direction 300 of the x-ray apparatus 3, along which the x-rays impinge on the object to be imaged and the detector 31. Said detector 31 may be an image intensifier 31 that is adapted to convert the x-rays into visible light from which an image can be generated by the x-ray apparatus (e.g. via a camera like a CCD camera). The x-ray generating means 30 and the image intensifier 31 are rigidly connected by a C-shaped arm 32 that can be tilted around a pitch axis A.sub.pitch and pivoted about a vertical axis A.sub.grav as indicated in FIG. 3 such that the viewing direction can be spatially adjusted (with respect to its orientation and position in space).

(14) In particular the method according to the invention is essentially based on

(15) (1) building a statistical model of the x-ray image distortion models before the calibration and using the resultant statistical model to correct the X-ray image distortion during the calibration (cf. flow chart according to FIG. 9), and/or

(16) (2) a phantom and a phantom detection method to obtain the projection parameters of the X-ray device for the acquired image as well as the position of the X-ray device with respect to the phantom at the time when the image is acquired (cf. flow chart according to FIG. 10).

(17) It is possible that for a newly introduced digital X-ray apparatus where a flat panel detector is used, there is no need to do the step (1).

(18) In such a situation, the phantom and the phantom detection method still can be used to obtain the projection parameters of the X-ray device for the acquired image as well as the position of the X-ray device at the time with respect to the phantom when the image is acquired. By tracking the phantom with an external pose tracking system, it is then possible to transform this position from the phantom coordinate system to any other reference coordinate system, although such a transformation is not necessary for certain applications.

(19) Further, the method for building the statistical model/image distortion correction can be performed independently with respect to the phantom detection method.

(20) Now, for building the statistical model, the fiducial device 10 is used. As all the fiducials 100 are arranged in one (extension) plane, the reduction of imaging space due to the distortion correction plate 10 is negligible compared to other calibration cages used in most existing systems, leading to a minimization of interference with the patient anatomy being imaged.

(21) As can be seen from FIGS. 1 and 2 the individual fiducials 100 may be cross-shaped and arranged on the lattice points of a (fictitious) square lattice. Such a grid allows one to easily track the distortions of the image due to the earth magnetic field and eventually other influences. Preferably, the fiducial device 10 comprises 77 fiducials 100 which are embedded into a particularly circular plate 11, 12 of the device 10, which plate 11, 12 extends along said extension plane. The plate 11, 12 may comprise an outer (circular) rigid frame 11 for fixing the fiducial device 10 to an image detector (intensifier) 31 of an x-ray apparatus 3, which frame 11 holds a central plate element 12, which carries the fiducials 100. A central (square) part 14 of the plate 11, 12, namely of the central plate element 12, is designed to be detachable from the fiducial device 10, so that a number of central fiducials 100 embedded into said part can be removed from the plate 10 leaving a recess (hole) 14 in the plate 11, 12. Thus, when an image is taken from the fiducial device 10 with recess 14, only projections of the boundary fiducials 101 are present in the image.

(22) Building the statistical model is only needed to be done off-line one time for a certain x-ray apparatus (C-arm) 3. At this stage, the central part 13 of the fiducial device 10 will be mounted to have a complete fiducial device as shown in FIG. 1 and the complete fiducial device 10 will be rigidly attached to the image intensifier 31 as shown in FIG. 3 such that the viewing direction 300 is particularly normal to the extension plane of the fiducial device 10. Thus if no distortion would be present, the projections of the fiducials would form the same square lattice (grid) as the fiducials (up to a scale factor).

(23) By placing now the C-arm 32 (or viewing direction 300) in different orientations and acquiring images of the fiducial device 10, i.e. of its fiducials extending along a single extension plane, a statistical model of the C-arm image distortion pattern can be established from fiducial projections detected from all acquired images. In detail, the individual imaging positions are selected by rotating the C-arm 32 (or viewing direction 300) around said axes A.sub.grav and A.sub.pitch while keeping the roll angle of the C-arm at 0 degree as shown in FIG. 3 in order to estimate the influence of the earth magnetic field on the image distortion. The first rotation Rot.sub.Grav about the Axis A.sub.grav is around the gravity direction, i.e. rotating the whole C-arm 32 horizontally to reveal the influence of the magnetic field of the earth, and the sampling interval is preferably 20°. The second rotation is the pitch rotation of the C-arm 32 and the sampling interval is also preferably 20°. Additionally three more images are generated resulting in 383 model building images in total.

(24) Given these model building images {I.sup.n}; n=1, 2, . . . , N=383, a cross-correlation based template matching is used to automatically detect the projections of all the fiducials from each image In and the detected fiducial positions {P.sub.i.sup.n}, i=1, 2, . . . , 77, n=1, 2, . . . , 383, are saved for the model building. The statistical model is then built as follows:

(25) For each model building image I.sup.n, the displacement vector D.sub.n=[D.sup.1.sub.n,x,D.sup.1.sub.n,y, . . . , D.sup.77.sub.n,x,D.sup.77.sub.n,y].sup.T is computed, where D.sup.i.sub.n=[D.sup.i.sub.n,x,D.sup.i.sub.n,y].sup.T are defined as the vectors starting from the ideal/nominal fiducial projections P.sup.i, which are determined from the physical geometry of the fiducial layout 100, to the detected fiducial projection P.sup.i.sub.n.

(26) Taking the displacement vectors {D.sub.n}, n=1, 2, . . . , 383 as the input, the mean displacement map M.sub.D is at first automatically computed as the average of all displacement vectors. The mean displacement map is then automatically subtracted from each displacement vector before a principal component analysis.sup.10 is performed on the variation of the displacement vectors to compute the eigen displacement maps {P.sup.i.sub.D} and the eigen values {λ.sup.i} of the variation. Then, any displacement map can be parameterized by the statistical model with much less number of parameters.

(27) One can show that over 98% of the variations can be explained by the first 3 eigen displacement maps P.sup.1.sub.D, P.sup.2.sub.D, P.sup.3.sub.D of the model as depicted in FIG. 4. The top left panel of FIG. 4 shows the mean displacement map M.sub.D.

(28) Intra-operatively, when image distortion correction is needed, the central part 13 of the image distortion correction plate (fiducial device) 10 will be removed before it is attached to the image intensifier 31. Thus, only the boundary fiducials 101 will be visible in an intra-operatively acquired image. The projections of these boundary fiducials 101 can be robustly automatically detected by using the same cross-correlation based template matching as mentioned before such that the displacement vectors of these fiducials 101 with respect to the corresponding nominal fiducial projections can be precisely estimated. The estimated displacement vectors of these boundary fiducials are then used to instantiate the statistical model of the x-ray apparatus (C-arm) image distortion pattern, leading to the extrapolation of the displacement vectors for those missing fiducials at the central part 13 of the plate 11, 12 of the fiducial device 10. Preferably, the extrapolation is done based on the statistical instantiation method proposed by Rajamani et al..sup.11, which describes a Mahalanobis distance weighted least square fit of a linear combination of the above described eigen displacement maps to the displacement vector of the boundary fiducials.

(29) Thus, the method according to the invention advantageously allows one to use the sparse fiducials around the border of the plate 11, 12 that are visible in each intra-operatively acquired image for estimating the projection locations of the missing fiducials by combining the detected boundary fiducial projections with the statistical model of the image distortion pattern.

(30) Based on the resulting displacement vector of all fiducials, the image is then undistorted using a fifth order polynomial-based approach. Moreover, the fiducial projections are preferably inpainted from the image with the method implemented in OpenCV library.sup.12. An example for the above-described statistical model based image distortion correction is shown in FIG. 5.

(31) In order to calibrate each intra-operatively acquired image, one needs to compute both the intrinsic and the extrinsic parameters of the image. This is particularly achieved in the present approach by means a mobile phantom as shown in FIG. 6

(32) There are totally 16 sphere-shaped calibration fiducials embedded in this phantom 20, namely 7 big first calibration fiducials {M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3}, i=1, 2, 3, having a diameter of preferably 5.0 mm and 9 small second and third calibration fiducials 201, 202 having a diameter of 2.5 mm. The 16 calibration fiducials {M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3}, 201, 202 are arranged in three different planes: all 7 big first calibration fiducials {M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3} are placed in a first plane and the remaining 9 small second and third calibration fiducials 201, 202 are distributed in a parallel second and a parallel third plane, wherein each of the three planes is different from the other planes.

(33) Furthermore, the 7 big calibration fiducials {M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3} are arranged to form three line patterns as shown in FIG. 6. Every line pattern consists of three first calibration fiducials {M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3}, i=1, 2, 3 with different ratios {r.sub.i=|M.sub.i.sup.1M.sub.i.sup.2|/|M.sub.i.sup.2M.sub.i.sup.3|}. The exact ratio for each line 21, 22, 23 is used below to identify which line pattern has been successfully detected.

(34) Particularly, the phantom 20 comprises three first calibration fiducials M.sub.1.sup.1, M.sub.1.sup.2, M.sub.1.sup.3 being arranged along a first line 21, wherein particularly the phantom comprises three further first calibration fiducials M.sub.3.sup.1, M.sub.3.sup.2, M.sub.3.sup.3 being arranged on a further (third) line 23 running parallel to the first line 21, and wherein particularly the phantom 20 comprises a further first calibration fiducial M.sub.2.sup.2 arranged on a further (second) line 22 running from a first calibration fiducial M.sub.1.sup.2 that is arranged between the two first calibration fiducials M.sub.1.sup.1, M.sub.1.sup.3 on the first line 21 to a first calibration fiducial M.sub.3.sup.2 that is arranged between two first calibration fiducials M.sub.3.sup.1, M.sub.3.sup.3 on the third line 23, wherein particularly the second line 22 runs perpendicular to the first 21 and the third line 23 such that also three first calibration fiducials M.sub.2.sup.1, M.sub.2.sup.2, M.sub.2.sup.3 are arranged on the second line 22.

(35) Preferably, the phantom 20 comprises a carrier 24 that is formed out of a (e.g. clear) plastic material and preferably comprises a plate 25 extending along the first plane, into which the first calibration fiducials are preferably embedded. Said plate 25 comprises two parallel boundary regions 26, wherein a side wall 27 protrudes from each of the boundary regions 26 via which the phantom 20 can rest on the target anatomy or the image detector (e.g. in the reference position) 31 such that the plate (first plane) 25 runs across the viewing direction 300 of the x-ray apparatus 3 (z-axis of the local coordinate system of the phantom 20, cf. FIG. 7).

(36) Each side wall 27 comprises a first step 28 into which the second calibration fiducials 201 are embedded, which first steps 28 extend along said second plane, as well as a second step 29 into which the third calibration fiducials 202 are embedded, which second steps 29 extend along said third plane.

(37) After the distortion of an intra-operatively acquired image is corrected, a sequence of automatic image processing operations may be applied to the image. As those calibration fiducials {M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3}, 201, 202 are made from a metal (steel), a simple threshold-based method is preferably first used to segment the image. Connected-component labeling is then applied to the binary image to extract a set of separated regions. Morphology analysis is further preferably applied to each label connected-component to extract two types of regions: candidate regions from (big) first calibration fiducial projections and candidate regions from small calibration fiducial projections. The centers of these candidate regions are regarded as projections of the center of a potential calibration fiducial. Due to background clutter, it is feasible that some of the candidate projections are outliers and that one may miss some of the true fiducial projections. Furthermore, to calculate both the intrinsic and the extrinsic parameters, we have to detect the phantom 20 from the image. Here, phantom detection means to establish the correspondences between the detected 2D fiducial projection centers and their associated 3D coordinates x,y,z in the local coordinate system of the phantom 20.

(38) For this purpose, a robust simulation-based approach is utilized. The pre-condition to use this method to build the correspondences is that one of the three line patterns of first calibration fiducials M.sub.i.sup.1, M.sub.i.sup.2, M.sub.i.sup.3, i=1, 2 or 3, has been successfully detected. Due to the fact that these line patterns are defined by (big) first calibration fiducials, the chance to missing all three line patterns is rare.

(39) We model the C-arm projection using a pin-hole camera.
α[I.sub.x, I.sub.y, 1].sup.T=K(R[x, y, z].sup.T+T)=P[x, y, z, 1].sup.T  (1)

(40) where α is a scale factor, K is the intrinsic calibration matrix, R and T are the extrinsic rotation matrix and translational vector, respectively. Both the intrinsic and the extrinsic projection parameters K, R and T can be combined into a 3-by-4 projection matrix P in the local coordinate system established on the mobile phantom.

(41) The idea behind the simulation-based method is to do a pre-calibration to compute both the intrinsic matrix K as well as the extrinsic parameters R.sub.0 and T.sub.0 of the C-arm in a reference position from Eq. (2) as shown in FIG. 7. According thereto, the viewing direction 300 of the x-ray apparatus 3 aligns with the z-axis of the local coordinate system of the phantom 20, which is fixed via its side walls 27 to the fiducial device 10 such that the first plane of the phantom 20 (i.e. the plate 25 of the carrier 24) extends perpendicular to the viewing direction 300 and to the plate 11, 12 of the fiducial device 10. Thus, in the reference position the plane spanned by the x- and y-axis of the local coordinate system of the phantom extends parallel to the plate 11, 12. The origin of the local coordinate system is preferably arranged in the center of the plate 25 of the carrier 24 of the phantom 20.

(42) Then, assuming that the intrinsic matrix K is not changed from one image to another (we only use this assumption for building the correspondences), the projection of the x-ray apparatus 3 (C-arm) at any other position with respect to the phantom 20 can be expressed as
α[I.sub.x, I.sub.y, 1].sup.T=K(R.sub.0(R.sup.xR.sup.yR.sup.z[x,y,z].sup.T+T)+T.sub.0)  (2)

(43) where R.sup.x, R.sup.y, R.sup.z and T are the rotation matrices around three axes (assuming the z-axis is in parallel with the view direction of the C-arm at the reference position, see FIG. 7) and the translation vector from an arbitrary acquisition position to the reference position, respectively, expressed in the local coordinate x, y, z of the mobile phantom 20. To detect the phantom projection when the x-ray apparatus (C-arm) 3 is placed in a new position (cf. FIG. 7), the simulation-based approach consists of two steps.

(44) Firstly, one wants to get rid of the influence of the parameters R.sup.z, α, and T on the phantom detection by normalizing the image acquired at the new position (arbitrary acquisition position in FIG. 7) as follows. Assuming that one knows the correspondences of fiducials on one line pattern, which is defined by 3 landmarks M.sup.1, M.sup.2, M.sup.3 with their correspondent projections at IM.sup.1, IM.sup.2, IM.sup.3, one can define a 2D coordinate system based on IM.sup.1, IM.sup.2, IM.sup.3, whose origin 0 is located at (IM.sup.1+IM.sup.2)/2 and the x-axis is defined along the direction origin 0.fwdarw.IM.sup.3. Accordingly a 2D affine transformation T.sub.normalize can be computed to automatically transform this line pattern based coordinate system to a standard 2D coordinate system with its origin at (0, 0) and x-axis along direction (1, 0) and at the same time to normalize the length of the vector IM.sup.1.fwdarw.IM.sup.3 to 1. By applying T.sub.normalize to all the fiducial projections, it can be observed that for a pair of fixed R.sup.x and R.sup.y, we can get the same normalized image no matter how the other parameters R.sup.z, α, and T are changed because the influence of these parameters is just to translate, rotate, and scale the fiducial projections, which can be compensated by the normalization operation. Therefore, the fiducial projections after the normalization will only depend on the rotational matrices R.sup.x and R.sup.y.

(45) Since the distribution of the fiducial projections in the normalized image only depends on the rotation matrices R.sup.x and R.sup.y, it is natural to build a look-up table which up to a certain precision (e.g., 1°) contains all the normalized fiducial projections with different combination of R.sup.x and R.sup.y. This is done off-line by simulating the projection operation using Eq. (2) based on the pre-calibrated projection model of the x-ray apparatus (C-arm) 3 at the reference position R.sub.0, T.sub.0 as shown in FIG. 7. For an image acquired at position other than the reference, said normalization operation as described above is automatically applied to all the detected candidate fiducial projections. The normalized candidate fiducial projections are then compared to those in the look-up table to find the best match. Since the items in the look-up table are generated by a simulation procedure, we know exactly the correspondence between the 2D fiducial projections and their corresponding 3D coordinates. Therefore, we can establish the correspondences between the candidate fiducial projections and the fiducials embedded in the phantom. An example of a phantom detection is shown in FIG. 8.

(46) As soon as the correspondences are established, we can further fine-tune the fiducial projection location by applying a cross-correlation based template matching. After that, preferably a direct linear transformation algorithm as described in detail in.sup.13 is used to compute the projection matrix P, i.e. the desired intrinsic and extrinsic parameters.

REFERENCES

(47) 1. Yaniv, Z., et al.: Fluroscopic image processing for computer-aided orthopaedic surgery. MICCAI'98, LNCS 1496 (1998) 325-334. 2. Hofstetter, R., et al.: Fluoroscopy as an imaging means for computer assisted surgical navigation. Computer Aided Surgery 4 (1999) 65-76. 3. Livyatan, H., et al.: Robust automatic C-arm calibration for fluoroscopy-based navigation: a practical approach. MICCAI'02, LNCS 2489 (2002) 60-68. 4. Kainz, B., et al.: Fast marker based C-arm pose estimation. MICCAI'08, LNCS 5242 (2008) 652-659. 5. Chintalapani, G., Jain A. K., Taylor, R. H.: Statistical characterization of C-arm distortion with application to intra-operative distortion correction. SPIE Medical Imaging 2007, Vol. 6509 (2007) 65092Y.1-65092Y.8. 6. Chintalapani, G., Taylor, R. H.: C-arm distortion correction using patient CT as fiducial. ISBI 2007, (2007) 1180-1183. 7. Yao, J., et al.: A C-arm fluoroscopy-guided progressive cut refinement strategy using a surgical robot. Computer Aided Surgery 5 (2000) 373-390. 8. Jain, A. and Fichtinger, G.: C-arm tracking and reconstruction without an external tracker. MICCAI'06, LNCS 4190 (2006) 494-502. 9. Otake, Y., et al.: An iterative framework for improving the accuracy of intraoperative intensity-based 2D/3D registration for image-guided orthopedic surgery. IPCAI'10, LNCS 6135 (2010) 23-33. 10. Jollie, I. T.: Principal component analysis (Springer series in statistics). 2nd Edition, 2002, Springer-Verlag, New York. 11. Rajamani, K. T., et al.: Statistical deformable bone models for robust 3D surface extrapolation from sparse data. Medical Image Analysis, 11 (2007) 99-109. 12. OpenCV Library:http://sourceforge.net/projects/opencvlibrary/files/opencvwin/2.2/ 13. Hartley, R. and Zisserman, A.: Multiple view geometry in computer vision. 2.sup.nd Edition, 2004, Cambridge Univeristy Press.