

ORIGINAL ARTICLE 



Year : 2019  Volume
: 44
 Issue : 3  Page : 191200 

A novel design of Nfiducial phantom for automatic ultrasound calibration
Maria Chatrasingh, Jackrit Suthakorn
Department of Biomedical Engineering, Center for Biomedical and Robotics Technology (BART LAB), Faculty of Engineering, Mahidol University, Salaya, Thailand
Date of Submission  11Sep2018 
Date of Decision  11Jun2019 
Date of Acceptance  16Jun2019 
Date of Web Publication  13Sep2019 
Correspondence Address: Dr. Jackrit Suthakorn Department of Biomedical Engineering, Center for Biomedical and Robotics Technology (BART LAB), Faculty of Engineering, Mahidol University, Salaya Thailand
Source of Support: None, Conflict of Interest: None  Check 
DOI: 10.4103/jmp.JMP_92_18
Abstract   
Background: Freehand ultrasound (US) is a technique used to acquire threedimensional (3D) US images using a tracked 2D US probe. Calibrating the probe with a proper calibration phantom improves the precision of the technique and allows several applications in computerassisted surgery. Nfiducial phantom is widely used due to the robustness of precise fabrication and convenience of use. In principle, the design supports singleframe calibration by providing at least three noncollinear points in 3D space at once. Due to this requirement, most designs contain multiple Nfiducials in unpatterned and noncollinear arrangements. The unpatterned multiple Nfiducials appearing as scattered dots in the US image are difficult to extract, and the extracted data are usually contaminated with noise. In practice, the extraction mostly relied on manual interventions, and calibration with Nfiducial phantom has not yet achieved high accuracy with single or few frame calibrations due to noise contamination. Aims: In this article, we propose a novel design of the Nfiducial US calibration phantom to enable automatic feature extraction with comparable accuracy to multiple frame calibration. Materials and Methods: Along with the design, the Random Sample Consensus (RANSAC) algorithm was used for feature extraction with both 2D and 3D models estimation. The RANSAC feature extraction algorithm was equipped with a closedform calibration method to achieve automatic calibration. Results: The accuracy, precision, and shape reconstruction errors of the calibration acquired from the experiment were significantly matched with the previous literature reports. Conclusions: The results showed that our proposed method has a high efficiency to perform automatic feature extraction compared to conventional extraction performed by humans.
Keywords: Freehand ultrasound, Nfiducial ultrasound calibration phantom, Random Sample Consensus
How to cite this article: Chatrasingh M, Suthakorn J. A novel design of Nfiducial phantom for automatic ultrasound calibration. J Med Phys 2019;44:191200 
How to cite this URL: Chatrasingh M, Suthakorn J. A novel design of Nfiducial phantom for automatic ultrasound calibration. J Med Phys [serial online] 2019 [cited 2019 Oct 18];44:191200. Available from: http://www.jmp.org.in/text.asp?2019/44/3/191/266852 
Introduction   
The threedimensional (3D) freehand ultrasound (US) technique involves the determination of positions and orientations of realtime US images in a 3D space. This technique transforms a conventional US into an advanced US that could determine the position of the scanned objects. In medical application, sweeping the 3D freehand US over the anatomy reconstructs the 3D volume of the anatomy.^{[1],[2],[3]} This technique provides an advantage over the conventional 3D US for the unlimited range of reconstruction. Furthermore, the conventional USguided intraoperative surgery, which surgeons perform surgery based on US image visualizing, could be enhanced by determining the realtime location of a cancer lesion or an organ of interest with respect to surgical tools.^{[4],[5],[6]}
If one wants the US to provide positioning a scanned object accurately, an US calibration process must be performed. The realtime positions and orientations of the US probe in 3D are determined by a tracking device. A tracking sensor is attached on the outer body of the US probe so that the tracking device could read the position of the sensor moving with the probe. The conventional US probes are attached with a tracking sensor so that the motion of the probe is known in real time. The US calibration process finds the rigid transformation relating the attached sensor to the US image. This rigid transformation is varied depending on the US probe of choice and the location where the sensor attached on the probe. The US calibration involves matching the US images with the known corresponding position in space. The calibrating structure with known geometry is called US calibration phantom. [Figure 1] shows the relationship of spatial transformations in the US calibration process which can be expressed as:  Figure 1: Frame of motion and the relation between each frame in ultrasound calibration
Click here to view 
where X_{Ph} = [x, y, z,1]^{T} and X_{Im} = [u, v,0,1]^{T} are the 4 × 1 position vectors of the corresponding point sets in the phantom frame and image frame, respectively. is the 4 × 4 matrix combining the rotation and translation of frame A which refers to frame B. is the transformation to determine the position and orientation of the tracking device which refers to the calibrated phantom estimated using a spatial localizer or a tracking sensor. is the transformation between the tracking sensor on the US probe and the tracking device. is the transformation between US images and the sensor on the probe. S is a 4 × 4 scale matrix.
During the calibration process, the tracked probe scans the US calibration phantom to capture the feature of the phantom display in US images, while the tracking device reads the corresponding transformation matrices: and The algorithm behind US calibration involves (1) extraction the point set displayed in US images: and matching to the corresponding phantom geometry: and (2) estimation of and S according to Equation 1.
Most stateoftheart US calibrations focus on the accuracy of the technique regardless of the contribution of other factors to its efficiency. Practical US calibration requires the ability to perform intraoperatively. Chen et al.^{[7]} addressed that the sterilization protocols might require US recalibration inside the operation room. During surgery, the surgeon might need to validate and recalibrate the accuracy of the freehand US. Therefore, the practical US calibration needs to be fast, automatic, and accurate, for the use intraoperatively. As per our literature review, we barely found any study in which high accuracy and fast and automatic functionality were simultaneously achieved in the US calibration. The choice of the US calibration phantom used in the calibration procedure is thought to be one of the pivotal factors responsible for achieving efficiency in the process.
Several US calibration phantom designs have been proposed in the last two decades each with particular benefits.^{[8]} The point phantom, crosswire phantom, and planetype phantom are known to provide insufficient dimensionality to state a corresponding 3D position of a phantom in a single scan. As a result, these phantoms require a number of input US images and the corresponding tracking data to iteratively define and S. On the other hand, the Nfiducial phantom (sometimes called Zfiducial phantom)^{[7],[9],[10],[11],[12],[13],[14]} is one of the few phantom designs that provide full X_{Ph} at a single scan. It used the construction of three wires oriented in a plane to form an N shape or Z shape. The intersection of the Nfiducial (three wires) with the scan plane determines one 3D position where the scan plane intersects the fiducial (mathematical explanations are in the next section). Whenever at least three noncollinear Nfiducials were intersected within a scan plane, the position and orientation of that plane in space could be determined. In other words, at least three noncollinear X_{Ph} in a single frame was considered sufficient to solve the least square problem in Equation 1 for the frame transformations, and S. Compared to other designs of the phantom that provide singleframe capability such as 2D alignment phantom,^{[15]} the Nfiducial phantom is more widely used due to the robustness of precise fabrication and convenience of use.^{[8]}
In theory, calibration using Nfiducial phantoms could be fast, demanding only oneframe input image and using a (noniterative) closedform solution ^{[16],[17],[18]} to solve for and S. However, the main drawback of using this phantom is a timeconsuming process of extracting the feature of the phantom from the US images.^{[14]}
Since at least three noncollinear points are required to determine a plane, multiple Nfiducials are placed in an unpatterned (noncollinear) multitude of designs. Hsu et al. pointed out that the main drawback of using this phantom is the difficulty to segment pixels that represent wireplane intersection in the US images.^{[14]} First of all, the unpatterned multiple Nfiducials appear as scattered dots in a US image which is difficult to be distinguished from the background speckle noise, signal artifact, or floating debris in the US medium. Second, an individual wire does not appear as a single dot or disc in the US scans as expected. Ameri et al.^{[19]} have demonstrated that a US scan of a rope or a rodlike component (as in the Nfiducial phantom) would produce various patterns of US hyperechoic noise including widespread or blurred. This not only compounds the segmentation difficulty but also introduces the inaccuracy to identify them.
Several works have expressed the need for performing the extraction manually.^{[9],[10],[11],[20]} The manual extraction process is tedious and timeconsuming, making the technique unlikely to be used in the operating room. According to Ameri et al.,^{[19]} the effect of US hyperechoic noise could result in local bias in the signal such that distorted displayed features. In such a case, manual extraction is considered an inaccurate extraction method since it relies on the intensity display in the image.
Lindseth et al.^{[11]} and Hsu et al.^{[14]} have successfully automatized the segmentation of images from an Nfiducial phantom given a predefined searching region. Their approaches do require manual identification of either the searching window or determination of the default scaling parameters. However, the accuracy of the extracted data has not been evaluated qualitatively. The concern on the effect of signal distortion has not been accounted or discussed.
The sophisticated approaches for automatic extraction include modifying the phantom so that the scan images display some linear patterns that could be picked up by the algorithmic image processing. Several attempts have been made in the recent past to align multiple Nfiducials in an obvious pattern. Comeau et al.^{[9]} proposed a design where three Nfiducials form a floor and two parallel side walls of a rectangular box that, in the scan, appear as dots forming corners of a rectangle. Although the design could easily be segmented, it lacks a repetitive pattern which enables one to recognize and correct the distortion of an individual dot. Plus toolkit ^{[21]} suggested the Nfiducial phantom design where the equivalentsized fiducials were flipped left to right layer by layer. ŚwiaątekNajwer et al.^{[22]} and Rao et al.^{[23]} also calibrated a US probe on a phantom with a similar concept to Plus toolkit. The lateral strands of all fiducials were vertically coaligned such that two parallel and straight dotted lines visible in scan images served as the landmark for detection. Nevertheless, these patterns constraint only on the lateral strands but not on the middle strands. The distortion of the middle strands could remain unrecognized.
A different approach includes constraining the extracted positions in three dimensions. Knowing the collinearity in 2D space constraints is insufficient. Chen et al.^{[7]} added the fiducial registration error (FRE) of the estimated Nfiducial in 3D space as feedback to iteratively exclude unwanted data in the next stage. The algorithm extracts US images in two phases: faster in the first extraction phase and then slowed down in the FRE iterative feedback phase (with an average of 12.5 s per calibration trial). Even though Chen's extraction has the 3D constraint, the iterative feedback method requires a lot of execution time and a large amount of input data without guarantee of convergence. This is in contrast to the benefits of Nfiducial designs which are capable of singleframe calibration.
Other than the automatic extraction, the step of computing the calibration results, and S, is accounted for fast and accurate calibration. Most research on US calibration used the iterative methods to minimize FRE of the calibration results. The iterative methods could determine the calibration in high accuracy with noisecontaminated data. The drawback of the methods is that it requires a large number of US image extracted data, uses a relatively high amount of time, and generally could not guarantee of the true minimization. The different method is using a closedform solution to compute the calibration results. Pagoulatos et al.^{[10]} used the closedform solution based on singular value decomposition (SVD)^{[16]} to calculate the calibration results and achieved lower accuracy compared to the iterative method. Boctor et al.^{[17]} used the closedform solution based on robot hand–eye calibration to determine the transformation between two moves of the probe using only two corresponding frames. Boctor et al. only presented the preliminary evaluation and effectiveness. The calibration involving the Nfiducial phantom with a closedform solution has yet to achieve high accuracy with few US frame. The low accuracy could result from the limitation of the closedform solutions which inaccuracy occurs with noisecontaminated input data.
As per our literature review, we barely found any study in which high accuracy and fast and automatic functionality were simultaneously achieved in the US calibration. The aforesaid automatic extraction method ^{[7],[11],[14]} relies on semiautomatic approaches and iterative methods for calibration to achieve highaccuracy results. The calibration using the closedform solutions was never performed by data taken from automatic extraction and evaluated with clear verification.
The objective of our research was to develop an automatic US calibration process by proposing a phantom design and the corresponding automatic extraction algorithm. The phantom and the algorithm should extract precise data so that the closedform solution ^{[16]} can accurately estimate the transformation matrix. The proposed phantom design is based on Nfiducial phantom. The phantom provides recognizable features on a US image while still retaining the noncollinearity of multiple Nfiducials. The automatic extraction algorithm corresponding to this design was proposed based on the Random Sample Consensus (RANSAC)^{[24]} model estimation. The constraints in both 2D and 3D positions are considered in the model. Thus, it should better extract the precise positions than the mentioned algorithms ^{[7],[21]} which only 2D constraints are used. The efficiencies of the proposed phantom design and the automatic extraction were evaluated based on the accuracy of the extraction itself and the accuracy of the overall calibration using the same closedform method.
Materials and Methods   
Phantom design
The Nfiducial phantom used the geometrical relation of at least three noncollinear Nfiducials to position the US scanning plane, using the intersection of the plane with the fiducial. For one Nfiducial, the US scan plane intersects the fiducial at three points, P _{r},P _{m}, and P _{l}, as illustrated in [Figure 2]a and [Figure 2]b.  Figure 2: The Nfiducial in the phantom. (a) Top view of an Nfiducial phantom that contains only one Nfiducial and (b) the relation of an Nfiducial intersects with the ultrasound scan plane
Click here to view 
The 2D vectors, i.e., P _{r},P _{m}, and P _{l}, represented the corresponding points in the US image frame reference. The position vector of the point where the plane intersected with the middle strand, i.e., was determined by the distance ratio between P _{r},P _{m}, and P _{l}:
where W and H were the width and height, respectively, of an Nshape fiducial. x and y were the vectors of directions along the width and height, respectively. z was the vector with the vertical direction. r was the scalar ratio of the distance between P _{r} and P _{m} and the distance between P _{r} and P _{l}, i.e.,
To attain singleframe calibration, at least three noncollinear P _{m} needs to be identified in the same US scan to determine a transformation of 6 degrees of freedom. The n numbers of fiducials scanned in the same US image frame were considered; the point sets {P _{m, j}}; j = 1,2…, n were the intersection points of the middle strand of the j^{th} Nfiducials with a scan plane. Note that the 3D position vectors of these point sets were linearly independent (noncollinear) if the matrix formed with these vectors as row vectors had at least a rank of 3.
Most studies varied x_{a, j} or z_{a, j} arbitrarily among the set of points and kept H_{j} and W_{j} constants so that the point set {P _{m, j}} was noncollinear [Figure 3]a. As a consequence, none of the sets of {P _{l, j}}, {P _{m, j}}, and {P _{r, j}} were seen as collinear in the scanning plane and made extraction difficult [Figure 3]b. In our proposed design [Figure 3]c, H_{j} was arbitrary while x_{a, j} was kept linearly dependent. Since x_{a, j} indicated the position of the left strand and H_{j} was constant, the left and right strands of multiple Nfiducials formed two parallel planes. Thus, {P _{l, j}} and {P _{r, j}} formed two dotted lines in the scanning plane. The lines were easily recognizable and established the lateral boundary of {P _{m, j}}, as depicted in [Figure 3]d.  Figure 3: Comparison between the convention Nfiducial phantom (a and b) design and our design with N fiducial arrangement (c and d)
Click here to view 
Considering an individual Nfiducial's height along with the constant width (or linear change of width), all the Nfiducials are arranged in a linear pattern avoiding the coalignment of all middle wires. For the prototype design in this article, the heights of Nfiducials were varied while the widths were constant.
[Figure 4] shows the prototype made of two lasercut acrylic plates for the front and back sides (90 mm × 70 mm × 5 mm) and two side plates (100 mm × 70 mm × 5 mm). The assembly of four plates formed an openended box of dimensions (110 mm × 90 mm × 70 mm). The front and back plates had 15 holes with diameter 1 mm, laser cutting corresponding to 5 Nfiducials stacked in the vertical direction. The front plate had an extension part for attaching a tracking sensor. It ensured that the origins of every fiducial (where its left strand and middle strand intersect) coincided in the x–y plane [Figure 2] and were placed outside of the box.  Figure 4: The prototype of our phantom design (a) Top view and (b) closeup view showing the arrangement of 5 Nfiducials
Click here to view 
Automatic extraction algorithm
RANSAC ^{[24]} is a method used to estimate the parameters of a certain model from a set of data contaminated by a large amount of outlier data. In the input US images, the extracted dots from the image processing technique may contain the outlier dots generated by an US echo shadow or floating debris up to 20% of the data. RANSAC was used to exclude the outliers and estimate the parameters of two models: one is the 2D model entitled by the pattern of Nfiducial alignment and another is the 3D plane, knowing that the scanning features are the results of the scanning plane intersection with the phantom.
First, thresholding and morphological operations were applied to US images. Then, RANSAC was used to estimate the parameters of the models. [Figure 5]a depicts the feature point sets (displayed as red dots) extracted from an example input US image where RANSAC was used to estimate five parameters of the model of a palette of evenly spaced rectangles: α, δ_{1}, δ_{2}, δ_{3}, and δ_{4}. The equations to define the model were based on the basic equation of 2D lines as follow:  Figure 5: Two estimated models for the Random Sample Consensus algorithm (a) the twodimensional model where red dots represent the position of the data extracted via initial image processing technique and (b) the threedimensional plane model where red dots represent threedimensional positions from N shape geometry
Click here to view 
two parallel lateral lines:
and
five horizontal lines:
where n is the number of Nfiducials stacked vertically in the design. α, δ_{1}, δ_{2}, δ_{3}, and δ_{4} are the parameters in the equation of lines.
With each estimation of the 2D model, the point sets were then categorized into three groups (the left strands, the right strands, and the middle strands) and assigned to each orthogonal line, such that P _{l, j}, P _{m, j}, and P _{r, j} are points in the US image representing the left, middle, and right strands of j^{th} Nfiducial, respectively. The assignment of the point set was defined using distance matrices.
For the next stage, the left and right strand 2D points, P _{l, j} and P _{r, j}, were replaced by the points, and , that defines the intersection in the palette (intersection between Equation 7 and Equations 56). and could be thought of as the projection of the points onto the estimated model.
and the point representing the middle strand, , was replaced by , its orthogonal projection on the corresponding j^{th} horizontal line (Equation 7), as follows:
where Then, each set of and was used to estimate a 3D point set using Equations 24, as shown in red dots in [Figure 5]b.
RANSAC was again used to estimate three parameters representing the scanning plane in 3D space from P _{m} point set,
Estimated plane equation:
and the 3D point representing the middle strand, , was replaced by , the orthogonal projection of P _{m, j} on the estimated plane where
The estimations of both models are continued in an iterative fashion where the summation of the estimated error is calculated to verify if both the estimated models are true.
On reaching the termination criterion of RANSAC, X_{Ph} and X_{Im} were defined as and . X_{Ph} and X_{Im} were substituted into Equation 1 such that and S can be determined via the calibration method.
Calibration method
To calculate the 4 × 4 matrix , the closedform solution for mapping two 3D point sets using the SVD given by Arun et al.^{[16]} was used. The twopoint sets were X_{Ph} and X_{Im} acquired from the last section. The closedform solution is described using the SVD of the point sets' covariance matrix to determine the 3 × 3 rotation matrix and the 3 × 1 position vector between the two sets. Although most of the closedform solutions ^{[16],[17],[18]} have a drawback in determining nonuniform scaling on US images, several US calibration works overlooked the issue and acquired decent efficiency.^{[10],[11],[25],[26],[27]}
The rotation, translation, and scale factors were derived as follows. Given K the covariance matrix of and which are normalized and about their centroids, , SVD of K was calculated by the following formula:
The rotation matrix was determined by
where
The translation between the two sets was
where and are the centroids of and , respectively.
Acquiring and gives according to the following equation:
was identified by solving the equation,
The scalar factor λ where S = λI is estimated according to Dosse and Ten Berge ^{[28]} as:
Acquisition system
The acquisition system consisted of a US scanner (iU22 xMatrix Ultrasound system, L125/50mm Linear Probe), an optical tracking system (Polaris Vicra ^{®}, Northern Digital Inc., Canada), a 2.40 GHz, Intel(R) Core(TM) i5 CPU computer with MATLAB 2017 (The MathWorks, Inc., Natick, Massachusetts, United States) image acquisition through a frame grabber, and our developed US calibration phantom. The optical tracking system was composed of a tracking camera, which was the world coordinate reference throughout the study. One tracking sensor was also attached to the phantom, and another sensor was attached to the US probe. The overall scenario of the experiment is shown in [Figure 6]a, [Figure 6]b, [Figure 6]c.  Figure 6: Experiment setup (a) experimental scenario, (b) probe's sensor attachment and breast phantom, (c) ultrasound calibration phantom, and (d) A series of ultrasound image displaying a target lesions
Click here to view 
A linear probe at 40mm scanning depth was calibrated and evaluated in this study. During the calibration, the calibration phantom and two evaluation phantoms were filled with cooking gelatin and distilled water at a ratio of 1:5 (w/v). The ratio was previously calibrated for isotropic scaling in the US image acquired from the probe. With the default US machine correctness of resolution, the gelatin solution should provide a sound speed of approximately 1540 m/s.
Evaluations
Simulation result
Fifty scanning planes were simulated giving the corresponding 3D coplanar point set, where the Nfiducials intersect with the plane. The 3D coplanar points were then transformed into the 2D image frame giving a random homogeneous scale: λ ∈ [0.5, 1]. Note that we synthesized the position data and inputted them to RANSAC extraction directly. The isotropic Gaussian noise with zero mean and variance (σ ∈ [0, 10]) was added to distort the 2D simulated data. Twenty percent of outliers were randomly added to contaminate the data. For each noise variance, the 2D data were extracted with the automatic RANSAC algorithm, and the calibration matrixes with the SVD solutions were calculated. The following were evaluated:
 The reduction in noise variance after extraction
 Errors in rotation estimation using extracted data. The deviation from the identity matrix ^{[29]} to measure the relative error between the estimated rotation matrix and the synthetic one was evaluated. Let be the synthetic rotation matrix while R being the estimated rotation matrix. . is the relative error between the true rotation and the estimated rotation. denotes the Frobenius norm of the matrix.
 Errors in translation estimation using extracted data
 Errors in scale estimation using extracted data
