Fast Scatter Artifacts Correction for Cone-Beam CT without System Modification and Repeat Scan

1.2 MB·English

by  Wei Zhao

1
Fast Scatter Artifacts Correction for Cone-Beam CT
without System Modification and Repeat Scan
Wei Zhao, Jun Zhu, Luyao Wang
Abstract—We provide a fast and accurate scatter artifacts applied to yield the final scatter estimation. By simply sub-
5
correction algorithm for cone beam CT (CBCT) imaging. The tractingthefinalscatterestimatefromtherawprojectiondata,
1
method starts with an estimation of coarse scatter profile for
0 we can obtain the corrected projection data. The corrected
a set of CBCT images. A total-variation denoising algorithm
2 image was then reconstructed using the corrected projection
designed specifically for Poisson signal is then applied to derive
n the final scatter distribution. Qualitatively and quantitatively data. MC simulations, experimental phantom data and in vivo
a evaluations using Monte Carlo (MC) simulations, experimental human data were used to validate and evaluate the proposed
J CBCT phantom data, and in vivo human data acquired for a scatter correction method.
9 clinicalimageguidedradiationtherapywereperformed.Results
1 showthattheproposedalgorithmcansignificantlyreducescatter II. MATERIALANDMETHODS
artifacts and recover the correct HU within either projection
] domainorimagedomain.Furthertestshowsthemethodisrobust The measured raw projection data can be modeled as the
h with respect to segmentation procedure. summation of the primary and scatter signals. A general
p scatter artifact correction strategy is to estimate the scattered
Index Terms—Cone beam computed tomography (CBCT),
-
radiation distribution and then subtract it from the measured
d scatter correction, denoising, polychromatic reprojection, spec-
e trum estimation. raw projection data, i.e.
m
I (α,(cid:126)x)=I(α,(cid:126)x)−I (α,(cid:126)x), (1)
p s
.
s I. INTRODUCTION whereI isthecorrectedprojectiondata,namely,theestimated
c p
si DUE to the increased axial coverage of multi-slice com- primary projection data, I is the raw data and Is is the scatter
y puted tomography (CT) and the introduction of flat- signal. Indices α and (cid:126)x stand for projection view angle and
h panel detector, the size of X-ray illumination fields has grown detector channel number respectively. For simplicity, we will
p dramatically,causinganincreaseinscatterradiation.Formost drop the index (cid:126)x in the following discussions.
[
of cases, artifacts are present when the model used in the

1 reconstruction algorithm is not consistent with the projection A. Coarse scatter estimation
v data acquisition model. Existing reconstruction algorithms Theproposedmethod,asshowninFig.1,startswithanim-
9 don’t model the scatter radiation, so scatter artifacts appear age segmentation of the scatter artifacts contaminated CBCT
0
in the CT images. Typical scatter artifacts show as shading, image volume, followed by a polychromatic reprojection of
4
4 streaks between high contrast objects, reduced contrast reso- the segmented image volume using a given x-ray spectrum
0 lution and inaccurate Hounsfield Units (HU) numbers. Scatter and known data acquisition geometry. The reprojected data is
1. correction has been extensively studied in the past serval thensubtractedfromthemeasuredrawprojectiondataateach
0 decades but a clinically sensible solution remains illusive. givenviewangletogenerateacoarseestimateofthenecessary
5 Current scatter correction methods can be roughly classified scatter profile used for denoising algorithm.
1 into five approaches: physical scatter rejection, analytical Based on the segmented image volume, the primary signal
:
v modeling, Monte Carlo (MC) simulation, primary modulation can be modeled as follows:
Xi and scatter measurements [1]. (cid:90) Emax (cid:34) (cid:90) l (cid:35)
To develop a scatter artifacts correction method that is Iˆ =N dEΩ(E)η(E)exp − µ(E,s)ds , (2)
p
r
fast, accurate, and, ideally, has no need for extra scans, extra 0 0
a
measurements or extra hardware support or modifications is whereN isthetotalnumberofphotons,Ω(E)isthepolychro-
thegoalofthisstudy.Weprovideafastandaccuratescatterar- maticX-rayspectrum,η(E)istheenergy-dependentefficiency
tifactscorrectionmethodforconebeamCT(CBCT)imaging. of the detector, which is simply considered as proportional
In the proposed method, the difference between the measured to photon energy E because most CT scanners use energy-
raw projection data and a polychromatic reprojection of a integrating detectors. E is the maximum photon energy
max
segmentedimagevolumeisusedtogenerateacoarseestimate of the spectrum. µ(E,s) is the energy-dependent linear at-
of scatter profile. Based on the coarse scatter estimation, a tenuation coefficient and l is the propagation path length for
denoising algorithm specifically designed for Poisson signal eachray,andcanbecalculatedusingaGPU-basedray-tracing
which employs a total variation regularization term was then algorithm [2], [3]. By using these notations, the flood field I
0
can be modeled as follows:
W.Zhao(e-mail:[email protected]),J.ZhuandLY.Wangarewiththe
DepartmentofBiomedicalEngineering,HuazhongUniversityofScienceand (cid:90) Emax
I =N dEΩ(E)η(E). (3)
Technology,Hubei,China. 0
0
2
yieldasmoothscatterdistribution[6].Thedenoisingalgorithm
Uncorrected Image Raw Data Denoised Scatter Corrected Data
which is based on the total-variation regularization is aimed
to solve the following optimization problem,
(cid:90) β (cid:90)
I =argmin d(cid:126)x(I −IˆlogI )+ d(cid:126)x|∇I |2. (6)
s Is s s s 2 s
Multi-segmentation Subtract
Segmented Image Denoising The first term of (6) is a data-fidelity term that is designed
Reprojection Data Coarse Scatter CCoorrrercetecdt eimd aIgme age specifically for Poisson statistical signal and keep I close
s
to the data Iˆ, while the second term is a total-variation
s
regularization term to keep the solution I smooth, namely,
s
being dominated with low frequency content. β is a constant
Normalized flux00000.....0000012345Spectrum tovoabrjdieaectttieiovrnemalifnuaenpcpthtrioeoanrcehlias[t6ivc].oenivmexporatnadncecaonf btheestowlvoedterumssin.gThains
00 En50ergy (ke1V00) 150
It has to note that the proposed method can also be imple-
Fig.1. Flowchartoftheproposedscattercorrectionmethod.Themethodcan mented in image-domain. For those cases when we only have
starts with either the raw projection data or an scatter contaminated volume
images or volume, we need to estimate scatter artifact error
as input. After an initial CT image reconstruction, the image is segmented
into different components, based on which a polychromatic reprojection is imageswhichcanbeaddeddirectlytotheuncorrectedimages.
performed with predetermined spectrum. The reprojection data is subtracted To achieve this goal, a scatter projection error ∆p which
s
from the raw projection data to yield a coarse scatter estimate. The coarse
can be added linearly in the logarithmic raw-data domain can
scatter is then applied to a denoising algorithm which was specifically
designedforPoissonsignaltoyieldthefinalscatterestimate.Thecorrected be calculated. Because the tomographic reconstruction is a
projection data is obtained by subtracting the final scatter from the raw linear process and the order of summation and backprojection
projectiondataandthecorrectedCTimagecanbegeneratedusingastandard
operations are interchangeable.
filteredbackprojectionreconstructionalgorithm.
C. Monte Carlo simulations
AfterreplacingN withI andsubstitutingitinto(2),theresult
0
describes the estimated primary signal with flood field I : To validate the proposed algorithm, an anthropomorphic1
0
thorax phantom and an abdomen phantom were used to
(cid:104) (cid:105)
Iˆp = I0(cid:82)0EmaxdE(cid:82)ΩE(mEax)dη(EEΩ)(eExp)η−(E(cid:82))0lµ(E,s)ds . (4) glaetnioenratpeacMkaCgesimGAulTaEtio[n7]d.aTtaowqiuthanatitGateivanelty4-banasdedquMalCitastiivmeuly-
0 investigate the effect of scatter artifacts reduction, both of
Subtracting the above estimated primary signal from the the phantoms include a water insert at the central area and
measured total projection data I, the coarse scatter signal Iˆ the water insert has three small low contrast inserts (adipose,
s
can be estimated as follows: breast and liver) and a bone insert.
(cid:104) (cid:105) Primary projection data, scatter only projection data, and
I (cid:82)EmaxdEΩ(E)η(E)exp −(cid:82)lµ(E,s)ds
Iˆ =I− 0 0 0 . primary plus scatter projection data were extracted indepen-
s (cid:82)EmaxdEΩ(E)η(E) dently.Forafastersimulation,aparallelgeometryandaplane
0 (5) X-ray source (2D 320 × 120 mm2 rectangle source) were
Note that the segmented image volume is generated from used in the MC simulation. The distance from the source
a scatter contaminated reconstruction and the polychromatic to the center of rotation is 750 mm and the distance from
x-ray spectrum used in the reprojection can be recovered the detector to center of rotation is 450 mm. A circular scan
using an indirect transmission measurement-based spectrum was simulated and a total of 360 projections per rotation are
estimation method [4]. acquired over an angular range of 360◦. The polychromatic
X-ray source spectrum is 125 kVp and it is generated using
B. Denoising the coarse scatter Spektr software [8] with 5 mm aluminum filtration. For each
Since the coarse scatter estimate Iˆ is highly dependent on of the simulations, a total number of 3×1010 photons were
s
emitted.
the segmentation procedure, it may yield inaccurate results,
In a first experiment, we investigate scatter correction for
especiallyforlowcontrastobjectsthathavesimilarattenuation
the anthropomorphic thorax phantom and abdomen phantom
properties as the background and for the edges of two neigh-
in both projection domain and image domain. The primary
boring materials. To compensate for the inaccuracy caused by
imagewasusedasreferenceimagesandthetotalimageswere
theinaccuratesegmentation,insteadofregularizingthecoarse
applied to the proposed scatter correction algorithm.
scatter using a convolution-based scatter model [5], in this
Totesttherobustnessofmethodwithrespecttothesegmen-
study, we directly denoise the coarse scatter estimate using a
tation, two components (lung and tissue), three components
statistical-based denoising algorithm.
(lung, tissue and bone) and four components (lung, water,
As we know, scatter radiation distribution is typically
tissue and bone) segmentation were performed based on the
predominantly of low frequency content not only in spatial
uncorrected thorax phantom images respectively.
domain but also in temporal domain. A denoising technique
specifically designed for Poisson signal was applied to Iˆs to 1http://www.qrm.de/content/pdf/QRM-Thorax.pdf.
3
True Scatter Coarse Scatter Denoised Scatter Thorax Phantom Water Insert Abdomen Phantom Water Insert
100 Breast
motnahP elgnA weiV 8600 egamI yramirP Ad2ip oseB 4o1 n e 3L iver
x
aro 40
h
T
20 eg
am
Position 0 I latoT
50
motnahP nemodb 432000 detcerroC nointciaemjoDPor
A
10
Fig. 2. Sinograms of true scatter, coarse scatter, and denoised scatter fo0r detcerro niamoD eg
an anthropomorphic thorax phantom and an abdomen phantom. The coarse Cam
I
scatterwhichisthedifferenceoftherawprojectiondataandthepolychromatic
reprojection data have global scatter profiles but also contain errors which
were partially generated from inaccurate segmentation. The denoised coarse
Fig. 3. Results of the scatter correction for the MC simulation data of the
scattersmooththeseerrorsandfitthetruescatterprofileswell.
thorax phantom and the abdomen phantom within both projection domain
andimagedomain.(Firstrow)Resultsfortheprimarysignalandtheycanbe
D. Physical phantom experiments and patient study regarded as ground truth. (Second row) Results for the primary plus scatter
signal(totalsignal)andscatterinducedshadingartifactsandstreaksareclearly
The proposed algorithm was also evaluated using experi- seen. (Third row) Results for the scatter corrected images within projection
data domain. (Fourth row) Results for the scatter corrected images within
mental data of a Catphan500 phantom (The Phantom Labo-
image domain. Shading artifacts and streaks are significantly reduced after
ratory, Salem, NY) scanned using a CBCT on-board imager correctionwithinbothprojectiondomainandimagedomain.Displaywindow:
(Varian 2100EX System, Varian Medical Systems, Palo Alto, [-1200HU,500HU]forthoraxphantomimagesand[-300HU,300HU]forthe
waterinsertimagesandtheabdomenphantom.
CA).Atotalof678projectionswereevenlyacquiredina360
degree rotation with 2×2 binning. Both wide collimation and
(a) (b)
narrow collimation modes were applied with the same scan 800 Primary 800 Primary
parameters. 600 Total 600 Total
A set of in vivo pelvis data was also used to evaluate the Corrected (PD) Corrected(PD)
400 Corrected(ID) 400 Corrected(ID)
proposed method. This data was acquired in half-fan mode
using state-of-the-art of the kV CBCT imaging system of UH 200 UH 200
the Varian TrueBeam system(Varian Medical Systems, Palo 0 0
Alto, USA). The tube potential was set to 125 kVp and -200 -200
678 projections were acquired on a flat-panel detector (Varian -400 -400
Breast Adipsoe Liver Bone Breast Adipose Liver Bone
4030CB imager) with 1024 columns and 768 rows.
Fig. 4. Results of reconstruction values in HU numbers at different ROIs
III. RESULTS ofthewaterinsertof(a)thethoraxphantomand(b)theabdomenphantom
beforeandaftercorrectionwithinbothprojectiondomainandimagedomain.
A. Monte Carlo simulations TheHUnumbersaresuccessfullyrecoveredafterscattercorrectionusingthe
proposedmethod.PDandIDstandforprojectiondomainandimagedomain,
Fig. 2 shows the sinograms of the MC reference scatter, respectively.
coarsescatterandthedenoisedscatterfortheanthropomorphic
thorax phantom and the abdomen phantom. Fig. 3 shows the B. Experimental phantom and patient study
results of the scatter correction using the proposed method
Fig. 6 shows the CT images and the corresponding differ-
for the MC simulation data of the thorax phantom and the
ence images of the Catphan(cid:13)R 500 phantom with and without
abdomen phantom within both projection domain and image
scatter correction. Fig. 7 shows CT images of an in vivo
domain.Theprimaryimageswerereconstructedusingprimary
pelvis scan with the kilovoltage CBCT imaging system of
projections and they were regarded as scatter-free images.
the Varian TrueBeam system (Varian Medical Systems, Palo
Fig. 4 depicts the HU numbers of breast, adipose, liver and
Alto, USA) without and with scatter correction. After scatter
boneinserts(showninFig.3)ofprimaryimages,totalimages
correction,theHUnumberofthedarkregionwassuccessfully
and scatter corrected images.
recoveredandthemissinganatomicalstructuresarerecovered.
The influence of different segmentation methods on the
It took about 1 min to correct a typical Varian clinical dataset
accuracy of the scatter correction for the thorax phantom is
(512×512×81) using a NVIDIA GeForce GTX 480 card.
depicted in Fig. 5.
4
Segmented Image Thorax Phantom Water Insert Thorax Phantom Water Insert Uncorrected Corrected
(PD) (PD) (ID) (ID)
detcerroC stnenopmoC Lung WBoTantiesesr u e
4
detcerroC stnenopmoC Lung BoTnises ue laixA
3
Lung Tissue
detcerroC stnenopmoC
2
Fig. 5. The influence of the segmentation methods on the accuracy of the m
scattercorrectionisshownforthethoraxphantomwithwaterinsertinboth o
projection domain (PD) and image domain (ID). (First row) Results for the oz
four components segmentation with lung, water, tissue and bone. (Second la
row) Results for the three components segmentation with lung, tissue and ix
A
bone. (Third row) Results for the two components segmentation with only
lungandtissue.Displaywindow:[-1200HU,500HU]forthethoraxphantom
imagesand[-300HU,300HU]forthewaterinsertimages.
Narrow Image Wide Image Corrected (PD) Corrected (ID) la
n
o
r
o
seg C
am
I la
ttig
ecn aS
erefeR
o
t ffiD Fofigt.h7e.VSacraiattnerTcrourerBecetaiomnfsoyrstaempe.lvAisxiCalBvCiTews,cacnoruosninalgvthieewkVanidmasagginitgtaslyvstieemw
aredepictedrespectively.Theuncorrectedimagesacquiredonthecommercial
systemshowshadingartifactsandpartoftheanatomicstructurewasmissing.
Fig. 6. Catphan(cid:13)R500 phantom with and without scatter correction. The After scatter correction, the shading artifacts were successfully reduced and
the missing anatomic structure can be clearly seen. Display window: [-200
narrowcollimationimageisconsideredtobethescatter-freeimage.Shading
HU,100HU]forallimages.
and reduced HU are shown in the wide collimation image.After correction,
scatter artifacts are greatly reduced and HU are successfully recovered for
bothprojectiondomain(PD)andimagedomainmethods(ID).Thedifference no modifications to existing system hardware, and can also
imagesshoweachimagesubtractedwiththenarrowcollimationimage.Note handle cross-scatter for dual source dual energy CT scanner.
that the wide collimation scan and the narrow collimation scan are two
Furthermore, the method can also be applied to spiral cone-
independent scans and the registration can not be perfect. Display window:
[-150HU,150HU]forboththeCBCTimagesandthedifferenceimages. beam CT where the anti-scatter grid can be uninstalled to
further reduce radiation dose.
IV. DISCUSSIONANDCONCLUSIONS
ACKNOWLEDGMENT
It has to note that the selection of the denoising parameters
The authors would like to thank Drs Jennifer Smilowitz,
depends on segmented image. Since the noise level of the
Guanghong Chen, Ke Li and Stephen Brunner for providing
MC simulation data is very high, a relatively large K and
the experimental data. The authors are also grateful to Dr Lei
β (compared to the patient data) were used to smooth the
Xing for his careful reading of the manuscript and his useful
segmentation error and to yield a denoised scatter that fit the
comments.
true scatter well enough. To our belief, most of the clinical
imageshavemuchlowernoiselevelandthesegmentedimages REFERENCES
would be superior to the segmented images of the MC data in
[1] L.Zhu,etal.,Medicalphysics,vol.36,p.2258,2009.
this study. Nevertheless, a larger K and β are recommended [2] G. Pratx and L. Xing, Medical physics, vol. 38, no. 5, pp. 2685–2697,
for badly segmented images to reinforce regularization. 2011.
[3] X.Jia,etal.,Phys.Med.Biol.,vol.59,no.4,p.R151,2014.
In summary, this work presents a novel strategy to estimate
[4] W.Zhao,etal.,Phys.Med.Biol.,vol.60,no.1,p.339,2015.
acoarsescatteredradiationprofile,basedonwhichadenoising [5] W. Zhao, et al., in SPIE Medical Imaging. International Society for
algorithm was applied to yield the final scatter estimation. we OpticsandPhotonics,2014,pp.903310–903310.
[6] X.Jia,etal.,Medicalphysics,vol.39,no.12,pp.7368–7378,2012.
demonstrated that a significant increase in image uniformity
[7] S.Jan,etal.,Phys.Med.Biol.,vol.56,no.4,p.881,2011.
and HU accuracy were achieved after correction. The correc- [8] J. Siewerdsen, et al., Medical physics, vol. 31, no. 11, pp. 3057–3067,
tion algorithm requires minimum computational cost, requires 2004.