Getting Started (Image Processing Toolbox)

Image Processing Toolbox

Example 2 -- Advanced Topics

In this exercise you will work with another intensity image, rice.tif, and explore some more advanced operations. The goals of this exercise are to remove the nonuniform background from rice.tif, convert the resulting image to a binary image by using thresholding, use components labeling to return the number of objects (grains or partial grains) in the image, and compute object statistics.

1. Read and Display an Image

Clear the MATLAB workspace of any variables and close open figure windows. Read and display the intensity image rice.tif.

clear, close all
I = imread('rice.tif');
imshow(I)

2. Use Morphological Opening to Estimate the Background

Notice that the background illumination is brighter in the center of the image than at the bottom. Use the imopen function to estimate the background illumination.

background = imopen(I,strel('disk',15));

To see the estimated background image, type

```
imshow(background)
```

Here's What Just Happened

Step 1. You used the toolbox functions imread and imshow to read and display an 8-bit intensity image. imread and imshow are discussed in Exercise 1, in 2. Check the Image in Memory, under the "Here's What Just Happened" discussion.
Step 2. You performed a morphological opening operation by calling imopen with the input image, I, and a disk-shaped structuring element with a radius of 15. The structuring element was created by the strel function. The morphological opening has the effect of removing objects that cannot completely contain a disk of radius 15. For more information about morphological opening, see Dilation- and Erosion-Based Functions.

3. Display the Background Approximation as a Surface

Use the surf command to create a surface display of the background approximation, background. The surf function requires data of class double, however, so you first need to convert background using the double command.

figure, surf(double(background(1:8:end,1:8:end))),zlim([0 255]);
set(gca,'ydir','reverse');

The example uses MATLAB indexing syntax to view only 1 out of 8 pixels in each direction; otherwise the surface plot would be too dense. The example also sets the scale of the plot to better match the range of the uint8 data and reverses the y-axis of the display to provide a better view of the data (the pixels at the bottom of the image appear at the front of the surface plot).

Here's What Just Happened

Step 3. You used the surf command to examine the background image. The surf command creates colored parametric surfaces that enable you to view mathematical functions over a rectangular region. In the surface display, [0, 0] represents the origin, or upper-left corner of the image. The highest part of the curve indicates that the highest pixel values of background (and consequently rice.tif) occur near the middle rows of the image. The lowest pixel values occur at the bottom of the image and are represented in the surface plot by the lowest part of the curve.
The surface plot is a Handle Graphics^® object, and you can therefore fine-tune its appearance by setting properties. For information on working with MATLAB graphics, see the MATLAB graphics documentation.

4. Subtract the Background Image from the Original Image

Now subtract the background image, background, from the original image, I, to create a more uniform background.

```
I2 = imsubtract(I,background); 
```

Now display the image with its more uniform background.

```
figure, imshow(I2)
```

Here's What Just Happened

Step 4. You subtracted a background approximation image from rice.tif. Because subtraction, like many of MATLAB mathematical operations, is only supported for data of class double, you must use the Image Processing Toolbox image arithmetic imsubtract function.

The Image Processing Toolbox has a demo, ipss003, that approximates and removes the background from an image. For information on how to run this (and other demos), see Image Processing Demos in the Preface.

5. Adjust the Image Contrast

The image is now a bit too dark. Use imadjust to adjust the contrast.

I3 = imadjust(I2, stretchlim(I2), [0 1]);

Display the newly adjusted image.

```
figure, imshow(I3);
```

Here's What Just Happened

Step 5. You used the imadjust command to increase the contrast in the image. The imadjust function takes an input image and can also take two vectors: [low high] and [bottom top]. The output image is created by mapping the value low in the input image to the value bottom in the output image, mapping the value high in the input image to the value top in the output image, and linearly scaling the values in between. See the reference pages for imadjust for more information.

You called imadjust with stretchlim(I2) as the second argument. The stretchlim function automatically computes the right [low high] values to make imadjust increase (stretch) the contrast of the image.

6. Apply Thresholding to the Image

Create a new binary thresholded image, bw, by using the functions graythresh and im2bw.

level = graythresh(I3);
bw = im2bw(I3,level); 
figure, imshow(bw)

