A New Weighting Scheme for Fan-beam and Circle Cone-beam CT Reconstructions
Wei Wang, Xiang-Gen Xia, Chuanjiang He, Zemin Ren, Jian Lu, Tianfu Wang, Baiying Lei
JJOURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 1
A New Weighting Scheme for Fan-beam and CircleCone-beam CT Reconstructions
Wei Wang, Xiang-Gen Xia, Chuanjiang He, Zemin Ren, Jian Lu, Tianfu Wang and Baiying Lei
Abstract —In this paper, we first present an arc based algorithmfor fan-beam computed tomography (CT) reconstruction viaapplying Katsevich’s helical CT formula to 2D fan-beam CTreconstruction. Then, we propose a new weighting functionto deal with the redundant projection data. By extending theweighted arc based fan-beam algorithm to circle cone-beamgeometry, we also obtain a new FDK-similar algorithm forcircle cone-beam CT reconstruction. Experiments show that ourmethods can obtain higher PSNR and SSIM compared to theParker-weighted conventional fan-beam algorithm and the FDKalgorithm for super-short-scan trajectories.
Index Terms —FBP algorithm, Parker’s weight, super-short-scan, fan-beam CT, circle cone-beam CT
I. I
NTRODUCTION C OMPUTED tomography (CT) has been widely used inclinical diagnosis and industrial applications since itsability of providing inner vision of an object without destruct-ing it. In classical tomography, the fan-beam CT reconstructionalgorithm is fundamental since it can be heuristically extendedto 3D helical cone-beam [1] [2] [3] [4] [5] [6] and circlecone-beam [7] [8] [9] [10] [11] [12] CT reconstructions. Thestandard fan-beam reconstruction method is the ramp filter-based filtered backprojection (FBP) algorithm [13] [14], whichcan be derived from the Radon inversion formula. For 2DCT image reconstruction, to exactly and stably reconstruct thewhole image, it requires to measure all the line integrals of theX-rays that diverge from all directions and pass through theobject. To make a fan-beam CT that samples data on a circulartrajectory meet this condition, the detector must be largeenough to cover the fan-angle of ± γ m = arcsin( R m /R o ) This work was supported partly by National Natural ScienceFoundation of China (Nos.12001381, 61871274 and 61801305),China Postdoctoral Science Foundation (2018M64081), PeacockPlan (No. KQTD2016053112051497), Shenzhen Key Basic ResearchProject (Nos. JCYJ20180507184647636, JCYJ20170412104656685,JCYJ20170818094109846, and JCYJ20190808155618806).Wei Wang, Tianfu Wang, and Baiying Lei are with the School of Biomed-ical Engineering, Shenzhen University, National-Regional Key TechnologyEngineering Laboratory for Medical Ultrasound, Guangdong Key Laboratoryfor Biomedical Measurements and Ultrasound Imaging, School of BiomedicalEngineering, Health Science Center, Shenzhen University, Shenzhen, China.(e-mail: [email protected], [email protected], [email protected]).Xiang-Gen Xia is with the Department of Electrical and ComputerEngineering, University of Delaware, Newark, DE 19716, USA. (e-mail:[email protected]).Chuanjiang He is with the College of Mathematics and Statistics,Chongqing University, Chongqing, China (e-mail: [email protected]).Zemin Ren is with the College of Mathematics and Physics, ChongqingUniversity of Science and Technology, Chongqing, China (e-mail: [email protected]).Jian Lu is with the Shenzhen Key Laboratory of Advanced MachineLearning and Applications, Shenzhen University, Shenzhen, China (e-mail:[email protected]). to avoid truncated projections and the X-ray source musttravel on a continuous arc of π + 2 γ m on the circle toensure that the line integrals of the X-rays diverging from alldirections in the 2D plane and passing through the object aremeasured, where R o is the radius of the scanning trajectoryand R m < R o is the radius of the object. In the literature,the range of π + 2 γ m is called as short-scan. In [15], Parkerproposed a weighting function to weight the short-scan fan-beam projection data before convoluting it to avoid completingthe projection data of the remain angles ([ π + 2 γ m : π ])and implement the reconstruction algorithm efficiently. Whenthe range of scanning angles is larger than the short-scan,redundant projection data is measured. In [16], Silver extendedParker’s weighting function by utilizing the virtual detectorsto tackle these redundant projection data.In 2002, Noo et al. [17] proposed another type of FBPalgorithm for fan-beam CT reconstruction by decomposing theconvolution of the ramp filter into a successive convolutionsof a Hilbert filter and a derivative filter. After that, many otheralgorithms [9] [18] [19] based on the Hilbert transform wereproposed. One advantage of these new algorithms is that theycan exactly reconstruct a part of the image even though therange of the scanning angles is less than the short scan. Thesenew algorithms can be seen as special cases that applying the3D exact helical cone-beam inversion formula [20] [21] to the2D CT reconstructions. In [22], You et al. derived a Hilberttransform based FBP algorithm for fan-beam full- and partial-scans, in which the backprojection does not include position-dependent weights. In [9] and [17], to deal with the redundantdata, a continuous weight function was also proposed to weightthe filtered sinograms.A specific view-dependent data differentiation was a com-mon processing step in these Hilbert transform based al-gorithms. In [3] [23], the implementation of this step wasresearched in order to improve the resolution and qualityof the reconstructed image. In [4], Zamyatin et al. used theTaylor series expansions to approximate the derivative of theHilbert transform and proposed new algorithms for fan-beamand helical cone-beam CT reconstruction.The fan-beam algorithms can be heuristically extended to3D cone-beam CT reconstructions [1] [2] [3] [5] [6] [7][9] [11] [12] [24], which are referred to as approximatingalgorithms. On the other hand, there exist a lot of exact cone-beam algorithms [10] [20] [21] [25] [26] [27] in the literature.The advantages of the approximating algorithms are that theyare easy to implement and flexible to modulate the ramp kernelfor different clinical applications, while those of the exactalgorithms are that they can reconstruct images with good a r X i v : . [ ee ss . I V ] J a n OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 2 resolution even when the cone angle is large. In [2] [6] [10][11] [12] [24] [27], the weighting functions are also proposedto tackle the redundant cone-beam projection data.In this paper, we first present an arc based algorithm forfan-beam CT reconstruction by applying Katsevich’s helicalCT [17] formula to 2D fan-beam CT reconstruction. Then,we propose a new weighting function to tackle the redundantprojection data. By extending the arc based algorithm to thecircle cone-beam geometry, we also obtain a new algorithm forcircle cone-beam CT reconstruction. Our weighting functionis different from the ones used in [9] [15] [17]. The weightingfunctions in [9] [15] [17] depend on the rotation angle of theX-ray source and the diverging direction of the X-ray, andare required to be continuous with respect to these argumentswhile our weighting function depends on the rotation angleof the X-ray source and the positions of the pixels of thereconstructed image, and has no continuity constraint. Thus,the CT images reconstructed by our method may have higherresolutions and less artifacts. Moreover, there exists a hyper-parameter in the weighting functions of [9] [17] that greatlyinfluences the performance of the algorithm while our weight-ing function has no hyper-parameter.The rest of the paper is organized as follows. In Sec-tion II we briefly introduce the related works. In SectionIII, we present an arc based algorithm for fan-beam CTreconstruction and derive a weighting function for tacklingthe redundant data, and then extend the algorithm with theweighting function to circular cone-beam CT reconstruction.Numerical evaluation is presented in Section IV, followed byconclusion in Section V.II. R
ELATED WORKS
In this section, we briefly describe some works related toour method.
A. Noo’s fan-beam formula
In [17], Noo et al. decomposed the convolution of the rampfilter into a successive convolutions of the Hilbert filter andthe derivative filter in the Fourier domain by observing that | σ | = (1 / π )( i πσ )( − i sign ( σ )) , (1)where | σ | , i πσ and − i sign ( σ ) are the Fourier transformof the ramp filter, the derivative filter and the Hilbert filter,respectively. By utilizing this relationship, they proposed aFBP reconstruction formula for fan-beam CT.Let g ( λ, γ ) = g ( λ, a ( γ )) (2)be the measured projections, where λ is a polar angle of theX-ray source, a ( λ ) = ( R o cos λ, R o sin λ ) is a position ofthe X-ray source, R o is the radius of the trajectory circleof the X-ray source, γ = cos γe + sin γe is a divergingdirection of the X-ray, γ ∈ [ − γ m , γ m ] , γ m is the half fan-angle, e = − (cos λ, sin λ ) , and e = ( − sin λ, cos λ ) . Thenthe reconstruction formula for fan-beam CT with equi-angulardetetcor reads: f ( x ) = 12 π (cid:90) Λ d λ (cid:107) x − a ( λ ) (cid:107) [ w ( λ, φ ) g F ( λ, φ )] φ = φ ∗ ( λ,x ) , (3) where Λ is the X-ray source trajectory satisfying the datacompleteness condition [17], φ ∗ ( λ, x ) = arctan x · e R o + x · e , | φ ∗ ( λ, x ) | < π/ (4)is the angle characterizing the ray that diverges from a ( λ ) andpasses x , g F ( λ, φ ) = (cid:90) γ m − γ m d γh H (sin( φ − γ )) Å ∂∂λ + ∂∂γ ã g ( λ, γ ) , (5) h H ( s ) = − (cid:90) + ∞−∞ d σ i sgn( σ )e i2 πσs = 1 πs (6)is the Hilbert filter and w ( λ, φ ) = c ( λ ) c ( λ ) + c ( λ + π − φ ) (7)is the weighting function with c ( λ ) = cos π ( λ − λ s − d )2 d if λ s < λ < λ s + d, if λ s + d < λ < λ e − d, cos π ( λ − λ e + d )2 d if λ e − d < λ < λ e (8)and d is an angular interval over which c ( λ ) smoothly dropsfrom 1 to 0. In the experiments, they set d = 10 π/ .Applying the changes of variables u = D tan γ and ˜ u = D tan φ in equations (3) and (5), they obtain thereconstruction formula for fan-beam CT with equally spacedcollinear detectors: f ( x ) = 12 π (cid:90) Λ d λ R o + x · e [ w ( λ, ˜ u ) g F ( λ, ˜ u )] ˜ u =˜ a ∗ ( λ,x ) (9)with g F ( λ, ˜ u ) = (cid:90) u m − u m d uh H (˜ u − u ) D √ D + u Å ∂∂λ + D + u D ∂∂u ã g ( λ, u ) , (10)where u m = D tan γ m , ˜ u ∗ ( λ, x ) = D tan φ ∗ ( λ, x ) = Dx · e R o + x · e (11)is the detector location of the line connecting x to a ( λ ) and w ( λ, ˜ u ) = c ( λ ) c ( λ ) + c ( λ + π − u/D )) . (12) B. Katsevich’s helical cone-beam formula
In [21], Katsevich proposed an exact reconstruction formulafor helical cone-beam CT and Noo et al. [28] researched howto efficiently and accurately implement it. The reconstructionformula can be written as f ( x ) = − π (cid:90) λ o ( x ) λ i ( x ) d λ (cid:107) x − a ( λ ) (cid:107) g F Å λ, x − a ( λ ) (cid:107) x − a ( λ ) (cid:107) ã (13)where λ i ( x ) and λ o ( x ) are the extremities of the π -line passingthrough x with λ i ( x ) < λ o ( x ) , g F ( λ, θ ) = (cid:90) π d γh H (sin γ ) g (cid:48) ( λ, cos γθ +sin γ ( θ × m ( λ, θ ))) , (14) OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 3 g (cid:48) ( λ, θ ) = lim ε → g ( λ + ε, θ ) − g ( λ, θ ) ε , (15) h H ( s ) is the Hilbert filter defined in equation (6), g ( λ, θ ) isthe measured projection data and vector m ( λ, θ ) is normal tothe κ -plane K ( λ, ψ ) of the smallest | ψ | value that contains theline of direction θ through a ( λ ) .III. P ROPOSED M ETHOD
A. Arc based fan-beam algorithm
Applying equation (13) to the 2D fan-beam reconstruction,the π -line becomes a chord of the circle of the scanningtrajectory, the κ -plane K ( λ, ψ ) coincides with the image planeto reconstruct and so θ × m ( λ, θ ) = θ ⊥ .Let a ( λ ) = ( R o cos λ, R o sin λ ) be the position of the X-raysource on the trajectory circle of radius R o and g ( λ, θ ) = (cid:90) + ∞ d tf ( a ( λ ) + tθ ) (16)be the measured projection data, where θ ∈ S is a divergingdirection of the X-ray and S is the unit circle in the 2D plane.Then, according to equation (13), the reconstruction formulafor fan-beam CT can be written as f ( x ) = − π (cid:90) chord ( x ) d λ (cid:107) x − a ( λ ) (cid:107) g F Å λ, x − a ( λ ) (cid:107) x − a ( λ ) (cid:107) ã , (17)where | x | < R m < R o , R m is the radius of the test object, g F ( λ, θ ) = (cid:90) π d¯ γh H (sin ¯ γ ) g (cid:48) ( λ, cos ¯ γθ + sin ¯ γθ ⊥ ) , (18) h H ( s ) and g (cid:48) ( λ, θ ) are defined by equations (6) and (15),respectively, chord ( x ) is a chord passing through x anddividing the trajectory circle into two arcs, and the integrald λ is done on any one of the two arcs.Fig. 1: An example of the X-ray source trajectory. Theprojection data measured along the three red arcs can be usedto reconstruct the intensity of x .From equation (17), we can argue that the intensity of anypoint x with | x | < R m < R o in the circle can be reconstructedif there exists a chord passing through x such that all theprojections g ( λ, θ ) on the corresponding arcs of the chordpassing through the neighborhood of x are measured. Notethat g ( λ, θ ) can be measured from another direction since g ( λ, θ ) = g ( λ , − θ ) , where θ = a ( λ ) − a ( λ ) | a ( λ ) − a ( λ ) | . See Fig. 1 forexample.For any fixed x in the trajectory circle, there may existmany available chords passing through x and so the integral arcs can be different. Moreover, when the range of scanningarc is larger than the short scan, redundant projection dataneeds to be processed. A simple method to tackle theseredundant projection data is to calculate equation (17) for allavailable chords and average them. However, averaging themis equivalent to filtering the reconstructed image by a low-passfilter and so may reduce the resolution of the reconstructedimage. In this paper we only consider two types of chords thatrespectively pass through the two endpoints of the scanningarcs.Therefore, we propose the following formula for fan-beamreconstruction: f ( x ) = − π (cid:90) λ P λ d λ (cid:36) ( x, λ ) (cid:107) x − a ( λ ) (cid:107) g F Å λ, x − a ( λ ) (cid:107) x − a ( λ ) (cid:107) ã , (19)where λ < λ P correspond to the two endpoints of thescanning arcs, (cid:36) ( x, λ ) = 12 ( (cid:36) ( x, λ ) + (cid:36) ( x, λ )) (20)is a weighting function, and (cid:36) ( x, λ ) and (cid:36) ( x, λ ) are definedby (cid:36) ( x, λ ) = ® , if a ( λ ) ∈ ˇ(cid:0) a ( λ ) , a ( λ ( x, λ )) , else , (21) (cid:36) ( x, λ ) = ® , if a ( λ ) ∈ ˇ(cid:0) a ( λ ( x, λ P )) , a ( λ P ) , else , (22)where a ( λ ( x, λ )) is another endpoint of the chord passingthrough a ( λ ) and x , a ( λ ( x, λ P )) is another endpoint of thechord passing through a ( λ P ) and x , and ˇ(cid:0) a ( λ ) , a ( λ ( x, λ )) denotes the arc starting from a ( λ ) and ending at a ( λ ( x, λ )) .
1) Implementation for the curved-line detector:
In thissubsection, we describe how to implement equation (19) whenthe fan-beam projections are measured by using an equi-angular curved-line detector.Fig. 2: Geometry of data acquisition by using an equi-angularcurved detector.Let a ( λ ) = ( R o cos λ, R o sin λ ) be the position of the X-ray source, θ ( λ, γ ) = cos γe ( λ )+sin γe ( λ ) be the divergingdirection, g c ( λ, γ ) = g ( λ, θ ( λ, γ )) be the measured projectiondata using the curved-line detector, where γ is the samplingcoordinate of the fan-angle, e ( λ ) = ( − sin( λ ) , cos( λ )) , e ( λ ) = ( − cos( λ ) , − sin( λ )) (See Fig. 2). By the chain rule[17], equation (15) can be implemented as g ( λ, γ ) := g (cid:48) ( λ, θ ( λ, γ )) = ∂g c ( λ, γ ) ∂λ + ∂g c ( λ, γ ) ∂γ . (23) OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 4
Applying the changes of variable ¯ γ = γ (cid:48) − γ , equation (18)can be rewritten as g ( λ, γ ) := g F ( λ, θ ( λ, γ ))= − (cid:90) π d γ (cid:48) h H (sin( γ − γ (cid:48) )) g ( λ, γ (cid:48) ) . (24)Let γ ∗ ( x, λ ) = arctan Ä x · e ( λ ) R o + x · e ( λ ) ä be the fan angle of themeasured X-ray that diverges from a ( λ ) and passes through x . Then, equation (19) can be rewritten as f ( x ) = − π (cid:90) λ P λ d λ (cid:36) ( x, λ ) (cid:107) x − a ( λ ) (cid:107) g ( λ, γ ∗ ) . (25)Note that the main differences of equation (25) and equation(3) are the weighting functions w ( λ, φ ) and w ( x, λ ) .
2) Implementation for straight-line detector:
In this sub-section, we show how to implement the algorithm for the fan-beam projections measured by using an equi-space straight-line detector.Fig. 3: Geometry of data acquisition by using an equi-spaceline detector.Let a ( λ ) = ( R o cos λ, R o sin λ ) be the position of the X-raysource, θ ( λ, u ) = ue ( λ )+ De ( λ ) √ u + D be the diverging direction, g l ( λ, u ) = g ( λ, θ ( λ, u )) be the measured projection data usingthe straight-line detector, where u is the sampling coordinateparallel to the detector, D is the distance between the X-raysource and the detector, e ( λ ) = ( − sin( λ ) , cos( λ )) , e ( λ ) =( − cos( λ ) , − sin( λ )) (See Fig. 3). By the chain rule, equation(15) can be implemented as g ( λ, u ) := g (cid:48) ( λ, θ ( λ, u )) = ∂g l ( λ, u ) ∂λ + u + D D ∂g l ( λ, u ) ∂u . (26)Applying the changes of variables u = D tan( γ ) and ¯ u = D tan( γ (cid:48) ) , equations (24) and (25) can be respectively rewrit-ten as g ( λ, u ) := g F ( λ, θ ( λ, u ))= − (cid:90) u m − u m d¯ u D √ ¯ u + D h H ( u − ¯ u ) g ( λ, ¯ u ) (27)and f ( x ) = − π (cid:90) λ P λ d λ (cid:36) ( x, λ ) R o + x · e ( λ ) g ( λ, u ∗ ) , (28)where u m = D tan γ m , γ m = arcsin( R m /R o ) and u ∗ ( x, λ ) = D x · e ( λ ) R o + x · e ( λ ) is the detector location of themeasured X-ray that diverges from a ( λ ) and passes through x . B. Circle cone-beam algorithm
The fan-beam algorithm (19) can be extended for circlecone-beam CT reconstruction heuristically. In the following,we give the detailed circle cone-beam reconstruction algo-rithms for flat-plane and curve-plane detectors, respectively.Fig. 4: Geometry of data acquisition by using a flat-planedetector.
1) Implementation for the flat-plane detector:
Let a ( λ ) =( R cos( λ ) , R cos( λ ) , be the position of the X-ray source, ( u, v, w ) be the local detector coordinates with unit vectors e u ( λ ) = ( − sin( λ ) , cos( λ ) , ,e v ( λ ) = ( − cos( λ ) , − sin( λ ) , ,e w ( λ ) = (0 , , , (29)and g f ( λ, u, w ) = g ( λ, θ ( λ, u, w )) be the measured projectiondata by using a flat-plane detector with diverging direction θ ( λ, u, w ) = ue u ( λ ) + De v ( λ ) + we w ( λ ) √ u + D + w , (30)where D is the distance between the detector plane and theX-ray source (See Fig. 4). To make formula (19) effective forcircle cone-beam CT reconstruction, we need to modify θ ⊥ inequation (18) and the weighting function (cid:36) ( x, λ ) in equation(20). In this paper, we heuristically set θ ⊥ ( λ, u, w ) = θ ( λ, u, w ) × ( De v ( λ ) + we w ( λ ) √ D + w × e u ( λ )) (31)and (cid:36) d ( x, λ ) = (cid:36) ( Proj ( x ) , λ ) , (32)where Proj ( x ) represents projecting x on the circle planeformed by the trajectory of the X-ray source. Then, as derivedin Appendix A, the circle cone-beam reconstruction algorithmfor the flat-plane detector can be implemented as follows. • step1: derivative at constant direction g ( λ, u, w ) := g (cid:48) ( λ, θ ( λ, u, w ))= g f ( λ, u, w ) ∂λ + u + D D g f ( λ, u, w ) ∂u + uwD g f ( λ, u, w ) ∂w . (33) • step2: convolution with Hilbert filter: g ( λ, u, w ) := g F ( λ, θ ( λ, u, w ))= (cid:90) u m − u m d¯ u D √ ¯ u + D + w h H ( u − ¯ u ) g ( λ, ¯ u, w ) , (34)where u m = D tan γ m , γ m = arcsin( R m /R o ) . OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 5 • step3: backprojection f ( x ) = 12 π (cid:90) λ P λ d λ (cid:36) d ( x, λ ) v ∗ ( x, λ ) g ( λ, u ∗ , w ∗ ) , (35)where ( u ∗ ( x, λ ) , w ∗ ( x, λ )) is the position in the detectorof the measured X-ray that diverges from a ( λ ) and passesthrough x , which can be calculated by v ∗ ( x, λ ) = R o + x · e v ( λ ) ,u ∗ ( x, λ ) = Dv ∗ ( x, λ ) ( x · e u ( λ )) ,w ∗ ( x, λ ) = Dv ∗ ( x, λ ) ( x · e w ( λ )) . (36)Fig. 5: Geometry of data acquisition by using a curve-planedetector.
2) Implementation for the curve-plane detector:
Thecurved-plane detector array consists of N row × N cols detectors.The detector columns are perpendicular to the trajectory circleof the X-ray source, while the detector rows form circle arcsparallel to each other and to the trajectory circle, where thecenter of the arc on the trajectory circle plane coincides withthe X-ray source. Let a ( λ ) = ( R cos( λ ) , R cos( λ ) , be theposition of the X-ray source, ( α, v, w ) be the local detectorcoordinates with unit vectors defined in equation (29) and g c ( λ, α, w ) = g ( λ, θ ( λ, α, w )) be the measured projection databy using a curve-plane detector with diverging direction θ ( λ, α, w ) = D sin( α ) e u ( λ ) + D cos( α ) e v ( λ ) + we w ( λ ) √ D + w , (37)where D is the distance between the detector plane and the X-ray source (See Fig. 5). The curve-plane detector coordinates ( α, v c , w c ) can be be converted to the flat detector coordinates ( u, v f , w f ) via u = D tan( α ) , v f = v c , w f = w c cos( α ) . (38)Applying the changes of variable in equation (38), we canobtain the circle cone-beam reconstruction algorithm for thecurve-plane detector as follows. • step1: derivative at constant direction g ( λ, α, w ) : = g (cid:48) ( λ, θ ( λ, α, w ))= g c ( λ, α, w ) ∂λ + g c ( λ, α, w ) ∂α . (39) • step2: convolution with Hilbert filter: g ( λ, α, w ) := g F ( λ, θ ( λ, α, w ))= (cid:90) π d¯ α D √ ¯ D + w h H (sin( α − ¯ α )) g ( λ, ¯ α, w ) . (40) • step3: backprojection f ( x ) = 12 π (cid:90) λ P λ d λ (cid:36) d ( x, λ ) v ∗ ( x, λ ) cos( α ∗ ) g ( λ, α ∗ , w ∗ ) , (41)where ( α ∗ ( x, λ ) , w ∗ ( x, λ )) is the position in the detectorof the measured X-ray that diverges from a ( λ ) and passesthrough x , which can be calculated by v ∗ ( x, λ ) = R o + x · e v ( λ ) ,α ∗ ( x, λ ) = arctan ( x · e u ( λ ) v ∗ ( x, λ ) ) ,w ∗ ( x, λ ) = D cos( α ∗ ) v ∗ ( x, λ ) ( x · e w ( λ )) . (42) C. Additional discrete schemes
In order to implement our CT reconstruction algorithms, weneed to give the discrete definitions of the derivatives ∂g c ( λ,γ ) ∂λ and ∂g c ( λ,γ ) ∂γ in equation (23), ∂g l ( λ,u ) ∂λ and ∂g l ( λ,u ) ∂u in equation(26), g f ( λ,u,w ) ∂λ , g f ( λ,u,w ) ∂u and g f ( λ,u,w ) ∂w in equation (33), and g c ( λ,α,w ) ∂λ and g c ( λ,α,w ) ∂α in equation (39), the Hilbert filters h H ( u ) in equations (27) and (34), h H (sin( γ )) in equation (24)and h H (sin( α )) in equation (40), and the weighting functions (cid:36) ( x, λ ) in equation (21) and (cid:36) ( x, λ ) in equation (22).Let γ i = i ∇ γ, i = − N, − N + 1 , ..., , ..., N − , N,α i = i ∇ α, i = − N, − N + 1 , ..., , ..., N − , N,u j = j ∇ u, j = − M, − M + 1 , ..., , ..., M − , M,w k = k ∇ w, k = − L, − L + 1 , ..., , ..., L − , L (43)be the discrete coordinates of the detectors and λ s = s ∇ λ, s = 0 , , ..., P (44)be the discrete sampling angles. To avoid the space shift in thereconstructed images, we use the centered difference schemes OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 6 (a) Original (b) CFA (c) ACE (d) Ours
Fig. 6: CT images reconstructed from short-scan fan-beam projection data.
200 210 220 230 240 250 260 270 280 290 300020406080100120140160
OriginalCFAACEOur
Fig. 7: 1D intensity profile passing through the red line in Fig.6a.to calculate the derivatives: ∂g c ( λ s , γ ) ∂λ = g c ( λ s +1 , γ ) − g c ( λ s − , γ )2 ∇ λ ,∂g c ( λ, γ i ) ∂γ = g c ( λ, γ i +1 ) − g c ( λ, γ i − )2 ∇ γ ,∂g l ( λ s , u ) ∂λ = g l ( λ s +1 , u ) − g l ( λ s − , u )2 ∇ λ ,∂g l ( λ, u j ) ∂λ = g l ( λ, u j +1 ) − g l ( λ, u j − )2 ∇ u ,g f ( λ s , u, w ) ∂λ = g f ( λ s +1 , u, w ) − g f ( λ s − , u, w )2 ∇ λ ,g f ( λ, u j , w ) ∂u = g f ( λ, u j +1 , w ) − g f ( λ, u j − , w )2 ∇ u ,g f ( λ, u, w k ) ∂u = g f ( λ, u, w k +1 ) − g f ( λ, u, w k − )2 ∇ w ,g c ( λ s , α, w ) ∂λ = g c ( λ s +1 , α, w ) − g c ( λ s − , α, w )2 ∇ λ ,g c ( λ, α i , w ) ∂α = g c ( λ, α i +1 , w ) − g c ( λ, α i − , w )2 ∇ α . (45) Let b denote a cut-off frequency for the Hilbert filter h H ( u ) .Then, we may write [29]: k H ( u ) ≈ − (cid:90) b − b i sgn( σ ) e i πσu dσ = (cid:90) − b ie i πσu dσ − (cid:90) b ie i πσu dσ = ï πu e i πσt ò σ = − b − ï πu e i πσt ò bσ =0 = 12 πu (1 − e − i πbu − e i πbu + 1)= 12 πu (2 − πbu ))= 1 πu (1 − cos(2 πbu )) . (46)For h H (sin( γ )) and h H (sin( γ )) , we have h H (sin( γ )) = γ sin( γ ) h H ( γ ) = 1 − cos(2 πbγ ) π sin( γ ) (47)and so h H (sin( α )) = 1 − cos(2 πbα ) π sin( α ) . (48)In this paper, we set b = 1 / (2 ∇ u ) , / (2 ∇ γ ) and / (2 ∇ α ) for h H ( u ) , h H (sin( γ )) and h H (sin( α )) , respectively. Therefore,the discrete definitions of the Hilbert filters are: k H ( u j ) = 1 − cos( πu j / ∇ u ) πu j ,h H (sin( γ i )) = 1 − cos( πγ i / ∇ γ ) π sin( γ i ) ,h H (sin( α i )) = 1 − cos( πα i / ∇ α ) π sin( α i ) . (49)For the weighting functions (cid:36) ( x, λ ) and (cid:36) ( x, λ ) , sincethe endpoints a ( λ ( x, λ )) and a ( λ ( x, λ P )) may not coincidewith the sampling points a ( λ s ) , we need to interpolate them at OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 7 (a) Original (b) CFA (c) ACE (d) Ours
Fig. 8: CT images reconstructed from super-short-scan fan-beam projection data.the endpoints of the sampling angles λ ( x, λ ) and λ ( x, λ P ) .Let s ( x, λ ) =( λ ( x, λ ) − λ ) / ∇ λ, (cid:98) s ( x, λ ) (cid:99) = floor ( s ( x, λ )) , (cid:100) s ( x, λ ) (cid:101) = floor ( s ( x, λ )) + 1 ,s ( x, λ P ) =( λ ( x, λ P ) − λ ) / ∇ λ, (cid:98) s ( x, λ P ) (cid:99) = floor ( s ( x, λ P )) , (cid:100) s ( x, λ P ) (cid:101) = floor ( s ( x, λ P )) + 1 . (50)Then, we set (cid:36) ( x, λ s ) = , if s ∈ [0 , (cid:98) s ( x, λ ) (cid:99) ] ,s − (cid:98) s ( x, λ ) (cid:99) , if s ∈ ( (cid:98) s ( x, λ ) (cid:99) , (cid:100) s ( x, λ ) (cid:101) ) , , else ,(cid:36) ( x, λ s ) = , if s ∈ [ (cid:100) s ( x, λ P ) (cid:101) , P ] ,s − (cid:98) s ( x, λ P ) (cid:99) , if s ∈ ( (cid:98) s ( x, λ P ) (cid:99) , (cid:100) s ( x, λ P ) (cid:101) ) , , else , (51)where λ ( x, λ ) and λ ( x, λ P ) can be calculated via AppendixB. IV. E XPERIMENTAL R ESULTS
In this section, we give some simulation results to verifythe effectiveness of our algorithms. For 2D fan-beam CTreconstructions, we do the experiments on the simulatedprojection data measured by equi-angular curved-line detectorswhile for 3D circle cone-beam, the simualted projection datameasured by flat-plane detectors are used to reconstruct CTimages. The codes for implementing our methods can bedownloaded from https://github.com/wangwei-cmd/arc-based-CT-reconstruction. We also provide the codes that implementthe compared methods coded by ourselves in this paper.
A. Fan-beam with curved-line detectors
In this subsection, we present some CT images recon-structed from the fan-beam projection data measured by equi-angular curved-line detectors under short-scan and super-short-scan trajectories, and compare the results with thoseof the conventional fan-beam algorithm with Parker-extendedweighting function (CFA) [16] and Noo’s algorithm (ACE)[17]. To test the performances of our method and the comparedalgorithms, we randomly choose 500 full dose CT images (ofsize × ) from “the 2016 NIH-AAPM-Mayo Clinic LowDose CT Grand Challenge” [30] as the original images. Theparameters for the fan-beam CT with the equi-angular curved-line detector are set as follows: R o = 500 , R m = √ × and so γ m = arctan( R m /R o ) = 35 . ◦ . We set γ = [ −
36 : 0 . × π/ for the samplingpositions on the γ coordinate. For the sampling positions onthe λ coordinate, we set λ = [0 : 1 : 180 + 2 × × π/ for short-scan and λ = [0 : 1 : 180] × π/ for super-short scan. The hyper-parameter d in ACE [17] isset as d = 6 × π/ .In Fig. 6, we show one set of CT images reconstructed byCFA, ACE and ours from the short-scan fan-beam projectiondata. We can observe that the visual effects of the reconstructedimages by ACE and ours are very similar while that by CFAhas a lot of artifacts. To objectively estimate the qualities ofthe reconstructed images by the three methods, the averagepeak signal to noise ratio (PSNR) and structural similarity(SSIM) of the 500 reconstructed CT images are listed in TableI. We can see that our method has a slightly higher averagePSNR and SSIM than those of ACE. The average PSNR andSSIM of CFA are far lower than those of ACE and ours. It’sbecause that the CFA algorithm needs higher sampling rate onthe λ coordinate. By experiments, we find that when setting λ = [0 : 0 .
25 : 180 + 2 × × π/ , the PSNR and SSIMof CFA can be raised to the same level of ACE and ours.In Fig. 7, we plot the 1D line intensity profile passingthrough the red line in Fig. 6a. We can observe that theintensity lines of ACE and ours are almost coincident andresemble more closely to the one of the original compared tothat of CFA.In Fig. 8, a set of CT images reconstructed by CFA,ACE and ours from the super-short-scan fan-beam projectionsare presented. We can observe that all the three images OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 8 (a) β =[0:210] (b) β =[0:220] (c) β =[0:230] (d) β =[0:240] Fig. 9: CT images reconstructed from fan-beam projections with different β . The images in the first row are reconstructed byACE and in the second row by ours. (a) var=1 (b) var=10 (c) var=100 (d) var=200 Fig. 10: CT images reconstructed from noisy short-scan fan-beam projection data with different variances.reconstructed by CFA, ACE and ours suffer from intensityinhomogeneity on the left part of the reconstructed imagesdue to the data incompleteness. Moreover, the stripe visualeffect can be easily observed in Fig. 8b. In Fig. 8c, a verticalline can be observed. This may be caused by data missingat the endpoints of the scanning arcs. Compared to CFA andACE, our method has the best visual effect as can be seen fromFig. 8d. In Table II, the average PSNR and SSIM of the CTimages reconstructed by CFA, ACE and ours are, respectively,listed, from which we can see that our method has the highestaverage PSNR and SSIM.To better demonstrate the artifacts of ACE caused by dataincompleteness, the CT images corresponding to λ = [0 : 1 :210] × π/ , λ = [0 : 1 : 220] × π/ , λ = [0 : 1 :230] × π/ and λ = [0 : 1 : 240] × π/ reconstructedby ACE and ours are shown in Fig. 9. We can observe that when the range of λ is lower than ◦ , some stripe visualeffects can be observed in the CT images reconstructed byACE. As the range of λ exceeds ◦ , the artifacts of ACEcaused by the data missing almost disappear. For our method,even λ = [0 : 1 : 210] × π/ , there exists no stripe artifactin the reconstructed CT images.To test the performances of our method when the projectiondata is corrupted by noise, we add Poisson and Gaussian noiseto the short-scan projections via the following formulas [31]: g = exp( − g/M ) ,g = g + I ∗ poission ( g ) + I ∗ Gaussin ( m, var/I ) ,g = log( I /g ) ∗ M, (52)where M is the maximal value of the projection data g , m and var are the mean and variance of the Gaussian noise,respectively, I is the average photon count for the Poisson OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 9 (a) Original (b) FDK (c) ACE (d) Ours
Fig. 11: Sliced CT images reconstructed from short-scan circle cone-beam projection data.noise. In this experiment, we set I = 10 , m = 0 , and var =1 , , , , respectively. The reconstructed images fromthe noisy projection data are shown in Fig. 10. We can observethat as the variance var increases, the reconstructed CT imagehas more noise. B. Circle cone-beam with flat-plane detectors
In this subsection, we give some CT images reconstructedfrom the circle cone-beam projection data measured by flat-plane detectors to very the effectiveness of our method, andcompare the results with those of the Feldkamp-Davis-Kress TABLE I: The Averaged PSNR and SSIM of CT ImagesReconstructed by CFA, ACE and Ours from Fan-Beam Short-Scan Projection Data.
PSNR SSIMCFA 18.74 ± ± ± ± ± ± (FDK) algorithm [7] and Noo’s algorithm (ACE) [17], wherewe extend the weighting function of ACE such that it can be OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 10 (a) Original (b) FDK (c) ACE (d) Ours
Fig. 12: Sliced CT images reconstructed from super-short-scan circle cone-beam projection data.TABLE II: The Averaged PSNR and SSIM of CT ImagesReconstructed by CFA, ACE and Ours from Fan Beam Super-Short-Scan Projection Data.
PSNR SSIMCFA 18.51 ± ± ± ± ± ± used to reconstruct circle cone-beam CT images.To test the performances of our method and the compared algorithms, we randomly choose 2500 full dose CT images (ofsize × ) from “the 2016 NIH-AAPM-Mayo Clinic LowDose CT Grand Challenge” [30] and use them to assemble 50objects of size × × as the original images.The parameters for the circle cone-beam CT with the flat-plane detector are set as follows: u = [ −
494 : 1 : 494] , v =[ −
54 : 1 : 54] , R o = 1000 , and D = ceil ( R o + √ × . For the sampling positions on the λ coordinate, we set λ = [0 : 1 : 180 + 2 × × π/ OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 11
TABLE III: The Averaged PSNR and SSIM of CT ImagesReconstructed by FDK, ACE and Ours from Circle Cone-Beam Short-Scan Projection Data.
PSNR SSIMCFA 29.90 ± ± ± ± ± ± TABLE IV: The Averaged PSNR and SSIM of CT ImagesReconstructed by FDK, ACE and Ours from Circle Cone-Beam Super-Short-Scan Projection Data.
PSNR SSIMCFA 24.69 ± ± ± ± ± ± for short-scan and λ = [0 : 1 : 180] × π/ for super-short-scan. The th slice of the images of the objectis assumed to lie in the plane z = 0 (i.e. the plane formed bythe trajectory of the X-ray source). The hyper-parameter d inACE [17] is set as d = 10 × π/ .In Fig. 11, we present some sliced CT images of one objectreconstructed by FDK, ACE and ours from the short-scancircle cone-beam projection data. We can observe that thereexist some stripe artifacts in the CT images reconstructedby FDK. The visual effects of the CT images reconstructedby ACE and ours are very similar. We also use the PSNRand SSIM to measure the similarities of the reconstructed CTimages and the original. From Table III, we can observe thatthe average PSNR and SSIM of our method are slightly higherthan that of ACE and are higher than that of CFA, whichcoincides with our observations.In Fig. 12, some sliced CT images of one object recon-structed by FDK, ACE and ours from the super-short-scancircle cone-beam projection data are shown. From Fig. 12b,we can see that some stripe visual artifacts exit in the CTimages reconstructed by FDK. Moreover, the left hand sideparts of the CT images reconstructed by FDK suffer fromsevere intensity inhomogeneity. From Fig. 12b, we can observethat there also exit some undesirable vertical lines in the CTimages reconstructed by ACE, which may be caused by thedata incompleteness. Compared to Fig. 12b and Fig. 12c,the CT images reconstructed by our method suffer from lessintensity inhomogeneity and their visual effects are the best asshown in Fig. 12d. The PSNR and SSIM are used to evaluatethe qualities of the reconstructed images and listed in TableIV. We can see that our method has the highest PSNR andSSIM compared to CFA and ACE.V. C ONCLUSION
In this paper, we proposed a new weighting function todeal with the redundant projection data for the arc basedfan-beam CT algorithm, which was obtained via applyingKatsevich’s helical CT formula [21] to 2D fan-beam CT reconstruction. By extending the arc based algorithm to circlecone-beam geometry with the proposed weighting function,we also obtained a new FDK-similar algorithm for circlecone-beam CT reconstruction. Experiments showed that ourmethods can obtained higher PSNR and SSIM compared tothe related algorithms when the scanning trajectories are super-short-scan. A
PPENDIX AA LGORITHM FOR CIRCLE CONE - BEAM CT WITHFLAT - PLANE DETECTOR
Fig. 13: Relation between γ and the flat-plane detector coor-dinates.As can be observed from Fig. 13, we have e u ( λ ) = −−−→ A A ||−−−→ A A || , e v ( λ ) = −−−→ A A −−−→ A A , e w ( λ ) = −−−→ A A −−−→ A A , ||−−−→ A A || = u, ||−−−→ A A || = ¯ u, ||−−−→ A A || = w, ||−−−→ A A || = D, (cid:54) A A A = α, (cid:54) A A A = ¯ γ. (53)We also have tan( α + ¯ γ ) = ||−−−→ A A ||||−−−→ A A || = ¯ u √ D + w and so d¯ γ = cos ( α + ¯ γ ) √ D + w d¯ u = √ D + w D + w + ¯ u d¯ u. (54)By the Law of Sines, we have sin ¯ γ = ||−−−→ A A || sin( (cid:54) A A A ) ||−−−→ A A || = √ D + w (¯ u − u ) √ D + w + ¯ u √ D + w + u . (55)Therefore, h H (sin ¯ γ ) = √ D + w + ¯ u √ D + w + u √ D + w h H (¯ u − u ) . (56)Note that θ ( λ, u, w ) = −−−→ A A ||−−−→ A A || and so cos ¯ γθ + sin ¯ γθ ⊥ = −−−→ A A ||−−−→ A A || . Substituting equations (54) and (56) into equation(18), we can obtain g F ( λ,θ ( λ, u, w )) = − √ u + D + w D × (cid:90) u m − u m d¯ u D √ ¯ u + D + w h H ( u − ¯ u ) g ( λ, ¯ u, w ) , (57) OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 12
From Fig. 4, we can observe that (cid:107) x − a ( λ ) (cid:107) = » ( u ∗ ) + ( w ∗ ) + D ( x − a ( λ )) · e v ( λ ) D = (cid:112) ( u ∗ ) + ( w ∗ ) + D D ( R o + x · e v ( λ )) . (58)Substituting equation (58) into equation (17), we can get f ( x ) = − D (cid:112) D + ( w ∗ ) + ( u ∗ ) × π (cid:90) λ P λ d λ (cid:36) d ( x, λ ) v ∗ ( x, λ ) g ( λ, u ∗ , w ∗ ) . (59)Canceling the factors − √ u + D + w D and − D √ D +( w ∗ ) +( u ∗ ) in equations (57) and (59), we obtainthe algorithm for circle cone-beam CT with the flat-planedetector. A PPENDIX BA LGORITHM FOR CALCULATING ONE ENDPOINT OF ACHORD
Let ( R o cos λ , R o sin λ ) and ( R o cos λ , R o sin λ ) be thetwo end-points of a chord on a circle with radius r = R o and ( x , x ) be a point on the chord. Then, we have R o (cos λ − cos λ ) = t ∗ ( x − R o cos λ ) ,R o (sin λ − sin λ ) = t ∗ ( x − R o sin λ ) , (60)where t ∈ [0 , . Solving equation set (60), we have cos λ = { R o ( − R o + x + x ) cos λ − x ( − R o + x R o cos λ + x R o sin λ ) } / { R o ( R o + x + x − x R o cos λ − x R o sin λ ) } , sin λ = { R o ( − R o + x + x ) sin λ − x ( − R o + x R o cos λ + x R o sin λ ) } / { R o ( R o + x + x − x R o cos λ − x R o sin λ ) } . (61)Thus, we can get λ by λ = arccos(cos λ ) , where λ needsto be changed by λ = 2 π − λ when sin λ < .R EFERENCES[1] G. Wang, T. H. Lin, P. C. Cheng, and D. M. Shinozaki, “A GeneralCone-Beam Reconstruction Algorithm,”
IEEE Transactions on MedicalImaging , vol. 12, no. 3, pp. 486–496, 1993.[2] K. Stierstorfer, A. Rauscher, J. Boese, H. Bruder, S. Schaller, andT. Flohr, “Weighted FBP - a Simple Approximate 3D FBP Algorithmfor Multislice Spiral CT With Good Dose Usage for Arbitrary Pitch,”
Physics in Medicine and Biology , vol. 49, no. 11, pp. 2209–2218, JUN7 2004.[3] F. Noo, S. Hoppe, F. Dennerlein, G. Lauritsch, and J. Hornegger, “ANew Scheme for View-Dependent Data Differentiation in Fan-Beam andCone-Beam Computed Tomography,”
Physics in Medicine and Biology ,vol. 52, no. 17, pp. 5393–5414, 2007.[4] A. A. Zamyatin, K. Taguchi, and M. D. Silver, “Practical Hybrid Con-volution Algorithm for Helical CT Reconstruction,”
IEEE Transactionson Nuclear Science , vol. 53, no. 1, pp. 167–174, 2006.[5] H. Kudo, T. Rodet, F. Noo, and M. Defrise, “Exact and ApproximateAlgorithms for Helical Cone-Beam CT,”
Physics in Medicine andBiology , vol. 49, no. 13, pp. 2913–2931, 2004.[6] X. Y. Tang, J. Hsieh, R. A. Nilsen, S. Dutta, D. Samsonov, andA. Hagiwara, “A Three-Dimensional-Weighted Cone Beam FilteredBackprojection (CB-FBP) Algorithm for Image Reconstruction in Volu-metric CT - Helical Scanning,”
Physics in Medicine and Biology , vol. 51,no. 4, pp. 855–874, 2006. [7] L. FELDKAMP, L. DAVIS, and J. KRESS, “Practical Cone-BeamAlgorithm,”
Journal of the Optical Society of America a-Optics ImageScience and Vision , vol. 1, no. 6, pp. 612–619, 1984.[8] L. F. Yu, X. C. Pan, and C. A. Pelizzari, “Image Reconstruction Witha Shift-Variant Filtration in Circular Cone-Beam CT,”
InternationalJournal of Imaging Systems and Technology , vol. 14, no. 5, pp. 213–221,2004.[9] H. Kudo, F. Noo, M. Defrise, and R. Clackdoyle, “New Super-Short-Scan Algorithms for Fan-Beam and Cone-Beam Reconstruction,” in
IEEE Nuclear Science Symposium and Medical Imaging Conference ,2003, pp. 902–906.[10] S. Tang and X. Tang, “Axial Cone-Beam Reconstruction by WeightedBPF/DBPF and Orthogonal Butterfly Filtering,”
IEEE Transactions onBiomedical Engineering , vol. 63, no. 9, pp. 1895–1903, 2016.[11] S. Tang, K. Huang, Y. Cheng, T. Niu, and X. Tang, “Three-DimensionalWeighting in Cone Beam FBP Reconstruction and Its Transforma-tion Over Geometries,”
IEEE Transactions on Biomedical Engineering ,vol. 65, no. 6, pp. 1235–1244, JUN 2018.[12] X. Y. Tang, J. Hsieh, A. Hagiwara, R. A. Nilsen, J. B. Thibault,and E. Drapkin, “A Three-Dimensional Weighted Cone Beam FilteredBackprojection (CB-FBP) Algorithm for Image Reconstruction in Vol-umetric CT Under a Circular Source Trajectory,”
Physics in Medicineand Biology , vol. 50, no. 16, pp. 3889–3905, 2005.[13] B. K. P. Horn, “Fan-Beam Reconstruction Methods,”
Proceedings of theIEEE , vol. 67, no. 12, pp. 1616–1623, 1979.[14] G. T. Herman and A. Naparstek, “Fast Image-Reconstruction Based ona Radon Inversion Formula Appropriate for Rapidly Collected Data,”
SIAM Journal on Applied Mathematics , vol. 33, no. 3, pp. 511–533,1977.[15] D. PARKER, “Optimal Short Scan Convolution Reconstruction forFanbeam CT,”
Medical Physics , vol. 9, no. 2, pp. 254–257, 1982.[16] M. D. Silver, “A Method for Including Redundant Data in ComputedTomography,”
Medical Physics , vol. 27, no. 4, pp. 773–774, 2000.[17] F. Noo, M. Defrise, R. Clackdoyle, and H. Kudo, “Image ReconstructionFrom Fan-Beam Projections on Less Than a Short Scan,”
Physics inMedicine and Biology , vol. 47, no. 14, pp. 2525–2546, JUL 21 2002.[18] G. H. Chen, “A New Framework of Image Reconstruction From FanBeam Projections,”
Medical Physics , vol. 30, no. 6, pp. 1151–1161,2003.[19] F. Noo, R. Clackdoyle, and J. Pack, “A Two-Step Hilbert TransformMethod for 2D Image Reconstruction,”
Physics in Medicine and Biology ,vol. 49, no. 17, pp. 3903–3923, SEP 7 2004.[20] Y. Zou and X. Pan, “Exact Image Reconstruction on PI-Lines FromMinimum Data in Helical Cone-Beam CT,”
Physics in Medicine andBiology , vol. 49, no. 6, pp. 941–959, MAR 21 2004.[21] A. Katsevich, “An Improved Exact Filtered Backprojection Algorithmfor Spiral Computed Tomography,”
Advances in Applied Mathematics ,vol. 32, no. 4, pp. 681–697, MAY 2004.[22] J. You and G. L. Zeng, “Hilbert Transform Based FBP Algorithm forFan-Beam CT Full and Partial Scans,”
IEEE Transactions on MedicalImaging , vol. 26, no. 2, pp. 190–199, 2007.[23] A. Katsevich, “A Note on Computing the Derivative at a ConstantDirection,”
Physics in Medicine and Biology , vol. 56, no. 4, pp. N53–N61, 2011.[24] K. Taguchi, B. S. S. Chiang, and M. D. Silver, “A New WeightingScheme for Cone-Beam Helical CT to Reduce the Image Noise,”
Physicsin Medicine and Biology , vol. 49, no. 11, pp. 2351–2364, 2004.[25] A. Katsevich, “Theoretically Exact Filtered Backprojection-Type Inver-sion Algorithm For Spiral CT,”
SIAM Journal on Applied Mathematics ,vol. 62, no. 6, pp. 2012–2026, AUG 2 2002.[26] L. Yu, D. Xia, Y. Zou, E. Y. Sidky, J. Bian, and X. Pan, “A RebinnedBackprojection-Filtration Algorithm for Image Reconstruction in HelicalCone-Beam CT,”
Physics in Medicine and Biology , vol. 52, no. 18, pp.5497–5508, 2007.[27] H. Schoendube, K. Stierstorfer, and F. Noo, “Accurate Helical Cone-Beam CT Reconstruction With Redundant Data,”
Physics in Medicineand Biology , vol. 54, no. 15, pp. 4625–4644, 2009.[28] F. Noo, J. Pack, and D. Heuscher, “Exact Helical Reconstruction UsingNative Cone-Beam Geometries,”
Physics in Medicine and Biology ,vol. 48, no. 23, pp. 3787–3818, 2003.[29] A. J. Wunderlich, “The Katsevich Inversion Formula for Cone-BeamComputed Tomography,” Ph.D. dissertation, Oregon State University,2006.[30] C. McCollough, “TU-FG-207A-04: Overview of the Low Dose CTGrand Challenge,”
Medical Physics , vol. 43, no. 3759-3760, pp. 3759–3760, 2016.
OURNAL OF L A TEX CLASS FILES, VOL. 14, NO. 8, AUGUST 2015 13 [31] Y. Liu, J. Ma, Y. Fan, and Z. Liang, “Adaptive-Weighted Total VariationMinimization for Sparse Data Toward Low-Dose X-Ray Computed Tomography Image Reconstruction,”