(2) to (4) were compared to the results using the data without model estimation.
Real data
The automatic extraction algorithm is evaluated by its ability of the algorithm to extract the features in US images accurately. Since the true locations of the feature were unknown, the current extraction ground was the extraction by humans. Five operators extracted the same set of 10 random US images. The extraction procedure included identifying the location of the dots with reference to the fiducials. The accuracy of the automatic extraction algorithms was validated using RMS errors between the locations extracted by the algorithm and by humans. On the other hand, the precision was measured using the standard deviation of the extracted features and compared between the algorithm and the manual method.
The data extracted by the proposed algorithm were used to determine the calibration results and compared with the manually extracted data. The accuracy of the calibration results was tested with a crosswire phantom where the center of the cross was the verifying position. A tracking sensor was attached to the crosswire phantom, and the tracking pointer was used to localize the related phantom geometry. The accuracy was evaluated by ensuring whether the 3D freehand US using results from calibration could localize the position where two wires crossed with respect to the tracker.
The precision of the calibration was verified using calibration reproducibility, which is the standard method to assess the precision ^{.[30]} Four corners of the US image and the center of the image were transformed into their 3D positions according to the tracker usage,
Other than point accuracy, the reconstruction accuracy was also taken into account. A halfsphere phantom made of gelatin was used as a test phantom to evaluate the US reconstruction. Five rubber balls embedded within the gelatin [Figure 7]a were used as the scan targets.  Figure 7: (a) The computeraided design of the test phantom and (b) the reconstruction of ultrasound images acquired when scanning the test phantom
Click here to view 
The balls were rigidly fixed to the plate and covered with a gelatinous medium with the ratio as described in the previous section. The operators used the freehand US to scan each knob until sufficient amounts of images were acquired. The US scans of the knob are shown in [Figure 6]d. The image and corresponding tracking data were recorded and calculated offline for reconstruction accuracy. The reconstruction accuracy for the sphere phantom was the estimate of reconstructed sphere diameter and the centroid location [Figure 7]a and [Figure 7]b.
Results   
Simulation data
The simulated positions as mentioned previously were evaluated for the process with and without the modelbased extraction algorithm at different noise levels. The ability of the extraction algorithm to reduce the noise in the data is shown in [Figure 8]a. The solid line represents the position data extracted by our RANSAC model. The noise contaminated in the data (σ input) decreases after the extraction (σ output). The dotted line represents the data without extraction which will be used as the baseline for the computation of calibration results.  Figure 8: The ability of the extraction algorithm to withstand the noised data. (a) The reduction of noise after extraction, (b) the ability to estimate the rotation matrix, (c) the ability to estimate the scale factor, and (d) the ability to estimate the translation
Click here to view 
After the extraction, the output data were used to calculate the calibration result with the SVD closedform solution. [Figure 8]b, [Figure 8]c, [Figure 8]d shows how the automatic extraction method could reduce errors in estimating the rotation, translation, and scale factors compared with the data contaminated by synthetic noise.
Experimental data
The extraction, precision, and accuracy are depicted in [Table 1], where comparison was made between the automatic extraction and the manual extraction. The manual extraction was the ground truth for accuracy; therefore, no accuracy is provided for itself. After calculation with the SVD closedform solution, the precision and accuracy of calibration results, with the amount of input in the US image, are shown in [Figure 9]a and [Figure 9]b.  Table 1: Comparison between the accuracies of automatic and manual extraction methods
Click here to view 
 Figure 9: Precision (a) and accuracy as verified by the crosswire phantom (b) of the algorithm versus the number of input scanned ultrasound images