Now call the whos command to see what type of array the thresholded image bw is.

```
whos
```

MATLAB responds with

Name                Size         Bytes  Class

I                 256x256        65536  uint8 array
I2                256x256        65536  uint8 array
I3                256x256        65536  uint8 array
background        256x256        65536  uint8 array
bw                256x256        65536  logical array
level               1x1              8  double array

Grand total is 327681 elements using 327688 bytes

Here's What Just Happened

Step 6. You called graythresh to automatically compute an appropriate threshold to use to convert the intensity image to binary. You then called im2bw to perform for thresholding, using the threshold, level, returned by graythresh.

Notice that when you call the whos command, you see the expression logical listed after the class for bw. This indicates the presence of a logical flag. The flag indicates that bw is a logical matrix, and the Image Processing Toolbox treats logical matrices as binary images. Thresholding using MATLAB logical operators always results in a logical image. For more information about binary images and the logical flag, see Binary Images.

7. Determine the Number of Objects in the Image

To determine the number of grains of rice in the image, use the bwlabel function. This function labels all of the connected components in the binary image bw and returns the number of objects it finds in the image in the output value, numobjects.

[labeled,numObjects] = bwlabel(bw,4);% Label components.
numObjects =

    80

The accuracy of your results depends on a number of factors, including:

The size of the objects
The accuracy of your approximated background
Whether you set the connectivity parameter to 4 or 8
Whether or not any objects are touching (in which case they may be labeled as one object) In the example, some grains of rice are touching, so bwlabel treats them as one object.

Here's What Just Happened

Step 7. You called bwlabel to search for connected components and label them with unique numbers. bwlabel takes a binary input image and a value specifying the connectivity of objects. The parameter 4, passed to the bwlabel function, means that pixels must touch along an edge to be considered connected. For more information about the connectivity of objects, see Pixel Connectivity.
You can also determine the number of objects in a label matrix by asking for the maximum pixel value in the image. For example,
max(labeled(:))
ans =
80

8. Examine the Label Matrix

You may find it helpful to take a closer look at labeled to see what bwlabel has created. Use the imcrop command to select and display pixels in a region of labeled that includes an object and some background.

To ensure that the output is displayed in the MATLAB window, do not end the line with a semicolon. In addition, choose a small rectangle for this exercise, so that the displayed pixel values don't wrap in the MATLAB command window.

The syntax shown below makes imcrop work interactively. Your mouse cursor becomes a cross-hair when placed over the image. Click at a position in labeled where you would like to select the upper left corner of a region. Drag the mouse to create the selection rectangle, and release the button when you are done.

grain = imcrop(labeled) % Crop a portion of labeled.

We chose the left edge of a grain and got the following results.

grain =

     0     0     0     0     0     0     0    60    60
     0     0     0     0     0    60    60    60    60
     0     0     0    60    60    60    60    60    60
     0     0     0    60    60    60    60    60    60
     0     0     0    60    60    60    60    60    60
     0     0     0    60    60    60    60    60    60
     0     0     0    60    60    60    60    60    60
     0     0     0     0     0    60    60    60    60
     0     0     0     0     0     0     0     0     0
     0     0     0     0     0     0     0     0     0

A good way to view a label matrix is to display it as a pseudo-color indexed image. In the pseudo-color image, the number that identifies each object in the label matrix maps to a different color in the associated colormap matrix. When you view a label matrix as an RGB image, the objects in the image are easier to distinguish.

To view a label matrix in this way, use the label2rgb function. Using this function, you can specify the colormap, the background color, and how objects in the label matrix map to colors in the colormap.

RGB_label = label2rgb(labeled, @spring, 'c', 'shuffle');
imshow(RGB_label);

Here's What Just Happened

Step 8. You called imcrop and selected a portion of the image that contained both some background and part of an object. The pixel values were returned in the MATLAB window. If you examine the results above, you can see the corner of an object labeled with 60's, which means that it was the 60th object labeled by bwlabel.

The imcrop function can also take a vector specifying the coordinates for the crop rectangle. In this case, it does not operate interactively. For example, this call specifies a crop rectangle whose upper-left corner begins at (15, 25) and has a height and width of 10.

rect = [15 25 10 10]; roi = imcrop(labeled, rect)

