inverts the Radon transform and can therefore be used to reconstruct
R = radon(I,theta);
iradon can then be called to
reconstruct the image
I from projection data.
IR = iradon(R,theta);
In the example above, projections are calculated from the original
Note, however, that in most application areas, there is no original
image from which projections are formed. For example, the inverse
Radon transform is commonly used in tomography applications. In X-ray
absorption tomography, projections are formed by measuring the attenuation
of radiation that passes through a physical specimen at different
angles. The original image can be thought of as a cross section through
the specimen, in which intensity values represent the density of the
specimen. Projections are collected using special purpose hardware,
and then an internal image of the specimen is reconstructed by
iradon. This allows for noninvasive imaging
of the inside of a living body or another opaque object.
iradon reconstructs an image from parallel-beam
projections. In parallel-beam geometry, each
projection is formed by combining a set of line integrals through
an image at a specific angle.
The following figure illustrates how parallel-beam geometry is applied in X-ray absorption tomography. Note that there is an equal number of n emitters and n sensors. Each sensor measures the radiation emitted from its corresponding emitter, and the attenuation in the radiation gives a measure of the integrated density, or mass, of the object. This corresponds to the line integral that is calculated in the Radon transform.
The parallel-beam geometry used in the figure is the same as the geometry that was described in Radon Transform. f(x,y) denotes the brightness of the image and is the projection at angle theta.
Parallel-Beam Projections Through an Object
Another geometry that is commonly used is fan-beam geometry,
in which there is one source and n sensors. For
more information, see Fan-Beam Projection. To convert parallel-beam projection
data into fan-beam projection data, use the
iradon uses the filtered backprojection algorithm
to compute the inverse Radon transform. This algorithm forms an approximation
of the image
I based on the projections in the
R. A more accurate result can be obtained
by using more projections in the reconstruction. As the number of
projections (the length of
theta) increases, the
IR more accurately approximates
the original image
I. The vector
contain monotonically increasing angular values with a constant incremental
Dtheta. When the scalar
known, it can be passed to
iradon instead of the
array of theta values. Here is an example.
IR = iradon(R,Dtheta);
The filtered backprojection algorithm filters the projections
R and then reconstructs the image using the
filtered projections. In some cases, noise can be present in the projections.
To remove high frequency noise, apply a window to the filter to attenuate
the noise. Many such windowed filters are available in
The example call to
iradon below applies a Hamming
window to the filter. See the
page for more information. To get unfiltered backprojection data,
'none' for the filter parameter.
IR = iradon(R,theta,'Hamming');
iradon also enables you to specify a normalized
D, above which the filter has zero response.
be a scalar in the range [0,1]. With this option, the frequency axis
is rescaled so that the whole filter is compressed to fit into the
[0,D]. This can be useful in cases
where the projections contain little high-frequency information but
there is high-frequency noise. In this case, the noise can be completely
suppressed without compromising the reconstruction. The following
iradon sets a normalized frequency value
IR = iradon(R,theta,0.85);
The commands below illustrate how to reconstruct an image from
parallel projection data. The test image is the Shepp-Logan head phantom,
which can be generated using the
The phantom image illustrates many of the qualities that are found
in real-world tomographic imaging of human heads. The bright elliptical
shell along the exterior is analogous to a skull, and the many ellipses
inside are analogous to brain features.
Create a Shepp-Logan head phantom image.
P = phantom(256); imshow(P)
Compute the Radon transform of the
phantom brain for three different sets of theta values.
R2 has 36 projections, and
theta1 = 0:10:170; [R1,xp] = radon(P,theta1); theta2 = 0:5:175; [R2,xp] = radon(P,theta2); theta3 = 0:2:178; [R3,xp] = radon(P,theta3);
a plot of one of the Radon transforms of the Shepp-Logan head phantom.
The following figure shows
R3, the transform with
figure, imagesc(theta3,xp,R3); colormap(hot); colorbar xlabel('\theta'); ylabel('x\prime');
Radon Transform of Head Phantom Using 90 Projections
Note how some of the features of the input image appear in this image of the transform. The first column in the Radon transform corresponds to a projection at 0º that is integrating in the vertical direction. The centermost column corresponds to a projection at 90º, which is integrating in the horizontal direction. The projection at 90º has a wider profile than the projection at 0º due to the larger vertical semi-axis of the outermost ellipse of the phantom.
Reconstruct the head phantom image from the projection data created in step 2 and display the results.
I1 = iradon(R1,10); I2 = iradon(R2,5); I3 = iradon(R3,2); imshow(I1) figure, imshow(I2) figure, imshow(I3)
following figure shows the results of all three reconstructions. Notice
I1, which was reconstructed from only
18 projections, is the least accurate reconstruction. Image
which was reconstructed from 36 projections, is better, but it is
still not clear enough to discern clearly the small ellipses in the
lower portion of the image.
I3, reconstructed using
90 projections, most closely resembles the original image. Notice
that when the number of projections is relatively small (as in
the reconstruction can include some artifacts from the back projection.
Inverse Radon Transforms of the Shepp-Logan Head Phantom