A dual-view digital tomosynthesis imaging technique for improved chest imaging Yuncheng Zhong, Chao-Jen Lai, Tianpeng Wang, and Chris C. Shawa) Department of Imaging Physics, The University of Texas MD Anderson Cancer Center, Houston, Texas 77054

(Received 23 January 2015; revised 11 June 2015; accepted for publication 22 July 2015; published 13 August 2015) Purpose: Digital tomosynthesis (DTS) has been shown to be useful for reducing the overlapping of abnormalities with anatomical structures at various depth levels along the posterior–anterior (PA) direction in chest radiography. However, DTS provides crude three-dimensional (3D) images that have poor resolution in the lateral view and can only be displayed with reasonable quality in the PA view. Furthermore, the spillover of high-contrast objects from off-fulcrum planes generates artifacts that may impede the diagnostic use of the DTS images. In this paper, the authors describe and demonstrate the use of a dual-view DTS technique to improve the accuracy of the reconstructed volume image data for more accurate rendition of the anatomy and slice images with improved resolution and reduced artifacts, thus allowing the 3D image data to be viewed in views other than the PA view. Methods: With the dual-view DTS technique, limited angle scans are performed and projection images are acquired in two orthogonal views: PA and lateral. The dual-view projection data are used together to reconstruct 3D images using the maximum likelihood expectation maximization iterative algorithm. In this study, projection images were simulated or experimentally acquired over 360◦ using the scanning geometry for cone beam computed tomography (CBCT). While all projections were used to reconstruct CBCT images, selected projections were extracted and used to reconstruct single- and dual-view DTS images for comparison with the CBCT images. For realistic demonstration and comparison, a digital chest phantom derived from clinical CT images was used for the simulation study. An anthropomorphic chest phantom was imaged for the experimental study. The resultant dual-view DTS images were visually compared with the single-view DTS images and CBCT images for the presence of image artifacts and accuracy of CT numbers and anatomy and quantitatively compared with root-mean-square-deviation (RMSD) values computed using the digital chest phantom or the CBCT images as the reference in the simulation and experimental study, respectively. High-contrast wires with vertical, oblique, and horizontal orientations in a PA view plane were also imaged to investigate the spatial resolutions and how the wire signals spread in the PA view and lateral view slice images. Results: Both the digital phantom images (simulated) and the anthropomorphic phantom images (experimentally generated) demonstrated that the dual-view DTS technique resulted in improved spatial resolution in the depth (PA) direction, more accurate representation of the anatomy, and significantly reduced artifacts. The RMSD values corroborate well with visual observations with substantially lower RMSD values measured for the dual-view DTS images as compared to those measured for the single-view DTS images. The imaging experiment with the high-contrast wires shows that while the vertical and oblique wires could be resolved in the lateral view in both singleand dual-view DTS images, the horizontal wire could only be resolved in the dual-view DTS images. This indicates that with single-view DTS, the wire signals spread liberally to off-fulcrum planes and generated wire shadow there. Conclusions: The authors have demonstrated both visually and quantitatively that the dual-view DTS technique can be used to achieve more accurate rendition of the anatomy and to obtain slice images with improved resolution and reduced artifacts as compared to the single-view DTS technique, thus allowing the 3D image data to be viewed in views other than the PA view. These advantages could make the dual-view DTS technique useful in situations where better separation of the objects-of-interest from the off-fulcrum structures or more accurate 3D rendition of the anatomy are required while a regular CT examination is undesirable due to radiation dose considerations. C 2015 American Association of Physicists in Medicine. [http://dx.doi.org/10.1118/1.4928214] Key words: digital tomosynthesis, chest imaging, chest radiology, tomosynthesis

5238

Med. Phys. 42 (9), September 2015

0094-2405/2015/42(9)/5238/14/$30.00

© 2015 Am. Assoc. Phys. Med.

5238

5239

Zhong et al.: Dual-view DTS imaging

1. INTRODUCTION Chest radiography is widely used as a preliminary imaging procedure for screening and diagnosing abnormalities in the thorax. However, chest radiography, as a projection imaging technique, produces two-dimensional images with overlapping anatomical structures that may impede the detection and visualization of subtle abnormalities. Although computed tomography (CT) has been widely considered as the gold standard for 3D chest x-ray imaging and can easily eliminate the overlapping problem, its significantly higher cost and radiation dose to the patient preclude its use for screening or preliminary diagnosis. An alternative technique, digital tomosynthesis (DTS), was developed to render crude 3D images that can be reviewed in the posterior–anterior (PA) view to reduce the effects of overlapping anatomy. Since its inception, chest DTS has continuously drawn interest in research and clinical applications. Clinical evaluations have shown that PA view DTS is superior to the twoview (PA and lateral) chest radiography with a comparable level of patient exposure.1–9 However, because limited angle

5239

projection data are used for reconstruction, the resolution of DTS images, although reasonably high in planes parallel to the detector (PA view), is intrinsically poor in the depth direction (perpendicular to the PA view planes). Thus, the images may only be displayed and viewed in the PA view while sagittal or axial views are essentially unusable. Furthermore, artifacts due to limited scanning angle and sampling appear in the fulcrum plane as high-contrast signals are spilled over from off-fulcrum planes with their strength decreasing with the distance from the fulcrum. These artifacts create another challenge in the interpretation of the DTS images. Increasing the angular range of DTS scans has been theoretically and experimentally shown to improve the resolution to a limited degree in the depth direction.10 Various configurations of DTS scans have been explored and their effects on the quality of reconstructed DTS images have been investigated.11,13 In this paper, we describe and demonstrate dual-view DTS technique with both simulation and experimental studies. With the dual-view DTS technique, projection images are acquired in two orthogonal views, PA and lateral, and then reconstructed with an iterative algorithm to improve the accuracy of the reconstructed volume image data for more accurate rendition of the anatomy and slice images with improved resolution and reduced artifacts, thus allowing the 3D image data to be viewed in views other than the PA view. To demonstrate these benefits, we measured and compared the signal spread of high-contrast wires from the fulcrum plane into the off-fulcrum planes in the dual-view DTS, single-view DTS, and cone beam CT (CBCT) images. Dual-view DTS images of a digital chest phantom and an anthropomorphic chest phantom were generated and compared with single-view DTS images and CBCT images for anatomical accuracy, image artifacts, and accuracy of CT numbers. 2. MATERIALS AND METHODS 2.A. DTS geometry

F. 1. Cone beam CT images were reconstructed from projections acquired over 360◦. (a) A subset of cone beam CT projections centered around the PA direction (solid lines) were extracted for reconstruction of the single-view DTS images. (b) For dual-view DTS images, an additional subset of projections centered around the lateral direction (dotted lines) were added for reconstruction. The angular range for scanning varies from 30◦ to 60◦. Medical Physics, Vol. 42, No. 9, September 2015

Various tomosynthesis geometries have been proposed and investigated.13 In addition to linear source motion with which the x-ray source moves in parallel to the detector plane, the x-ray source may follow partial isocentric motion along an arc with the center of rotation placed inside or close to the object while the detector either remains stationary or moves linearly in direction perpendicular to the rotation axis.13 We implemented DTS on a bench-top experimental cone beam CT system (Fig. 1), in which the phantom rotates while the x-ray source and the detector remain stationary, simulating the scanning geometry used by clinical CBCT or CT systems with which the x-ray source and detector gantry rotates around the object. This scanning geometry differs from typical implementation of DTS techniques but allows us to compare various DTS techniques and the CBCT technique using the projections obtained from the same scan. 2.B. Reconstruction

For CBCT with complete projection image data (300 views over 360◦), we used the well-known Feldkamp-Davis-Kress

5240

Zhong et al.: Dual-view DTS imaging

5240

(FDK) filtered back-projection algorithm for reconstruction to keep the computation task more manageable.14,15 For DTS images, we used the maximum likelihood expectation maximization (MLEM) iterative algorithm with a priori constraints.16–20 With this algorithm, the formation of x-ray images is modeled as a statistical process in which the signal on each detector pixel fluctuates independently following the Poisson probability distribution. The likelihood function is represented by the product of probability values of all pixels as follows: L(µ) =

 Ii (µ) Iim e−Ii (µ) Iim !

i

,

(1)

where L is the likelihood function Iim is the x-ray quanta measured on the ith detector pixel, and Ii (µ) is the expected xray quanta on the ith detector pixel due to the object function of µ. Because x-ray quanta are largely proportional to the image signals, they may be linearly converted from the image signals. For actual implementation, the logarithm of the likelihood function, ln L(µ), is maximized. This leads to the derivation of the following iteration process:21  µ(n+1) = k

i µ(n) k

 (n) − ai j µ j

aik * I0e j , Iim aik

+ -,

F. 2. Flow chart of iterative reconstruction used in this study.

(2)

i

where µk and µk are the kth voxel for the estimated object from the nth and (n +1)th iteration, respectively, and aij is the weighting factor that measures the contribution of the jth voxel, µ j , to x-ray attenuation in computing the projection signal for the ith pixel. The ratio of the weighted sum of the estimated projection values through the kth voxel to that of the measured values is used as the correction factor to adjust the value of µk (n). Notice that the correction described by Eq. (2) is only for a specific view. For complete iteration, this correction process is repeated for all views. To allow the iteration process to converge sooner, we incorporated the generalized Gaussian Markov random field constraint, computed as the logarithm of the sum of the absolute differences between any two neighboring voxel values and added to ln L(µ) for maximization.22 This was incorporated to help minimize the variations of the attenuation values of neighboring voxels and to suppress noise in the reconstructed images while preserving the edge details, (n)

(n+1)

 (n) − ai j µ

j + aik * I0e j i (n+1) (n) , µk = µk   r −1 ,  m r Ii aik + I0r λ w pq µ p − µq



i

(3)

{ p,q} ∈C

where the additional term in the denominator represents the constraint. The {p, q} are counted for the neighbors pairs of voxels. The weighting parameter w pq = 1 for the nearest neighbors and 1.4 for the second nearest neighbors. The parameters were chosen as λ = 1 and r = 1.1 in the reconstruction. This iterative process is repeated until the differences between the estimated object and that from previous iteration reaches a preset low threshold level. To determine the number of iterations used, we computed and plotted the root-mean-squaredeviation (RMSDs) between images from two successive iterMedical Physics, Vol. 42, No. 9, September 2015

ations and plotted them as a function of the iteration number. It was found that the RMSDs decreased with each iteration but stayed at the 5%–8% level with 30 or more iterations. We therefore chose to use 30 iterations for reconstruction of all DTS images in the study. A flow chart for implementation of the above described reconstruction algorithm is shown in Fig. 2. An arbitrary value was assigned to all voxels as the initial estimate of the object from which forward projections were computed. The ratios of the weighted sum of the computed projection values to that of measured projection values were calculated and used to correct the object on a voxel-by-voxel basis for all views. 2.C. Simulation with a digital chest phantom

A digital chest phantom was used to simulate the dual-view DTS imaging. The digital chest phantom used was constructed from three-dimensional (3D) chest images obtained using a clinical CT system (LightSpeed16, GE Medical Systems, Inc., Milwaukee, WI) with 120 kVp x-rays. To simulate CBCT and DTS imaging, projection images through the digital chest phantom were first computed for 360 evenly spaced angular views for 360◦ scan around the chest. Projection image signals were computed by integrating x-ray attenuation along various x-ray paths through the object as follows: Ii (µ) = I0e

 − ai j µ j j

,

(4)

where Ii is the calculated x-ray quanta detected by the ith pixel in a projection image, I0 is the unattenuated x-ray quanta detected, µ j is the estimated linear attenuation coefficient for the jth voxel in the 3D matrix representing the object, and aij is the projection matrix that weighs the contribution of µ j

5241

Zhong et al.: Dual-view DTS imaging

to the total attenuation of x-rays going through the jth voxel in the object and falling on the ith pixel of the detector for a specific view. Algorithms have been developed to calculate the x-ray projections with pixel-driven,23 voxel-driven, or distance-driven algorithms.24 We adopted an algorithm with which the x-ray attenuation by a voxel is first projected onto pixels covered by the x-rays going through the voxel. The next step is to combine x-ray attenuations from all contributing voxels for each different pixel to form the projection image.25,26 This algorithm was also used in the iterative algorithm for reconstruction of DTS images in our study. The imaging geometry for simulation was configured to be similar to that used in the imaging experiments (Fig. 1) with an isocenter-to-image distance and an x-ray source-to-image distance of 20 and 180 cm, respectively. Notice that the scanning geometry shown in Fig. 1 differs from those commonly used in the implementation of DTS with commercial x-ray imaging systems. The geometry allowed our bench-top CBCT system to be used to acquire projection images for use in reconstruction of both CBCT and DTS images. The detector was assumed to have a 43 × 43 cm2 image area covered by a 736 × 736 image matrix with a pixel size of 0.585 mm chosen as follows: The voxel size of the clinical CT images (0.5 mm) was first projected to the detector plane with a magnification factor of 1.125 (180/160 cm) to be 0.5625 mm, which turned out to be close to the detector pixel size for 3 × 3 binning mode (3 × 0.2 mm = 0.6 mm). However, we chose the pixel size of 0.585 mm to obtain an image matrix size (736×736 for 43×43 mm) in multiples of 16×16 to allow GPU to be used to speed up reconstruction. Single-view DTS was implemented with a 60◦ scan centered around the depth (PA) direction and an angular step of 2.0◦, resulting in 31 projections. An additional lateral scan was added for simulation of dual-view DTS. The first set of dual-view DTS images were simulated using an angular range of 60◦ for each view and an angular increment of 2.0◦ or 4.0◦, resulting in 62 or 32 projections, respectively. The second set of images were simulated using an angular range of 30◦ for each view and an angular increment of 1.0◦ or 2.0◦, also resulting in 62 and 32 projections, respectively. This allowed us to compare dual-view DTS techniques with

5241

different scanning angles (30◦ or 60◦) and different numbers of projections (32 or 62). Both visual inspection and quantitative evaluations were performed to compare DTS and CBCT images with the digital chest phantom which was used as the reference images in the simulation study. The reconstructed images were first visually compared with the reference images for the presence of artifacts, accuracy in anatomical appearance and voxel values (CT numbers). For quantitative evaluation, the RMSD of the voxel values from those of the reference images was calculated over four representative regions-of-interest (ROIs) as follows:     (HU j − HU0j )2 RMSD =

j=1, N

N

,

(5)

where HU j and HU0j were the Hounsfield units (CT numbers) for the jth voxel in the reconstructed image and the reference image, respectively, and N was the number of the voxels included for comparison. The ROIs were 3D volumes in four anatomical structures: heart, lung, subdiaphragm, and mediastinum (Fig. 3). A smaller RMSD value would indicate that the reconstructed images were a better representation of the object. 2.D. Imaging experiment with an anthropomorphic chest phantom

For the experimental study, the CBCT images were used as the reference images as they are expected to provide a reasonably good 3D representation of the phantom. Various subsets of the projection image set for constructing the CBCT images were extracted to construct the single-view DTS images and four different sets of dual-view DTS images, corresponding to four different combinations of angular ranges and increments. Both visual inspection and RMSD calculations were conducted to evaluate and compare the two types of DTS images to CBCT images. A bench-top CBCT system was constructed with an xray tube (G1592, Varian Medical Systems, Inc., Palo Alto, CA) and a flat panel detector (XRD 1621, PerkinElmer, Inc.,

F. 3. Regions-of-interest for computing the RMSD values for quantitatively comparing the accuracy of CT numbers of the DTS images with the CBCT images and the digital phantom images (for the simulation study only). (a) and (b) are from different depth in depth (PA) direction. Medical Physics, Vol. 42, No. 9, September 2015

5242

Zhong et al.: Dual-view DTS imaging

Waltham, MA) mounted on an optical bench with a sourceto-image distance of 180 cm. Images were acquired in the 2 × 2 binning mode, resulting in 1024 × 1024 images with a pixel size of 400 µm. This acquisition mode was widely used in operating clinical CBCT systems. The resulting smaller image size also helped reduce the computation time for image reconstruction. A stepper motor-driven rotating stage was mounted between the x-ray tube and the detector with an isocenter-to-image distance of 23 cm to hold and rotate the phantom to simulate CBCT chest imaging with continuous gantry motion. The scan was performed with continuous x-ray exposures at 120 kVp and 36 mA. Thus, the total exposure used was 1.2 mAs per projection image, which was comparable to that of clinical CBCT level of 1.2 mAs per projection view.27 Three hundred projection images of an anthropomorphic chest phantom (Radiology Support Devices, Long Beach, CA) were acquired over a complete rotation (360◦) with an angular increment of 1.2◦. Single-view DTS projections were extracted from the projections every 2.4◦ within a 60◦ angular range in the depth (PA) direction, resulting in 26 projections. For dual-view DTS, same number of additional projections were extracted over the same angular range in the lateral direction and added to the projection image set. For dual-view DTS, a sampling angle of 2.4◦ or 1.2◦ was used to extract 26 projection views over a scanning angle of 60◦ or 30◦, respectively, centered around both PA and lateral directions. Alternatively, a sampling angle of 4.8◦ or 2.4◦ was used to extract 13 projection views over a scanning angle of 57.6◦ or 28.8◦, respectively, centered around both PA and lateral directions. The latter two

5242

scanning angles are close enough to the first two so that they are referred to as 60◦ and 30◦ as well. 2.E. Imaging experiment with high-contrast wires

The main advantage of dual-view DTS, as envisioned in this paper, is its ability to improve the accuracy of the reconstructed volume image data, to achieve more accurate rendition of the anatomy and to obtain slice images with improved resolution and reduced artifacts. This would allow 3D images to be viewed in the lateral (sagittal) or axial (superior to inferior) views. In addition to demonstrating this advantage with simulation and experimental imaging studies as described in Secs. 2.C–2.D, we imaged three high-contrast linear objects and studied how well their shapes were retained in CBCT and DTS imaging with various techniques. The objects imaged were three straight segments of metal wires with a diameter of 0.254 mm embedded in a Lucite block.12 All wires lied in a PA view plane, one oriented horizontally, in parallel to the circular trajectory of the focal spot; the second one oriented vertically, transverse to the trajectory; and the third one at oblique angle (45◦) to the trajectory. The projection images were acquired in the same way as the projection images of the anthropomorphic chest phantom were acquired in the experimental imaging study. Signal profiles across the wire object were extracted in the PA view plane in which the wires lied and in a lateral view plane to compare how well the shape of the wires were retained and how the wire signals spread to offfulcrum planes.

F. 4. Coronal view slice images of a digital chest phantom: (a) original CT image representing the phantom, (b) simulated CBCT image, (c) simulated dual-view DTS image (31/60◦ + 31/60◦), and (d) simulated single-view DTS image (31/60◦). Medical Physics, Vol. 42, No. 9, September 2015

5243

Zhong et al.: Dual-view DTS imaging

3. RESULTS 3.A. Simulation study

A coronal view image representing a slice of the digital chest phantom and those obtained with the CBCT, the singleview DTS, and the dual-view DTS techniques are shown for comparison in Fig. 4. The single-view DTS images were obtained with 31 views over 60◦ (referred to as 31/60◦) centered around the depth (PA) direction while the dual-view DTS images were obtained with additional 31 views over 60◦ centered around the lateral direction (referred to as 31/60◦ + 31/60◦). The original phantom image [Fig. 4(a)] appeared to be accurately represented by the CBCT image [Fig. 4(b)]. The dualview DTS image [Fig. 4(c)] is less accurate than the CBCT image, but significantly more accurate than the single-view DTS image [Fig. 4(d)] in depicting the anatomical details. The overlapping of the off-fulcrum structures was more severe in the single-view DTS image [Fig. 4(d)], as more off-fulcrum objects appeared as artifacts, thereby biasing image signals and degrading the visibility of details in the fulcrum. For instance, image signals were darker in some regions [Fig. 4(d)], reflecting the spillover of off-fulcrum air signals into the fulcrum. Overall, the dual-view DTS images were a significantly more accurate representation of the original phantom than the single-view DTS images. A sagittal view image representing a slice of the digital chest phantom and those obtained with the CBCT, the singleview DTS, and the dual-view DTS techniques are shown for

F. 5. Sagittal view slice images of a digital chest phantom: (a) original CT image representing the phantom, (b) simulated CBCT image, (c) simulated dual-view DTS image (31/60◦ + 31/60◦, 31 projections for each view), and (d) simulated single-view DTS image (31/60◦). Medical Physics, Vol. 42, No. 9, September 2015

5243

comparison in Fig. 5. Again, the original phantom image [Fig. 5(a)] appeared to be accurately represented by the CBCT image [Fig. 5(b)]. The dual-view DTS image [Fig. 5(c)] appeared to be less accurate than the CBCT image, but significantly more accurate than the single-view DTS image [Fig. 5(d)] in depicting the anatomical details. Like the coronal view single-view DTS image, there was much loss of details and artifacts from off-fulcrum structures. Unlike the coronal view single-view DTS image, the sagittal view image appeared to be severely distorted and blurred at such a severe level that it appeared unusable for regular imaging tasks. The coronal view single-view DTS image, on the other hand, was less blurry and more informative than the sagittal view images. Thus, current DTS chest images [based on single-view scan in the depth (PA) direction] are usually displayed in the coronal view only. To study the effects of scanning and sampling angles, dualview DTS images were constructed with four different combinations of scanning angle and angular sampling: 31 or 16 views evenly spaced over 60◦ and 30◦ centered around both PA and lateral directions. Example coronal view and sagittal view slice images from the results are shown in Figs. 6 and 7, respectively. For the coronal view, the 31/60◦ + 31/60◦ image in Fig. 6(a) appears to have the most details and the best contrast. The 16/60◦ + 16/60◦ image in Fig. 6(b) has lower contrast and slightly less details. Compared to these two 60◦ images, the details in both 30◦ images [Figs. 6(c) and 6(d)] were subject to much more spillover from off-fulcrum structures. For instance, shadows of the rib bones are visible in the lung area in the 30◦ images but not in the 60◦ images. The vessels are shown as short segments in the lung area in the 30◦ images but only cross sections of the vessels are shown in the 60◦ images. For the sagittal view, the 31/60◦ + 31/60◦ image in Fig. 7(a) appears to have more and better defined details than the 16/60◦ + 16/60◦ image in Fig. 7(b). Compared to these two 60◦ images, the two 30◦ images Figs. 7(c) and 7(d) were missing more details and subject to more spillover from off-fulcrum structures. For instance, the lower lung borders in them look distinctively different from those in the two 60◦ images which look more similar to those in the digital phantom [Fig. 5(a)] and CBCT image [Fig. 5(b)]. The vessels in the lung area are also shown as significantly longer segments with lower contrast compared to those in the 60◦ images. Both observations indicate more spillovers of off-fulcrum objects into the slice in fulcrum. Table I lists the scanning and sampling parameters used in simulation of the five different DTS techniques and the associated RMSD values evaluated over different ROIs in the heart, lung, subdiaphragm, and mediastinum (Fig. 3). Also listed are the average RMSD values over all the ROIs. The first observation is that the average RMSD value for single-view DTS (258.8) was 13.1 times greater than that for CBCT (19.8) while those for dual-view DTS were all lower (84.8–156.5). It appears that the sampling angle did not significantly affect the average RMSD values (84.8 versus 90.0 and 156.5 versus 154.8 for 60◦ and 30◦, respectively) while the scanning angle did (84.8 and 90.0 for 60◦ versus 156.5 and 154.8 versus 156.5 for 30◦) for dual-view DTS. Using the dual-view DTS techniques, the average RMSD value with the single-view

5244

Zhong et al.: Dual-view DTS imaging

5244

F. 6. Effects of scanning angle and number of projections on coronal view dual-view DTS slice images: (a) 31/60◦ + 31/60◦, (b) 16/60◦ + 16/60◦, (c) 31/30◦ + 31/30◦, and (d) 16/30◦ + 16/30◦.

DTS technique was reduced by a factor of 2.9–3.1 with a scanning angle of 60◦ and 1.7 with a scanning angle of 30◦. The reduction of the RMSD values apply to all ROIs and all dualview DTS techniques except a slight increase of the RMSD values with the two 30◦ techniques in the lung regions. 3.B. Experimental study with an anthropomorphic chest phantom

The reconstructed coronal view CBCT and DTS images of an anthropomorphic chest phantom are shown in Fig. 8 for comparison. The CBCT image [Fig. 8(a)], used as the reference for comparison, appears to have the best clarity and well defined anatomical details with only faintly visible shadows of ribs. The shape and dimension of anatomical structures in the two 60◦ dual-view DTS images [Figs. 8(c) and 8(d)] closely match those in the CBCT image. However, many small vessels are missing in these images. There are also faint shadows of the ribs and heart, which are slightly more visible than those in the CBCT image. There also appear to be substantial artifacts, including a dark shadow covering the top 20% in all images and streaking artifacts present in the top third and right side in almost all DTS images. The former may have been generated because the detector was not wide enough to cover the entire phantom in views close to the depth (PA) direction. The latter may have been caused by attempting to use widely separated projection views to reconstruct high-contrast objects: sharp edges of the shoulder parts in the phantom for artifacts in the top third of the image and the rib bone near the bottom of the left lung for artifacts in the lower right part of the image. While the artifacts due to truncated coverage of the phantom may be alleviated if padding extension was applied, the streaking artifacts can only be reduced by using smaller sampling angle Medical Physics, Vol. 42, No. 9, September 2015

or totally eliminated by using a complete set of projections in the case of CBCT. The shapes and dimensions of the single-view [Fig. 8(b)] and the two 30◦ dual-view DTS images [Figs. 8(e) and 8(f)] appear to be distorted as compared to those in the CBCT image. For instance, the lungs appear to be wider and taller while the mediastinum narrower. With the same window and level setting used for printing, the contrast of the lungs appears to be lower. Furthermore, the shadows of the heart and ribs are more prominent while the diaphragm is distorted and blurry as compared to that in the CBCT image. All these observations indicate that more off-fulcrum structures were spilled over and overlap with those in the fulcrum plane. The differences in preserving the anatomical details between the 30◦ and 60◦ dual-view DTS images demonstrate the important role the scanning angle plays in achieving accurate image reconstruction. Comparing all DTS images for presence of artifacts, the 26/30◦ + 26/30◦ image is almost artifact free except that there are some faint curvilinear artifacts over the top third of the image caused by the failure to cover the entire shoulder in views close to the depth (PA) direction. In contrast, the corresponding artifacts in the single-view DTS image and in the 13/30◦ + 13/30◦ and 26/60◦ + 26/60◦ dual-view DTS images, all obtained with the same sampling angle, are more widely spaced and those in the 13/60◦ + 13/60◦ dual-view DTS image have the widest spacing. The effects of angular spacing on the look of the artifacts are clearly reflected in these observations. The reconstructed sagittal view CBCT and DTS images are shown in Fig. 9 for comparison. Again, the CBCT image [Fig. 9(a)] appears to have the best clarity and well defined anatomical details, with only slight shadows of ribs near the spine. The shape and dimension of anatomical structures in

5245

Zhong et al.: Dual-view DTS imaging

5245

F. 7. Effects of scanning angle and number of projections on sagittal view dual-view DTS slice images: (a) 31/60◦ + 31/60◦, (b) 16/60◦ + 16/60◦, (c) 31/30◦ + 31/30◦, and (d) 16/30◦ + 16/30◦.

the two 60◦ dual-view DTS images [Figs. 9(c) and 9(d)] closely resemble those in the CBCT image. However, many small vessels are missing in these images. There are also faint shadows of the ribs (near the spine) and heart. There are also artifacts in the region where short sections of arms branch out from the torso in the phantom, forming high-contrast edges on both sides. These shadows are clearly visible but hardly recognizable as cross sections of arms. As with the two coronal view

60◦ dual-view DTS images [Figs. 8(c) and 8(d)], there appear to be substantial streaking artifacts in the upper third of the image and to a lesser extent in the lower right part of the image. These artifacts were generated due to the attempt to use a small number of widely spaced projection views for reconstruction of high-contrast objects: sharp cutoff edges on the two sides of the shoulders and the bone structure of the spine. Comparing the two 60◦ dual-view DTS images [Figs. 9(c) and 9(d)], the

T I. RMSD values for regions-of-interest in four representative regions (heart, lung, subdiaphragm, and mediastinum) in the reconstructed chest images. RMSD values Technique

Scanning technique

Sampling angle

Heart

Lung

Subdiaph.

Mediastinum

Average

Simulation

CBCT Single-view Dual-view Dual-view Dual-view Dual-view

360◦ 60◦ ◦ 60 + 60◦ 60◦ + 60◦ 30◦ + 30◦ 30◦ + 30◦

1.2◦ 2.0◦ 2.0◦ 4.0◦ 1.0◦ 2.0◦

16 262 84 94 76 78

25 185 87 102 210 205

16 164 58 56 67 65

22 424 110 108 273 271

19.8 258.8 84.8 90.0 156.5 154.8

Experiments

Single-view Dual-view Dual-view Dual-view Dual-view

60◦ 60◦ + 60◦ 57.6◦ + 57.6◦ 30◦ + 30◦ 28.8◦ + 28.8◦

2.4◦ 2.4◦ 4.8◦ 1.2◦ 2.4◦

209 82 68 105 103

165 70 76 152 157

180 76 84 105 125

181 82 82 183 179

183.8 77.5 77.5 136.3 141.0

Medical Physics, Vol. 42, No. 9, September 2015

5246

Zhong et al.: Dual-view DTS imaging

5246

F. 8. Coronal view slice images of an anthropomorphic chest phantom obtained using (a) CBCT technique, (b) single-view DTS technique (26/60◦), and dual-view DTS techniques with various scanning angles and number of projections: (c) 26/60◦ + 26/60◦, (d) 13/57.6◦ + 13/57.6◦, (e) 26/30◦ + 26/30◦, and (f) 13/28.8◦ + 13/28.8◦.

26 view image [Fig. 9(d)] is more blurry resulting in reduced contrast and more details lost. The use of a larger sampling angle (4.8◦ versus 2.4◦) also resulted in artifacts with more widely spaced curvilinear patterns. Again, this demonstrates the important role angular sampling played in tomosynthesis imaging. Although subjected to greater distortion and spillover of anatomical structures from off-fulcrum planes, the 26/30◦ + 26/30◦ dual-view DTS image [Fig. 9(e)] has the cleanest appearance (with least streaking artifacts) among all DTS images, followed by the 26/60◦ +26/60◦ dual-view DTS image [Figs. 9(c)], 13/30◦ +13/30◦ dual-view DTS image [Fig. 9(f)], and the single-view DTS image [Fig. 9(b)] while the 13/60◦ +13/60◦ dual-view DTS image [Fig. 9(d)] has the most prominent streaking artifacts. This indicates that an increase of the scanning angle without corresponding increase of the number of projections would result in serious streaking artifacts in both single- or dual-view DTS imaging. However, for the accuracy of the reconstructed 3D images, the two 60◦ dual-view DTS images [Figs. 9(c) and 9(d)] outperformed the two 30◦ dualview DTS images [Figs. 9(e) and 9(f)] and the single-view Medical Physics, Vol. 42, No. 9, September 2015

DTS image [Fig. 9(b)] in their ability to isolate the in-fulcrum structures from the off-fulcrum structures, as indicated by the less prominent presence of the shadows from the heart, ribs, the arms, and diaphragm outside the off-fulcrum planes. This may be attributed to the use of the larger scanning angle (60◦ versus 30◦) and the use of an additional scan in the lateral direction when comparing to the single-view DTS image. Table I lists the scanning and sampling parameters used in acquiring the projection images for the five different DTS techniques and the associated RMSD values evaluated over four different ROIs in the heart, lung, subdiaphragm, and mediastinum (Fig. 3) and their average. For the experimentally obtained phantom images, the CBCT images of the phantom were used as the reference in computing the RMSD values. This was justifiable as the RMSD values computed by comparing the CBCT images with the original phantom images in the simulation study were very low, indicating their close resemblance to each other. With the projection images for CBCT acquired with a sampling angle of 1.2◦, a sampling angle of 2.4◦ was used to extract 26 projection views over a

5247

Zhong et al.: Dual-view DTS imaging

5247

F. 9. Sagittal view slice images of an anthropomorphic chest phantom obtained using (a) CBCT technique, (b) single-view DTS technique (26/60◦), and dual-view DTS techniques with various scanning angles and number of projections: (c) 26/60◦ + 26/60◦, (d) 13/57.6◦ + 13/57.6◦, (e) 26/30◦ + 26/30◦, and (f) 13/28.8◦ + 13/28.8◦.

scanning angle of 60◦ in the depth (PA) direction for singleview DTS. As with the simulation study, the average RMSD values for all dual-view DTS techniques (77.5–141.0) were lower than those for the single-view DTS technique (183.8). It appears that angular sampling did not significantly affect the average RMSD values (77.5 versus 77.5 and 136.3 versus 141.0 for 60◦ and 30◦, respectively) while the scanning angle did (77.5 for 60◦ versus 136.3 and 141.0 for 30◦) for dualview DTS. Using the dual-view DTS techniques, the average RMSD value with the single-view DTS was reduced by a factor of 2.4 with a scanning angle of 60◦ and 1.3 with a scanning angle of 30◦. The reduction of the RMSD values applies to all ROIs and dual-view DTS techniques except that in lung and mediastinum regions, the reduction factors were only 1.1 with the two 30◦ techniques. 3.C. Experimental study with high-contrast wires

To evaluate how well the shapes of the wires are retained, we measured the linear profiles across the wires in the PA Medical Physics, Vol. 42, No. 9, September 2015

view plane containing all wires studied. To evaluate how much the wire signals spread to off-fulcrum planes, we measured the profiles in a lateral view plane containing the wires. Figures 11(a)–11(c) show the linear profiles in the PA view images across the vertically oriented wire, the obliquely oriented wire, and horizontally oriented wire, respectively. The wire profiles were measured at the locations shown by Fig. 10(a) for the PA view images and by Fig. 10(b) for the lateral view images. All profiles are displayed as unprocessed CT numbers without any rescaling or normalization to show the differences in the accuracy of CT numbers and image contrast between the different techniques used. All profiles from either single- or dual-view DTS images have the shape of peaked distributions that reflect x-ray attenuation through the wire cross sections and have similar full-widths at half maximum (FWHMs). Since the voxel size was 348 µm, the ideal diameter of the wire would be of one voxel, which was smaller than the measured FWHMs. Interestingly, the profiles from the CBCT image actually have slightly wider FWHMs. Although a high-contrast narrow wire object is relatively easy

5248

Zhong et al.: Dual-view DTS imaging

5248

F. 10. (a) PA view of the wire phantom illustrating how line profiles were measured. (b) Sagittal view illustrating how the peak line profiles were extracted.

to reconstruct, requiring only as few as two projections, the use of a large number of projections reconstructed may result in more alignment errors to be superimposed together, degrading the sharpness of the reconstructed wire object. Thus, although CBCT results in much more accurate reconstruction than DTS, it uses many more (300 or more) projections than DTS (26–52) thus the alignment error may accumulate to a higher level, reducing the sharpness of the reconstructed wire. The base signals (of the Lucite background) of all profiles are at a similar level (−300 to −200). The only exceptions are those around the horizontally oriented wire in the singleview DTS image [−580, Fig. 11(c)], which might be caused by the combined effects of less coverage, wire orientation,

and the high-contrast object. The main differences among the profiles lie in the wire contrasts, with the CBCT image showing the highest contrast, followed closely by the two 26+26 dual-view DTS images and then by the two 13+13 dual-view DTS images. The contrast generated by the single view DTS technique was similar to those generated by the two 13+13 dual-view DTS images for both the vertical and oblique wires but substantially lower for the horizontal wire. Figures 11(d)–11(f) show the linear profiles in the lateral view images across the vertically oriented wire, the obliquely oriented wire, and horizontally oriented wire, respectively. The plane corresponding to these images were chosen to contain the entire vertical wire and align with its central axis. The plane

F. 11. Original (un-normalized) peak line profiles across the (a) vertical, (b) oblique, and (c) horizontal wires as measured in the PA view plane containing all wires. And original (un-normalized) peak line profiles across the (d) vertical, (e) oblique, and (f) horizontal wires as measured in a lateral view plane containing or intersecting with all wires. Medical Physics, Vol. 42, No. 9, September 2015

5249

Zhong et al.: Dual-view DTS imaging

intersects with the oblique and horizontal wires at 45◦ and 90◦, respectively. Thus, the lateral view images contain only the cross sections of these wires, through whose centers horizontal signal profiles were measured along the depth (PA) direction. Thus, all profiles measure the maximum spread of the peak wire contrast from the fulcrum plane to off-fulcrum planes. All profiles, except that for horizontal wire from the singleview DTS image, have the same shape as although slightly wider than those from the PA view images, indicating that the wires could be resolved in the corresponding lateral view images. The single-view DTS image resulted in an almost flat profile across the horizontal wire, indicating that the wire could not be resolved at all in the depth (PA) direction. The CBCT image resulted in highest contrast for all three wires, closely followed by those from the two 26+26 dual-view DTS images and then the two 13+13 dual-view DTS images. For the vertical and oblique wires, the single-view DTS image resulted in about the same contrast as those in the two 13+13 dual-view DTS images, but no contrast at all for the horizontal wire.

4. DISCUSSIONS Our results show that compared with the single-view DTS technique with the same scanning and sampling angles, the dual-view DTS technique substantially improved the accuracy of reconstructed volume image data, resulting in more accurate rendition of the anatomy and slice images of better resolution and reduced artifacts, thus allowing the 3D image data to be displayed in views other than the PA view. This improvement can be explained by using the coverage of Fourier space determined by the angular range of the projection views used for reconstruction. In single-view DTS, the Fourier space along the depth (PA) direction is covered much less than that along the lateral direction, resulting in poor resolution along the depth (PA) direction [Figs. 5(d) and 9(b)]. In dual-view DTS, the additional scan data in the lateral view increase the coverage of the Fourier space thus resulting in more accurate 3D volume image data. This finding is corroborated by comparison of the single-view and dual-view DTS images and their corresponding RMSD values. As with the single-view DTS technique, the scanning angle and the sampling angle are the two most important parameters in implementing the dual-view DTS technique. The scanning angle is the major factor that impacts the resulting quality of the reconstructed images in single-view DTS imaging. While a larger scanning angle is desirable, it may require more projections to be acquired to keep the view angle reasonably sampled, which may also increase the radiation dose required. These considerations apply to dual-view DTS imaging as well. In addition, scanning in two orthogonal views would double the number of projections if the sampling angular is to be kept the same. This would double the radiation dose unless the exposure for each projection is cut in half. Alternatively, the sampling angle may be doubled or the scanning angle for each view cut in half to keep the number of projections the same. This way, both the total radiation dose and the expoMedical Physics, Vol. 42, No. 9, September 2015

5249

sure per projection may be kept the same. Our comparison of four different combinations of the scanning angle and number of projections provides a perspective on the image quality differences among these approaches. The comparison results show that, as with single-view DTS, the use of larger scanning angle has significantly higher impact on image quality, as judged by the accuracy and sharpness of anatomical details and the amount of artifacts present. Effects of angular step, which determines the number of projections, on image quality are hardly noticeable. However, in the upper half of the phantom images in Figs. 8 and 9, the artifacts generated by the incomplete coverage of the shoulder area are obviously significantly affected by the sampling angle. These artifacts are substantially reduced by using a smaller sampling angle, decreasing from 4.8◦ in Figs. 8(d) and 9(d) to 2.4◦ in Figs. 8(b), 8(c), 8(f), 9(b), 9(c), and 9(f), followed by the smallest angle, 1.2◦, in Figs. 8(e) and 9(e), in which the artifacts are barely noticeable. The experimental study with high-contrast wires was by no means a rigorous study on the imaging characteristics of the dual-view DTS imaging technique. However, we have demonstrated that linearly shaped objects lying in a PA view plane horizontally (parallel to the scanning direction) were unresolvable with the single-view DTS technique [Fig. 11(f)]. This resulted in severe degradation of image quality and image distortion in the lateral view display of the reconstructed 3D image data [Figs. 5(d) and 9(b)]. With the dual-view DTS technique, such objects were well resolved and allowed the reconstructed 3D image data to be displayed in the lateral view without distortion and with good clarity [Figs. 5(c) and 9(c)]. Although we have demonstrated the potential advantages of the dual-view DTS technique, this study is limited in that the effects of radiation dose were not considered. In our experiment, the CBCT scan was performed at an exposure level of 1.2 mAs per projection. Since 26 projections were used for singleview DTS, the total exposure was 31.2 mAs. The exposures for dual-view DTS were 31.2 or 62.4 mAs, depending on the configuration used. They are 5 or 10 times of typical exposure (6 mAs) for clinical DTS chest imaging or a lateral chest xray. Given the fact that both the artifacts and spillover of offfulcrum anatomy into fulcrum plane are a visible problem at typical exposure levels of clinical DTS, the advantages of dual-view DTS demonstrated at higher exposure levels should also hold at the lower exposure level of current clinical DTS. However, the reduction of artifacts and spillover of off-fulcrum structures may open the door to the visualization of more subtle details and higher exposure levels may become desirable because of this enhanced benefit. The lack of equalizing the radiation doses in generating the dual-view DTS images may not affect the comparison of the gross artifacts such as the spillover of the rib bones, heart chambers, or lung boarders from the off-fulcrum planes but may significantly affect the detection and visualization of subtle objects like lung nodules. Thus, further investigation of the characteristics and advantages of the dual-view DTS technique should include the radiation dose considerations. A more rigorous simulation study would require noise fluctuations be added to the projection images

5250

Zhong et al.: Dual-view DTS imaging

compute from the digital chest phantom. A more rigorous experimental study would require the dual-view DTS images be generated with the same radiation dose as the single-view DTS images. The dual-view DTS technique may be implemented with either single- or the dual-gantry approach. With the former, an x-ray source-detector gantry must be rotated around the patient from −30◦ to 30◦ [with 0◦ being the depth (PA) direction] for the PA scan and then rotated to a new position for the lateral scan from 60◦ to 120◦. Thus, a total rotation of 150◦ is required for the two scans. In the dual-gantry approach, the two x-ray source-detector gantries are oriented orthogonal to each other and rotated together for 60◦ to complete image acquisition. This concurrent rotation of the two gantries allows the scanning time to be substantially reduced and patient motion minimized. However, the dual-gantry approach would cost more to implement as it requires the use of two sets of x-ray source and detectors supported and rotated by two gantries. Furthermore, it also requires the complex task of aligning two gantries and registering projection images acquired from them. Like single-view DTS, dual-view DTS may be implemented by keeping the detector stationary while moving the x-ray source or by rotating both the x-ray source and detector mounted on a gantry. Alternatively, the object can be rotated while keeping the x-ray source-detector gantry stationary during the scan. Although this approach may be acceptable for imaging animals, it may be unacceptable for most patient imaging applications. The dual-view DTS technique may also be implemented with multiple stationary xray sources such as the carbon nano tubes (CNTs). With this approach, the tubes at different locations are turned alternately fired to acquire projection images at various view angles. Because there is no need to move a heavy and bulky conventional x-ray tube, the scanning time may be substantially reduced and there is no image blurring due to tube motion.28 Despite the advantages of the dual-view DTS technique demonstrated in this study, the clinical usefulness of the technique depends on several considerations. First of all, the dualview DTS technique may be useful in obtaining a more accurate and isotropic 3D rendition of the anatomy in situations that call for more flexible viewing especially in the lateral view but a regular CT examination would be undesirable due to cost or radiation dose considerations. Second, the dual-view DTS technique may be used when the separation of the infulcrum structures from the off-fulcrum structures needs to be improved to reduce the artifacts generated from the offfulcrum structures for a more clear view of structures-ofinterest.

5. CONCLUSIONS We have described and demonstrated a dual-view DTS technique for chest imaging. We have demonstrated with both simulation and image experiments that the dual-view DTS technique may be used to improve the accuracy of the reconstructed volume image data, to achieve more accurate Medical Physics, Vol. 42, No. 9, September 2015

5250

rendition of the anatomy, and to obtain slice images with improved resolution and reduced artifacts, thus allowing the 3D image data to be viewed in views other than the PA view. The dual-view DTS technique may find use in applications of DTS imaging where better separation of the fulcrum objects from the off-fulcrum objects, and more accurate 3D rendition of the anatomy and/or display with views other than PA view is desired.

ACKNOWLEDGMENTS This work was supported in part by Grant Nos. CA104759, CA124585, and CA138502 from the NIH-NCI, a research Grant No. EB00117 from the NIH-NIBIB, and a subcontract from NIST-ATP. This research was also supported in part by NIH through MD Anderson’s Cancer Center Support Grant No. CA016672.

a)Author

