Math
584,
Mathematics of Medical Imaging
Matlab Lesson 2:
Images and the Radon Transform
The Radon and
inverse Radon transforms
This is a simple
tutorial for using
Matlab for Math 584, Mathematics of Medical Imaging.
The black text in bold is what you should type at the Matlab command prompt. The non-bold black text which follows is the expected Matlab output. The red text will explain what you should be learning. Please experiment by changing the commands so you can get a feel for the language of Matlab.
In Matlab,
images are real-valued matrices. The
matrix entries represent the intensities or densities at the pixels of
the
image. One way to work with images is to load one from a file using the
imread function.
Download the file circle.jpg
into your Matlab working directory and
then type
>> X=imread('circle.jpg','jpeg');
Matlab can read several different image file formats. Type >>
help imread for more information. To display the image, type
>> figure; imagesc(X); colormap(gray);
The image X is a
real-valued
matrix. You can view the entries just
like any other matrix:
>> size(X), X(100,100)
ans =
256 256 This means that the image is 256x256 pixels.
ans =
255 The value of the
pixel at location (100,100) is 255.
Experiment by defining
your own
images. For example, try:
>> figure; imagesc(eye(150)); colormap(gray);
A phantom is an image
which can be
used to test the numerical accuracy of a 2-D image reconstruction
algorithm. Matlab
has several built-in phantoms. The
default phantom is a ‘modified’ Shepp-Logan
phantom:
>> P=phantom(256);
figure; imagesc(P); colormap(gray);
You can also easily
generate your
own ellipse-based phantoms. Type
>> help phantom for more
information.
For more information on
these
topics, see the following help topics: image, imagesc,
imread, imwrite,
phantom.
Begin with a phantom:
>> P=phantom(256);
To compute the radon
transform of
this matrix, type
>> R=radon(P,0:1:359);
The vector 0:1:359
represents the
discrete set of angles (in degrees) for the Radon transform. Plot the resulting matrix, and spend some time
trying to understand why the image appears as it does.
What are the horizontal and vertical
axes? What symmetries or approximate
symmetries exist? What does the dark
region of the image represent?
>> figure; imagesc(R); colormap(gray);
Now reconstruct the
image using
filtered back projection (this may take several minutes if you have a
slower
computer) :
>> P2=iradon(R,0:1:359); figure; imagesc(P2);
colormap(gray);
You can experiment with
different
options and filters for iradon. See >> help iradon.
We can compare the reconstructed image to the original by displaying
the difference:
>> figure; imagesc(P2-P);
colormap(gray);
This command produces the following error:
Matrix dimensions must agree.
We check
>> size(P2)
ans =
258 258
The
reconstruction has produced a slightly larger image. We can select a
256X256 part to compare to the original image.
>>figure; imagesc(P2(2:257,2:257)-P); colormap(gray);
1. Define an image (matrix) consisting of all zeros except for a 1 at the (m,n)th pixel. What do you expect the Radon transform of this data to look like? Now test your conjecture using Matlab. Try this for a few different pairs (m,n).
2.
Let R be the Radon transform of the 'modified' Shepp-Logan phantom as above. What is
the effect of the Matlab statement
>> R(:,
3. In the final step above we
selected P2(2:257,2:257) to compare to P. There are two other
reasonable possibilities. Find them, compare them to P and state which
of the three possibilities is most like the original image.