You are not restricted to rectangular regions of interest. The toolbox also has a roipoly command that enables you to select polygonal regions of interest. Many image processing operations can be performed on regions of interest, including filtering and filling. See Region-Based Processing for more information.

The call to label2rgb illustrates a good way to visualize label matrices. The pixel values in the label matrix are used as indices into a colormap. Using label2rgb, you can specify your own colormap or use one of the MATLAB colormap-creating functions, including gray, pink, spring, and hsv. For information on these functions, see colormap in the MATLAB Function Reference.

9. Measure Object Properties in the Image

The regionprops command measures object or region properties in an image and returns them in a structure array. When applied to an image with labeled components, it creates one structure element for each component. Use regionprops to create a structure array containing some basic properties for labeled.

graindata = regionprops(labeled,'basic')

MATLAB responds with

graindata = 

80x1 struct array with fields:
    Area
    Centroid
    BoundingBox

To find the area of the component labeled with 51's, use dot notation to access the Area field in the 51st element in the graindata structure array. Note that structure field names are case sensitive, so you need to capitalize the name as shown.

```
graindata(51).Area
```

returns the following results

```
ans =

   296
```

To find the smallest possible bounding box and the centroid (center of mass) for the same component, use this code:

graindata(51).BoundingBox, graindata(51).Centroid
ans =

  142.5000   89.5000   24.0000   26.0000
ans =

  155.3953  102.1791

To create a new vector, allgrains, which holds just the area measurement for each grain, use this code:

allgrains = [graindata.Area];
whos allgrains

Call the whos command to see how MATLAB allocated the allgrains variable.

  Name            Size         Bytes  Class

  allgrains       1x80           640  double array

Grand total is 80 elements using 640 bytes

allgrains is a one-row array of 80 elements, where each element contains the area measurement of a grain. Check the area of the 51st element of allgrains.

```
allgrains(51)
```

returns

```
ans =

   296
```

which is the same result that you received when using dot notation to access the Area field of graindata(51).

Here's What Just Happened

Step 9. You called regionprops to return a structure of basic property measurements for each thresholded grain of rice. The regionprops function supports many different property measurements, but setting the properties parameter to 'basic' is a convenient way to return three of the most commonly used measurements: the area, the centroid (or center of mass), and the bounding box. The bounding box represents the smallest rectangle that can contain a region, or in this case, a grain. The four-element vector returned by the BoundingBox field,

[142.5000 89.5000 24.0000 26.0000]

shows that the upper left corner of the bounding box is positioned at [142.5 89.5], and the box has a width of 24.0 and a height of 26.0. (The position is defined in spatial coordinates, hence the decimal values. For more information on the spatial coordinate system, see Spatial Coordinates.) For more information about working with MATLAB structure arrays, see Structures in the MATLAB programming and data types documentation.

You used dot notation to access the Area field of all of the elements of graindata and stored this data to a new vector allgrains. This step simplifies analysis made on area measurements because you do not have to use field names to access the area.

10. Compute Statistical Properties of Objects in the Image

Now use MATLAB functions to calculate some statistical properties of the thresholded objects. First use max to find the size of the largest grain. (If you have followed all of the steps in this exercise, the "largest grain" is actually two grains that are touching and have been labeled as one object).

```
max(allgrains)
```

returns

```
ans =

   695
```

Use the find command to return the component label of this large-sized grain.

```
biggrain = find(allgrains==695)
```

returns

```
biggrain =

    68
```

Find the mean grain size.

```
mean(allgrains)
```

returns

```
ans =

   249
```

Make a histogram containing 20 bins that show the distribution of rice grain sizes.

hist(allgrains,20)


Here's What Just Happened
Step 10. You used some of the MATLAB statistical functions, `max`, `mean`, and `hist` to return the statistical properties for the thresholded objects in `rice.tif`. The Image Processing Toolbox also has some statistical functions, such as `mean2` and `std2`, which are well suited to image data because they return a single value for two-dimensional data. The functions `mean` and `std` were suitable here because the data in `allgrains` was one dimensional. The histogram shows that the most common sizes for rice grains in this image are in the range of 300 to 400 pixels.

Example 1 -- Some Basic Topics Where to Go from Here