to whom correspondence should be addressed. Electronic mail: [email protected] 1Y. Yamada et al., “Fast scanning tomosynthesis for the detection of pulmonary nodules diagnostic performance compared with chest radiography, using multidetector-row computed tomography as the reference,” Invest. Radiol. 46(8), 471–477 (2011). 2A. Tingberg et al., “Breast cancer screening with tomosynthesis-initial experiences,” Radiat. Prot. Dosim. 147(1-2), 180–183 (2011). 3A. Tingberg, “X-ray tomosynthesis: A review of its use for breast and chest imaging,” Radiat. Prot. Dosim. 139(1-3), 100–107 (2010). 4D. J. Godfrey, H. P. McAdams, and J. T. Dobbins, “Stochastic noise characteristics in matrix inversion tomosynthesis (MITS),” Med. Phys. 36(5), 1521–1532 (2009). 5A. S. Chawla et al., “Design and development of a new multi-projection xray system for chest imaging,” IEEE Trans. Nucl. Sci. 56(1), 36–45 (2009). 6J. T. Dobbins et al., “Digital tomosynthesis of the chest for lung nodule detection: Interim sensitivity results from an ongoing NIH-sponsored trial,” Med. Phys. 35(6), 2554–2557 (2008). 7J. T. Dobbins et al., “Digital tomosynthesis of the chest,” J. Thorac. Imaging 23(2), 86–92 (2008). 8J. Vikgren et al., “Comparison of chest tomosynthesis and chest radiography for detection of pulmonary nodules: Human observer study of clinical cases,” Radiology 249(3), 1034–1041 (2008). 9A. Terzi et al., “Lung cancer detection with digital chest tomosynthesis: Baseline results from the observational study SOS,” J. Thorac. Oncol. 8(6), 685–692 (2013). 10B. Zhao et al., “Experimental validation of a three-dimensional linear system model for breast tomosynthesis,” Med. Phys. 36(1), 240–251 (2009). 11J. A. Gersh et al., “Improved volumetric imaging in tomosynthesis using combined multiaxial sweeps,” J. Appl. Clin. Med. Phys. 11(4), 181–195 (2010). 12T. T. Littleton, “A phantom method to evaluate the clinical effectiveness of a tomographic device,” Am. J. Roentgenol. 108(4), 847–856 (1970). 13J. T. Dobbins and D. J. Godfrey, “Digital x-ray tomosynthesis: Current state of the art and clinical potential,” Phys. Med. Biol. 48(19), R65–R106 (2003). 14A. C. Kak and M. Slaney, Principles of Computerized Tomographic Imaging (IEEE, New York, NY, 1988). 15L. A. Feldkamp, L. C. Davis, and J. W. Kress, “Practical cone-beam algorithm,” J. Opt. Soc. Am. A 1(6), 612–619 (1984). 16A. P. Dempster, N. M. Laird, and D. B. Rubin, “Maximum likelihood from incomplete data via the EM algorithm,” J. R. Soc. Stat. B39, 1–38 (1977). 17L. A. Shepp and Y. Vardi, “Maximum likelihood reconstruction for emission tomography,” IEEE Trans. Med. Imaging 1(2), 113–122 (1982). 18K. Lange, M. Bahn, and R. Little, “A theoretical study of some maximum likelihood algorithms for emission and transmission tomography,” IEEE Trans. Med. Imaging 6(2), 106–114 (1987). 19J. Kay, “The EM algorithm in medical imaging,” Stat. Methods Med. Res. 6, 55–75 (1997).

