Transcript ch6_Tomography

Objectives to understand:
1. Formation of tomographic sections using geometrical blurring.
Dependence of section thickness on angle of rotation.
2. Acquisition of CT data as projections. Mathematical description of
a projection.
3. Iterative reconstruction.
4. The Projection Theorem and associated mathematics.
Tomography objectives (continued):
5. Filtered Back Projection and associated mathematics.
6. CT exposure requirements.
7. Definition of CT image values.
8. Spiral CT and CT angiography.
In conventional radiographic projection imaging, the attenuation
effects from material situated throughout the x-ray path are
superimposed at the detector plane, as shown below.
Detector
Source
•
This can lead to considerable ambiguity, for example when a low
contrast object is superimposed on a dense or anatomically complex
object.
The object of tomography is to obtain images corresponding to the
attenuation distribution within a planar section through the body.
Another problem with conventional radiographic projection imaging
is the presence of x-ray scatter, which makes it difficult to make
quantitative assessment of attenuation values. The first of these
problems was addressed by conventional film tomography
sometimes called blurring tomography. Eventually, both were
addressed by computed tomography (CT).
First let us consider the blurring approach which was introduced
circa 1930 by Ziedses Des Plantes in Holland (in the same PhD
thesis which introduced film subtraction angiography--try to beat
that when you choose your thesis topic!)
The basic geometry for blurring tomography is shown in below.
Source motion
A
B
•
•
Detector motion
A B
A
B
A
A and B are objects at two different depths within the patient. By
coordinating the motion of the detector and the source, objects such
as A in the plane of rotation remain stationary on the detector.
Objects such as B, which are outside the plane of rotation, are
smeared out as the source and detector are rotated. The degree of
blurring increases with the angle of rotation. Therefore the
thickness of the in-focus plane decreases with increased rotation.
Several different types of blurring may be used. For example,
linear, circular or hypocycloidal motions are used and result in
different characteristic blurring artifacts. The choice of motion is
made by considering which type of artifact is less likely to mimic
the anatomy which may be under study.
Figure 3A below shows a frontal conventional radiograph of the head.
In Figure 3B, blurring tomography has been used to select a plane
through the maxillary sinuses. The arrows point to a mass in one of
the maxillary sinuses.
A
B
In computed tomography, x-rays pass in a direction approximately
transverse to the longitudinal axis of the body and pass only through
the slice of interest, thus overcoming the problem of overlapping
information from adjacent slices. The goal is to obtain an image
representing the two dimensional distribution of attenuation
coefficients at each picture element in the image. Various
mathematical schemes may be used to reconstruct the image from
the projection information.
An advantage of the slice geometry, analogous to that of the line
scanned x-ray systems discussed in the section on scatter, is that
scatter is greatly reduced. Additionally, in CT, electronic detectors
capable of quantitative recording of the transmitted information are
used. These two factors were essential in providing the
revolutionary contrast resolution made available by CT in the mid
1970’s.
Figure 4 illustrates the concept of projection data. At several angles
around the edge of the slice several ray paths, collectively called a
projection, are passed through the slice of interest. Labeled in
Figure 4 are the first projection P 1 and three of its rays and the ith
projection P i and its jth ray pij.
In Sir Godfrey Hounsfield’s original scanner, there were eighty rays
per projection and the image was represented by an 80 x 80 matrix.
One might think this was the result of some advance mathematical
optimization procedure. The fact is he was advancing his source
and detector using a screw mechanism which had eighty turns.
Figure 4
Hounsfield’s scanner, for which he earned the Nobel prize in
Medicine was called a rotate-translate scanner. Following rotation
to a new projection angle, a linear ray scan was done to collect all
the ray data for the projection.
Modern scanners use a continuously rotating gantry and arrays of
detectors so that at each projection angle, all ray transmission
values are simultaneously collected. There was an evolution of
source and detector configurations leading to several generations
of scanners.
Figure 5 below shows a third generation scanner, in which an array of
several hundred individual detectors are rotated along with the x-ray
tube in order to acquire the projection data.
Also shown is a fourth generation scanner in which only the x-ray
tube undergoes circular motion and the radiation is detected by a
stationary ring of detectors.
Finally, Figure 5 shows a more recently developed electron beam
scanner configuration in which an electron beam is magnetically
deflected along the surface of a cone and impinges on a circular
target. X-rays are transmitted through the patient and detected by a
slightly offset detector array. This type of scanner is particularly
well suited to cardiac applications where scan times on the order
of tens of milliseconds are desired.
In each case the x-ray beam is collimated to define the desired slice
and reference detectors just outside the slice are used to monitor
beam intensity fluctuations.
Use of the reference detector to normalize the projection data is
illustrated below in Figure 6.
Reference detector
IR
Iij
I0
r
The detected intensity along one ray in the projection is given by
(1)
where r is measured along the ray. The overall attenuation is due to
all attenuation elements along the ray.
The reference detector provides a signal proportional to the incident
intensity according to
The basic objective of CT is to solve for the attenuation values uij
at each point in the slice given the complete set of line integrals Pij.
We can get an approximate estimate of how many projections are
required by referring to Figure 7 below which shows rays incident
on a subject of radius D which is to be represented as an image
with pixel size d.
d
D
The number n of pixels per diameter, which is equal to the number
of required rays per projection is given by
(5)
The total number of pixels, which is equal to the number of
unknown attenuation coefficients in the image is then given by,
Therefore the number of projections with n rays per projection is
given by
This analysis ignores the fact that pixels far from the axis of
rotation are probed by a smaller number of rays than those near the
center. When a more careful analysis is done in terms of the
sampling of spatial frequencies, the required number is twice that
calculated.
It should also be noted that in present fan beam scanners the rays
are divergent and algorithms are required to obtain the information
provided by an equivalent set of parallel rays which lend
themselves more easily to mathematical analysis.
Figure 8, from Krestel-Imaging Systems for Medical Diagnosis,
shows two scans of a CT phantom, one obtained with 180
projections (A) and the other (B) with 720 projections. Clearly the
180 projection scan is inadequate while the 720 projection scan is
significantly better.
willi angles movie
Kalender et al
Kalender et al
Iterative Reconstruction
The early CT scanners used iterative image reconstruction schemes.
These are time consuming, especially with the 512 x 512 image
matrices used on modern scanners. Nevertheless, they are somewhat
interesting. The main idea will be illustrated below for the simple
case of a 2 x 2 matrix. Suppose the attenuation coefficient image is
as shown in Figure 9a.
5
7
12
6
6
6
2
8
4
4
1st estimate
6
4
6
4
10 10
5
6
7
2
11
9
6.5 5.5
6.5 5.5
4.5 3.5
4.5 3.5
2nd estimate
13
7
10
10
Even in the original Hounsfield scanner there were 6400 pixels. One
can imagine that this process would be time consuming even for a
reasonably fast computer using modern matrices of over 250,000
pixels.
We will now examine the reconstruction of the attenuation
coefficient image using two dimensional Fourier transform
methods.
A plane wave is a two dimensional function of the form
(8)
characterized by sinusoidal variation in the direction of the
propagation vector
, where
This is illustrated in Figure 10 for the case of a wave
We have already seen that a one dimensional intensity distribution
can be represented as a Fourier integral representing a continuous
summation over sinusoidal functions of various frequencies. The
extension of this idea is that a two dimensional function, such as an
attenuation coefficient distribution (x,y) can be represented as a
sum of a continuous distribution of two dimensional sinusoidal
plane wave functions as
where
is the two dimensional Fourier transform of
u(x,y) representing the weightings of the plane waves having spatial
frequency coordinate (kx,ky)
If it were possible by means of x-ray projection data to obtain the
weighting coefficients
, the attenuation distribution
could be found from equation 11. This is made possible by a
mathematical relationship called the Projection Theorem.
Consider a coordinate system with the x’ direction perpendicular to
the direction of the x-ray projection as shown in Figure 11.
For a given projection anglejthe detector array measures a
projection consisting of detected rays
(13)
where the y’ integral is over all attenuation contributions in the
patient.
Suppose that we now take the one dimensional Fourier transform of
these detected values, obtaining
where the integrals cover the entire area of the slice.
If we compare equation 14 with equation 12, we obtain,
(15)
Therefore by performing the projection measurement and
transforming the detected values, we obtain the Fourier expansion
coefficients along the kx’ axis in frequency space with the kx’ axis
making an angle relative to the kx axis.
The number of samples in k space will be equal to the number of
elements in the detector array. However it should be realized that
each data point in k space is computed from an integral over all of
the detected ray values.
By changing the projection angle, the expansion coefficient values
throughout k-space can be acquired. These are arranged in radial
fashion as shown in Figure 12.
The Nyquist Theorem requires that the sampling intervals
Dkand Dkr
be equal. Therefore
Dkq = kmax *p/Np
Dkr = 2*kmax/Nr
Np = pNr/2
Where Np is the number of required projections and Nr is the number
of samples along with each projection. This exceeds our previous
simple estimate by a factor of two.
Once k-space has been filled up the data can be interpolated into a
rectangular grid so that the attenuation image can be obtained from
equation 11.
(11)
One problem with the use of the projection theorem to obtain the
data needed for the calculation of the image from equation 11 is that
all of the projection data must be obtained prior to image
reconstruction.
This means that there is a delay in the appearance of the image
following the scan. We will now discuss a technique which permits
image reconstruction to occur simultaneously with data acquisition.
First we will discuss simple back projection, and then add a refinement
required to remove image artifacts. Consider an object in the patient
slice and the attenuation it produces in various projections as shown in
Figure 13A.
If it is assumed that the attenuation associated with each detected
projection is due to a uniform distribution of attenuation along the
projection direction, this attenuation value may be projected back
along the projection direction as shown in Figure 13B.
The mathematical representation of the back projected image B(x,y)
is given by
(16)
where
(17)
and
lies along the x’ projection axis as shown in Figure 14.
This back projected information will reinforce at the location which
produced the attenuation signal recorded in each projection.
However this simple scheme leads to the so called “Star Artifact”
caused by the residual back-projected intensity outside the object.
We can do better by considering how information from an object is
propagated to remote pixels in the image, in other words by
considering the point spread function of the back projection
operation.
Figure 15 shows the process by which an attenuation contribution
from pixel A is communicated to a distant pixel B.
B
∆xB
RAB
A
The total fraction of the attenuation value of pixel A erroneously
communicated to pixel B is proportional to the fraction F of back
projected rays through pixel A which intercept pixel B. This is given
by
The contribution of one pixel to another falls off like 1/ R where R is
the distance between a pixel at r and a pixel at r’ as shown in Figure
16.
y
r
r-r’ = R
r’
x
This spreading of information from one point to another is
equivalent to a convolution of the true image with a 1/R point spread
function. The back projected image B(x,y) is related to the true
image (x,y) by
The solution for (x,y) knowing B(x,y) and the point spread
function is called Filtered Back Projection.
Recall that the convolution theorem states that if
That the Fourier transform of 1 / R is 1 / k is presented without
proof. However, Figure 17, which represents results obtained with an
image processing program called Alice (Hayden Image Processing)
shows the 1 / R point spread function (A) and its Fourier transform
(B).
A
B
The profile through A, shown in C indicates the 1 / R behavior, while
the profile through B, shown in D, indicates the 1 / k behavior.
Similar behavior is found for all rays through these radially
symmetric distributions.
C
D
Equation 21 states that the k space representation of the back projected
image is wrong by a factor of 1/k.
(21)
One solution to obtain the correct attenuation distribution would be
to;
1
Back project to obtain B(x,y)
2
Take the Fourier transform to get
3
Multiply by k to get
4
Take Fourier transform to get (x,y)
Although this would work, a more efficient way which allows image
reconstruction to proceed before all of the data is collected is to
correct the projections P (x’,) prior to back projecting in order to
remove the effects of the point spread function.
In order to do this we note from equation 21 that the k-space
representation of the back projected image is wrong by a factor of
1 / k. The procedure usually used is to multiply this image by a filter
function
designed to remove the 1 / k mistake, i.e.,
Although the temptation is to just use
, the filter function
is usually truncated at the maximum k value which can be adequately
sampled by the pixel matrix.
The filter function and its Fourier transform are shown in Figure 18
along single lines in each space. The width and shape of R(x’) depend
on the details of the pixel size and the chosen cutoff frequency.
Figure 18
The image space version of equation 22 is, according to the
convolution theorem,
(23)
Substituting for the back projected image from equation 16, we have
Although the convolution integral in equation 24 is nominally a two
dimensional integral, the projection data exists only at y’ =0, and the
operation amounts to a linear convolution of each projection with the
radial component of the filter function, which itself has only radial
dependence.
The filtered back projection procedure then amounts to;
1
Measuring projections P(x’,)
2
Convolving with R(x’)
3
Back projecting according to equation 24.
The effect of the filtering is shown in Figure 19 where the projection
data adds constructively at the location of a point object at point A, but
tends to cancel at a point B just away from the actual point.
An estimate of the exposure required to see an attenuation difference
Dwe can make a couple of assumption and apply the signal to noise
equations previously derived, namely
(25)
Where N0 is the required input fluence (assuming unit detector
efficiency),DS/n is the required signal to noise ratio, C is the contrast,
A is the projected area, and T is the object transmission.
To get an approximate estimate for CT we assume that the required
exposure is the same as would be required if all of the exposure were
made in one projection and that the reconstruction process does not
add significant noise to the image(which is usually the case). Two
volume elements (voxels) to be distinguished are shown below.
The pixel size is d on each side, the slice thickness is t. The area is
given by
The contrast is given by
Equation 25 then becomes,
(26)
Note the strong dependence on pixel size d. This is due to the fact
that it enters into the projected area and the contrast in this
geometry.
Suppose we wish to see a 0.6% difference in attenuation with signal to
noise ratio of one. Assume that the average tissue coefficient is 0.19
cm- 1, t= 1cm, d = 1.5 mm and a patient thickness of 30 cm. This gives
and a required input fluence N of
Converting to exposure E in Roentgens we get,
An example of an abdominal CT image (Medicamundi 34/3) is shown
in Figure 22.
The liver, spine and kidneys are clearly evident. The bone is
brightest due to its high attenuation.
CT image values are usually stated in relative attenuation units rel
which are a scaled percentage of the attenuation difference relative
to water as
The relative attenuation coefficients for k = 1000 are shown in
Figure 23 taken from Krestel-Imaging Systems for Medical
Diagnostics
CT images are usually displayed by selecting a range (window) of
attenuation values centered at a chosen attenuation value (window
level). Values within this window are digitally enhanced to fill the
dynamic range of the CT video display. This permits detailed viewing
of structures which are finely separated in attenuation.
window
Level
In conventional computed tomography the slices are obtained by
incrementally moving the patient couch and repeating the scan. In
this approach the spatial resolution in the slice direction is
determined by detector size and collimation.
Another approach, which was proposed and developed by a former
UW Medical Physics graduate, Professor Willi Kalender, now at the
University of Erlangen, is to acquire data as the patient is
continuously moved through the scanner.
The data required to reconstruct a single slice can be chosen from
any portion of the continuous data set providing a continuous
distribution of average slice position.
Interpolation of data is used to produce a consistent set of data
corresponding to a single perpendicular slice through the patient.
Lieber Chuck:
• Spiral acquisition is great!
• The future is MRI!
• But some ignorants will
continue on multislice and
cone-beam spiral CT ....
• Thank you for your support!
Dein Willi
Spiral Galaxy M88
The net result of this process is an improvement in resolution in the
slice direction and more rapid acquisition of volume data. The
scanning geometry is shown in Figure 24.
One of the most promising applications of this approach is the
generation of angiographic images following the intravenous
administration of iodinated contrast material.
Contrast is injected through an arm vein or through a catheter
placed in the superior vena cava. Following a time delay to allow
the contrast to go through the heart, the spiral scan is begun.
Arterial images are generated from the 3D data set using a projection algorithm such
as the Maximum Intensity Pixel (MIP) algorithm. Shaded surface displays, such as in
the image below are also used. The presence of calcium can be a source of image
degradation.
CTA
Rubin et al, Dx Imaging Jan. 1993
XRA
The spiral CT technique has a distinct advantage over previously
tried intravenous injection angiographic techniques in that it produces
a full 3 dimensional data set which can be used to reconstruct views
from arbitrary directions. Previously tried techniques such as
intravenous DSA required additional injections in order to obtain new
views.
A disadvantage which the spiral CT approach shares with
intravenous DSA (discussed in the next section) is the uncertainty in
the timing of the imaging and the arrival of the contrast material.
Cardiac CT Imaging
Technical development & clinical results
Drs J.L. Sablayrolles, F.Besse, C.Jardin,
J.C. Roy, Q.Sénéchal
Centre Cardiologique du Nord - Saint Denis - France
P.Giat, C.Coric GE Medical Systems
Cine-CT of prosthetic mitral valve
in axial view
Cine-CT of prosthetic aortic valve
in coronal view
Electron Beam CT Coronaries
Moshage et al.
Radiology 1995
196:707-714
Since CT acquires a three dimensional data set, the data can be
viewed in easily interpretable surface and volume rendered displays,
especially when spiral CT is used to improve resolution in the slice
direction. Such displays have been very useful in planning surgical
approaches. An example of such a display is shown in Figure 26
which shows a surface display of the bones in the region of the hip.
The 3D data can be used to digitally dissect the body enabling the
surgeon to look inside various cavities and ascertain with precision
the relationships between anatomical objects.
Computed Rotational Angiography
CRA
2D
R Fahrig, AJ Fox, S Lownie and DW Holdsworth
AJNR 18:1507-1514 Sept 1997
Micro CT
Holdsworth and Thornton
Microfil MV-122 contrast agent
Limitations of Current Contrast-Enhanced CT Angiography
Signal
Vein
Artery
Time
Z
Conventional Whole Body CTA Geometry
Time-resolved MR Angiography
45 sec
59 sec
81 sec
Asymmetric enhancement: A-V transit time
T. Carroll, J. Du et al
Signal
Vein
Artery
Time
Z
Whole Body Time-Resolved Acquisition
Current Multiple Rotation Coronary CTA
Best
Case
Worst
Case
Z-Scan CT

Z
Direction of
conventional
table motion
Conventional gantry
Rotation angle.
Rotating detector array
+ focused grid
Z
Patient table
Linear X-ray source
produced by
• scanned electron beam or
• discrete pulsed source
distribution
Cone Beam CT
Cone Beam Detector