Click here to view 
The overall precision and accuracy (root mean square error) report for the proposed automatic calibration is shown in [Table 2].  Table 2: The precision and accuracy of the proposed automatic calibration at a number of ultrasound image input (n=3)
Click here to view 
The outcome of the reconstruction had the error of diameter estimation at 2.34%, resulting in 5.48% of volume estimation. The accuracy of estimating the reconstruction location was 1.04 mm.
Discussion   
The design of the arranged Nfiducial phantom facilitates automatic extraction by providing a specific pattern for the left and right lateral boundaries. This pattern could easily recognize and exclude outliers, which might have been produced from US speckle noise, dust, bubble, etc. The proposed conceptual design changes the height of the fiducials while retaining the width, thus creating the linear patterns in a US image. The design has advantages over Comeau et al.^{[9]} since in that design only three points were used to define a line and local distortion could easily corrupt the estimation of the line. Furthermore, their designs are not space saving in forming a line from several dots and would be difficult to use in a limited fieldofview linear array probe. For Lasso et al.,^{[21]} ŚwiaątekNajwer et al.,^{[22]} Rao et al.'s ^{[23]} designs, and Plus toolkit,^{[21]} the stack of samesized Nfiducials and mirrored fiducials would cause near linearity of point sets in the center of the design where it should be the most scanned area. The data acquired from that region would be rank deficient for estimating a transformation matrix. Our design places the origin of individual Nfiducials outside the phantom and interchanges the width of the Nfiducials. Therefore, nonlinearity in the phantom is guaranteed.
Using the proposed design with automatic extraction based on the RANSAC algorithm, the pattern provides a predictable model, both as a 2D model (array of rectangles) and a 3D plane model (scanning plane estimation). After the estimation of the model's parameters from the input US image, the outlier of datasets such as those from shadow or debris was excluded and the projection of the inlier dots was used instead of using the intersection of line models. This constraint extraction process ensures that the distorted data and outlier data would have less effect on the accuracy of the calibration. This facilitates data recovery, thereby reducing the noise level to approximately 22%, as shown in the results obtained from the simulation technique.
The projection of the point set in the estimated plane was ideal for estimating the transformation matrix. Carbajal et al.^{[31]} reported similar improved outcomes using the projection data (inplane error) instead of using raw data as the input of optimization for the transformation matrix.
Chen's extraction algorithm ^{[7]} relies on the FRE as feedback to exclude the outlier. The algorithm requires multiple input frames of an image to calculate the calibration result and then group the calibration results to exclude the outlier. As a result, the algorithm would use about 12.5 s per calibration and more than 60 s to accomplish the calibration of all US depth options. Plus toolkit ^{[21]} provides an automatic extraction algorithm for the designed Nfiducial phantom. Their algorithm uses the constraint of linearity for the lateral strand feature. To conserve the noncollinearity of the scanned data, the middle strands could not form a linear pattern in a US image. The algorithm of Plus toolkit has no constraint on the middle strands and could lead to extraction error. Our extraction, on the other hand, adds the 3D constraint (plane estimation) as the guide to extract the middle strands.
The automatic extraction technique and modelbased extraction algorithm were devised using the new phantom design which introduced notable improvements in the rotational and translational components acquired from SVD closedform solutions.^{[16]} The automatic method not only reduced the time of manual extraction but also improved the overall calibration accuracy.
In summary, the automatic calibration technique discussed in this article has an accuracy of 0.94 mm and a precision of 2.67 mm using threeimage input. The first US image input delivered accuracy comparable with that of the other proposed methods. Carbajal et al.^{[31]} achieved mean accuracy at 1.02 mm with the iterative method. Pagoulatos et al.^{[10]} used the same closedform solution for calibrating a phasearray probe and obtained 1.54mm mean accuracy at 9cm US depth.
Since temporal calibration was not performed in this study, the efficiency provided is not final. However, the results are promising, and the proposed method seems to have potential for automatic and realtime calibration.
In this proposed design, the middle wires can move more toward the lateral direction compared to the coaligned Nfiducial's strands. The SVD closedform solutions ^{[16]} have the constraint of homogenous scaling in 2D data. The US machine could provide anisotropic scaling. Improved accuracy was achieved by the algorithm, which provided nonconstrained scale homogeneity. Further studies aimed at improving the design and optimizing the calibration techniques for easier automatic extraction are currently in progress.
Note that the evaluation of probes other than linear array probe was not performed in this article since it is out of the focus of our aimed application. As a result, our prototype is designed to be suitable for small fieldofview probes and fit practical calibration with probe depth selection. However, there are two main concepts proposed in this article:first, alternative ways to align the multiple Nfiducials while changing either width, height, or shifting the placement of individual Nfiducial, and second, the RANSAC model estimation to fit both 2D designated patterns and 3D planar alignment of the scanned point set.
Conclusions   
In this study, we have proposed a new design of the Nfiducial phantom and a corresponding extraction algorithm which facilitates an automatic and fast US calibration. Varying sizes of Nfiducials were used in this design so that, in US scan, they align in an obvious pattern for detection algorithms. The corresponding modelbased automatic extraction algorithm was also developed to accurately extract the feature of the phantom in US images. The proposed design and algorithm enable automatic calibration with Nfiducial phantom with significant efficiency.
Acknowledgment
N.T. would like to thank Ramathibodi Hospital staff and radiologist to give me a lot of knowledge and recommendation to develop this system. Another group N.T. would like to thank BART LAB members who give me much knowledge in mechanical design, programming, and statistical analysis.
This project is supported by Thailand's National Research Universities Grant through Mahidol University. The first author is supported by the Mahidol University PhDMD program and FY2016 Thesis Grant for Doctoral Degree Student under NRCT through Mahidol University. The US machine and experimental space are supported by Ramathibodi Hospital, Mahidol University. The authors would like to thank Cholatip Wirakapun, MD, for providing medical aspects of this research.
Financial support and sponsorship
This research has been funded by the ComputerIntegrated Intelligent Medical System Project under the National Research University Grant through Mahidol University. Another fund resource is the Integration of Surgical Navigation and Surgical Robotics for Breast Biopsy in Breast Cancer using Mammogram and US Images on Breast Mathematical Model Project under the Government Research Budget through Mahidol University (Grant No. 1112558). The first author is supported by the Mahidol University PhDMD program and FY2016 Thesis Grant for Doctoral Degree Student under the National Research Council of Thailand (NRCT) through Mahidol University.
Conflicts of interest
There are no conflicts of interest.
References   
1.  Chung SW, Shih CC, Huang CC. Freehand threedimensional ultrasound imaging of carotid artery using motion tracking technology. Ultrasonics 2017;74:1120. 
2.  Brambati B, Oldrini A, Lanzani A. Transabdominal chorionic villus sampling: A freehand ultrasoundguided technique. Am J Obstet Gynecol 1987;157:1347. 
3.  Yagel S, Cohen SM, Shapiro I, Valsky DV. 3D and 4D ultrasound in fetal cardiac scanning: A new look at the fetal heart. Ultrasound Obstet Gynecol 2007;29:8195. 
4.  Zhang H, Banovac F, White A, Cleary K. Freehand 3D Ultrasound Calibration Using an Electromagnetically Tracked Needle. Vol. 6141; 2006. p. 61412M. 
5.  Gebhard RE, Eubanks TN, Meeks R. Threedimensional ultrasound imaging. Curr Opin Anaesthesiol 2015;28:5837. 
6.  Cool D, Sherebrin S, Izawa J, Chin J, Fenster A. Design and evaluation of a 3D transrectal ultrasound prostate biopsy system. Med Phys 2008;35:4695707. 
7.  Chen TK, Thurston AD, Ellis RE, Abolmaesumi P. A realtime freehand ultrasound calibration system with automatic accuracy feedback and control. Ultrasound Med Biol 2009;35:7993. 
8.  Hsu PW, Prager RW, Gee AH, Treece GM, Freehand 3D Ultrasound Calibration: A Review. In: Sensen C.W., Hallgrímsson B. (eds) Advanced Imaging in Biology and Medicine. Springer, Berlin, Heidelberg; 2009. 
9.  Comeau RM, Fenster A, Peters TM. Integrated MR and ultrasound imaging for improved image guidance in neurosurgery. In: Proceedings of the SPIE; 1998. p. 74754. 
10.  Pagoulatos N, Haynor DR, Kim Y. A fast calibration method for 3D tracking of ultrasound images using a spatial localizer. Ultrasound Med Biol 2001;27:121929. 
11.  Lindseth F, Tangen GA, Langø T, Bang J. Probe calibration for freehand 3D ultrasound. Ultrasound Med Biol 2003;29:160723. 
12.  Zhang WY, Rohling RN, Pai DK. Surface extraction with a threedimensional freehand ultrasound system. Ultrasound Med Biol 2004;30:146173. 
13.  Chen TK, Abolmaesumi P, Pichora DR, Ellis RE. A system for ultrasoundguided computerassisted orthopaedic surgery. Comput Aided Surg 2005;10:28192. 
14.  Hsu PW, Prager RW, Gee AH, Treece GM. Realtime freehand 3D ultrasound calibration. Ultrasound Med Biol 2008;34:23951. 
15.  Sato Y, Nakamoto M, Tamaki Y, Sasama T, Sakita I, Nakajima Y, et al. Image guidance of breast cancer surgery using 3D ultrasound images and augmented reality visualization. IEEE Trans Med Imaging 1998;17:68193. 
16.  Arun KS, Huang TS, Blostein SD. Leastsquares fitting of two 3d point sets. IEEE Trans Pattern Anal Mach Intell 1987;9:698700. 
17.  Boctor E, Viswanathan A, Choti M, Taylor RH, Fichtinger G, Hager G. A novel closed form solution for ultrasound calibration. 2004 2 ^{nd} IEEE International Symposium on Biomedical Imaging: Nano to Macro (IEEE Cat No 04EX821). Vol. 2; 2004. p. 52730. 
18.  Boctor EM, Iordachita I, Choti MA, Hager G, Fichtinger G. Bootstrapped ultrasound calibration. Stud Health Technol Inform 2006;119:616. 
19.  Ameri G, McLeod AJ, Baxter JS, Chen EC, Peters TM. Line fiducial material and thickness considerations for ultrasound calibration. In: Proceedings of SPIE; 2015. p. 941529. 
20.  Bouchet LG, Meeks SL, Goodchild G, Bova FJ, Buatti JM, Friedman WA. Calibration of threedimensional ultrasound images for imageguided radiation therapy. Phys Med Biol 2001;46:55977. 
21.  Lasso A, Heffter T, Rankin A, Pinter C, Ungi T, Fichtinger G. PLUS: Opensource toolkit for ultrasoundguided intervention systems. IEEE Trans Biomed Eng 2014;61:252737. 
22.  ŚwiaątekNajwer E, Krysztoforski K, Dragan SŁ, Będziński R. The investigation of the lower limb geometry using 3D sonography and magnetic resonance. Meas J Int Meas Confed 2012;45:70210. 
23.  Rao Y, Wenhua L, Xingyuan Z, Peng Y. 3D Calibration algorithm for freehand ultrasound image. Inf Technol J 2013;12:255460. 
24.  Fischler MA, Bolles RC. Random sample consensus: A paradigm for model fitting with applications to image analysis and automated cartography. Commun ACM 1981;24:38195. 
25.  Comeau RM, Sadikot AF, Fenster A, Peters TM. Intraoperative ultrasound for guidance and tissue shift correction in imageguided neurosurgery. Med Phys 2000;27:787800. 
26.  Wei Z, Wan G, Gardi L, Mills G, Downey D, Fenster A. Robotassisted 3DTRUS guided prostate brachytherapy: System integration and validation. Med Phys 2004;31:53948. 
27.  Amin D V, Kanade T, Jaramaz B, DiGioia AM, Nikou C, LaBarca RS, et al. Calibration method for determining the physical location of the ultrasound image plane BT. In: Niessen WJ, Viergever MA, editors. Medical Image Computing and ComputerAssisted Intervention – MICCAI 2001. Berlin, Heidelberg: Springer Berlin Heidelberg; 2001. p. 9407. 
28.  Dosse MB, Ten Berge J. Anisotropic orthogonal procrustes analysis. J Classif 2010;27:11128. 
29.  Huynh DQ. Metrics for 3D rotations: Comparison and analysis. J Math Imaging Vis 2009;35:15564. 
30.  Mercier L, Langø T, Lindseth F, Collins DL. A review of calibration techniques for freehand 3D ultrasound systems. Ultrasound Med Biol 2005;31:44971. 
31.  Carbajal G, Lasso A, Gómez A, Fichtinger G. Improving Nwire phantombased freehand ultrasound calibration. Int J Comput Assist Radiol Surg 2013;8:106372. 
[Figure 1], [Figure 2], [Figure 3], [Figure 4], [Figure 5], [Figure 6], [Figure 7], [Figure 8], [Figure 9]
[Table 1], [Table 2]