5251

Zhong et al.: Dual-view DTS imaging

20T. M. Benson and J. Gregor, “Framework for iterative cone-beam micro-CT

reconstruction,” IEEE Trans. Nucl. Sci. 52(5), 1335–1340 (2005). 21T. M. Buzug, Computed Tomography: From Photon Statistics to Modern Cone-Beam CT (Springer-Verlag, Berlin, Heidelberg, 2008). 22C. A. Bouman and K. Sauer, “A generalized gaussian image model for edgepreserving MAP estimation,” IEEE Trans. Image Process. 2(3), 296–310 (1993). 23R. L. Siddon, “Fast calculation of the exact radiological path for a 3dimensional CT array,” Med. Phys. 12(2), 252–255 (1985). 24B. De Man and S. Basu, “Distance-driven projection and backprojection in three dimensions,” Phys. Med. Biol. 49(11), 2463–2475 (2004).

Medical Physics, Vol. 42, No. 9, September 2015

5251 25Y.

Long, J. A. Fessler, and J. M. Balter, “3D forward and back-projection for x-ray CT using separable footprints,” IEEE Trans. Med. Imaging 29(11), 1839–1850 (2010). 26A. S. Chawla et al., “Optimized image acquisition for breast tomosynthesis in projection and reconstruction space,” Med. Phys. 36(11), 4859–4869 (2009). 27M. J. Murphy et al., “The management of imaging dose during imageguided radiotherapy: Report of the AAPM Task Group 75,” Med. Phys. 34(10), 4041–4063 (2007). 28J. Shan et al., “Stationary chest tomosynthesis using a carbon nanotube x-ray source array: A feasibility study,” Phys. Med. Biol. 60(1), 81–100 (2015).