# Math 584, Mathematics of Medical Imaging

Math 584, Mathematics of Medical Imaging. 1

# Matlab Lessons

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.

## Matlab Lesson 2: Images and the Radon Transform

### Images and Phantoms

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

>> 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);

Begin with a phantom:

>> P=phantom(256);

To compute the radon transform of this matrix, type

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) :

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:

??? Error using ==> minus
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);

### Exercises (Lesson 2)

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(:,10:20)=0.  How does this change affect the reconstructed image?

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.

## Matlab tips

• Always end Matlab statements with a semicolon, unless you want Matlab to display the resulting value of the statement (or expression).  You can get the value of any variable by simply typing its name (press enter without a semicolon).
• You can use the up and down arrows on your keyboard to retrieve previously executed Matlab statements.  For more efficient retrieval, type the first few characters of the statement before using the up arrow.
• From the Matlab prompt, type help [function name] to get help on any built-in Matlab function.
• It is sometimes tempting to use "for loops" to define and manipulate matrices and arrays.  However this should be avoided whenever possible because Matlab is very inefficient when implementing such loops.  There are certain Matlab operators and functions which can often be used in place of "for loops".  See the following help topics: times, divide, colon, find.
• You can combine a sequence of Matlab statements into a single Matlab-executable file called an m-file.  See help function.
• You can save large data objects in files (for future use).  See the help topics: save, load.
• You can have more than one statement on a single line: just separate the statements by semicolons.
• Have you discovered other useful tips and tricks?  Send them to cle@math.upenn.edu .