CS425 Lab: Image Processing Toolbox and Histograms


1. Image Processing Toolbox

"The Image Processing Toolbox is a a collection of functions that extend the capability of the MATLAB numeric computing environment. The toolbox supports a wide range of image processing operations." (from the online Image Processing Toolbox, User's Guide - "Getting Started" - "What is Image Processing Toolbox?"

1.1 Reading and Displaying Images

1.2 What is stored in an image?

Don't panic. When we get to images, we are still working with matrices that can be accessed through indices. The upper left-hand corner of the image is indexed by (1,1). The general format of indexing is (row, col) or (y,x). Be careful here because you are used to thinking in Cartesian (x, y) coordinates. The following diagram shows the indices to all four corners of a small image:

Coordinates for corners of smile picture

All the basic matrix operations that you have learned still apply. What you might be asking yourself is what is stored in an image matrix? Let's use a grayscale image ("pout.tif") as an example. This image is stored in a 291 x 240 matrix. That means there are 291 rows and 240 columns. Each element corresponds to a single pixel value.

In the case of "pout.tif", the values are between 0 and 255. Zero is full black and 255 is white. Values between zero and 255 are shades of gray.

How to inspect Pixel Values

Note: The "Pixel info" at the bottom displays the value at (column, row) or (x,y). This is different from the MATLAB subscripting syntax, which is in (row, column) or (y,x) format.

1.3 Displaying Image Information

The following table summarizes some ways to get information about an image. These are not specific to the Image Processing Toolbox.

Command Description
whos To get information about size, type, and bytes, of all variables.
whos I For information about an image stored in I.
size(I) To get the size of the image stored in I.
class(I) To get type of data stored in I

Color Formats

There can be different types of data stored in the image array (or matrix). The following table, based on the table on page 2-3 of the Image Processing Toolbox User's Guide, summarizes the different image types and values.

Image Type Intensity Values/Ranges
Binary
  • logical (boolean) array of 0's and 1's
  • 0 is black and 1 is white
Indexed (pseudocolor)
  • values in image are indices into a colormap
  • range of values where p is the length of the colormap:
    • single or double arrays [1, p]
    • logical, uint8, or uint16 [0, p-1]
Grayscale (intensity)
  • Values have ranges as follows:
    • single or double arrays [0, 1]
    • uint8 [0, 255]
    • uint16 [0, 65 535]
    • int16 [-32 768, 32 767]
Truecolor (RGB)
  • Values have the following range:
    • single or double arrays [0, 1]
    • uint8 [0, 255]
    • uint16 [0, 65 535]

Truecolor

Until now you have worked with grayscale images and you have sees that that there are different image types. Most of them work similarly. For each pixel there's a single value, from a range of values appropriate to a data type, stored in a two dimensional matrix. To represent color we have to add another dimension. For instance, if you read image "autumn.tif" into I2 and view the details of I2 using:

I2=imread('autumn.tif');
whos I2

The image size is: 206x345x3

The extra dimension is the color component. For each pixel, we now have a representation of the intensity of red, green, and blue. In autumn.tif, 255 represents the highest intensity, and 0 represents the lowest intensity. For instance, to make a pure red, red would be 255, green would be 0, and blue would be 0.

At pixel (93,180) in autumn.tif,

The following diagram represents the three color components at pixel values. Often, they are represented as channels as is shown in the diagram below. Usually, there is a red channel, a green channel and a blue channel.

Colors in Autumn image

If you only wanted to view the red channel, which of the following statements would you write?

  1. imshow( I2(1, :, :) );
  2. imshow( I2(:, :, 1) );
  3. imshow( I2(:, 1, :) );

What does the image look like?
What happens if you try another color channel?

1.4 Converting Between Different Formats

The following table lists commands provided by the Image Processing Toolbox that convert between the different formats given above.


Image format conversion
(Within the parenthesis you type the name of the image you wish to convert.)

MATLAB command:

Operation:

gray2ind()

Convert from intensity format to indexed format.

ind2gray()

Convert from indexed format to intensity format.

ind2rgb()

Convert from indexed format to RGB format.

mat2gray()

Convert a regular matrix to intensity format by scaling.

rgb2gray()

Convert from RGB format to intensity format.

rgb2ind()

Convert from RGB format to indexed format.

im2bw()

Convert an intensity/indexed/RGB image to a thresholded black and white binary image.

dither()

Convert from intensity/indexed/RGB format to dithered indexed format with custom colormap.
Defaults to black and white if no colormap is provided,

The command mat2gray is useful if you have a matrix representing an image but the values representing the gray scale range between, let's say, 0 and 1000. The command mat2gray automatically re-scales all entries so that they fall between 0 and 255 (if you use the uint8 class) or 0 and 1 (if you use the double class).

As you read through the lab notes, you will see examples of these functions, often to change an indexed or RGB image to grayscale to simplify processing.

Here are some examples of how to use some conversion commands:

%Start with a truecolor image
img=imread('onion.png');
imshow('onion.png')
original onion picture
%Convert to binary/Black and White
img_bin=im2bw(img); 
imshow(img_bin)
thresholded onion
%Convert to binary/Black and White 
%with custom threshold (default is 0.5)
img_bin2=im2bw(img, 0.3);
imshow(img_bin2)
thresholded onion
%Convert to grayscale
img_gray=rgb2gray(img);
imshow(img_gray)
grayscale onion
%reduce to indexed color
[img_indexed, map]=rgb2ind(img,8);
imshow(img_indexed, map)
8-bit dithered
%reduce to indexed color 
%without dithering
[img_nodither, map2]=rgb2ind(img,8, 'nodither');
imshow(img_nodither, map2)
8-bit undithered
%reduce to indexed color 
%using a custom colormap
map_custom=[
   1 1 1;  %white
   1 0 0;  %red
   0 0 0]; %black
img_custommap=dither(img, map_custom);
imshow(img_custommap, map_custom)
dithered to custom color map
%Of course you can swap colormaps...
map_custom2=[
   1 1 1;  %white
   0 0 1;  %blue
   0 0 0]; %black
imshow(img_custommap, map_custom2)
color map substitution

1.5 Saving Images

To save an image, use the imwrite function. The parameters you send to the function depend on the file format and image type.

Saving a Binary, Grayscale or Truecolor Image

When saving a binary, grayscale or truecolor image, all you need to specify are the matrix and filename. For example:

 imwrite(img_bin, 'onion_bw.png')

In this example, the image would be stored as a 1-bit image in png format. The format is automatically guessed from the file extension.

Saving an Indexed Image

When you save an indexed image, it is important to provide the colormap like this:

imwrite(img_indexed, map, 'onion.gif');

Notes

Options

There are various options that you can send to imwrite to control how different file types save. For example you can change the compression level of a JPEG like this:

imwrite(img, 'onion.jpg', 'Quality', 25);

2. Histograms

See Chapter 4 in your textbook.

Histograms are a way of visualizing the predominant intensities of an image. As a definition, image histograms are a count of the number of pixels that are at a certain intensity. When represented as a plot, the x-axis is the intensity value, and the y-axis is the number of pixels with that intensity value. The following are examples of histograms with predominately low, mid and high range intensities.

Which image would you expect to be darkest? lightest?

index_tire_hist

index_pout_hist

index_eight_hist

Notice that the x-axis is the intensity value from 0 to 255 (these images are uint8). The y-axis varies depending on the number of the pixels in the image and how their intensities are distributed.

MATLAB easily displays image histograms using the function imhist(I). The above plots were created with the following commands:

    tire=imread('tire.tif');
    figure,imhist(tire);
    title('Histogram of tire.tif--Predominately Low Intensity');
    pout=imread('pout.tif');
    figure,imhist(pout);
    title('Histogram of pout.tif--Predominately Mid-range Intensity');
    eight=imread('eight.tif');
    figure, imhist(eight);
    title('Histogram of eight.tif--Predominately High Intensity');
  

2.1 Histogram Equalization

See Chapter 5.5 in your textbook.

The idea behind Histogram Equalization is that we try to evenly distribute the occurrence of pixel intensities so that the entire range of intensities is used more fully. We are trying to give each pixel intensity equal opportunity; thus, equalization. Especially for images with a wide range of values with detail clustered around a few intensities, histograms will improve the contrast in the image.

In MATLAB, the function to perform Histogram Equalization is histeq(I).

Following are some examples of Histogram Equalization and its results on the images that have predominately low and mid range intensities:

2.1.1 Low Range Intensities

  Image Histogram
Original index_tire_orig index_tire2_hist
With Histogram Equalization index_tire_eq index_tire_eq_hist

The above diagrams were created using the following calls:

    tire=imread('tire.tif'); 
    imshow(tire);
    figure, imhist(tire); 
    tire_eq=histeq(tire); 
    figure, imshow(tire_eq);
    figure, imhist(tire_eq);
  

2.1.2 Mid Range Intensities

  Image Histogram
Original index_pout_orig index_pout2_hist
With Histogram Equalization index_pout_eq index_pout_eq_hist

The above diagrams were created using the following calls:

pout=imread('pout.tif'); 
imshow(pout);
pout=rgb2gray(pout);
figure, imhist(pout)
pout_eq=histeq(pout); 
figure, imshow(pout_eq);
figure, imhist(pout_eq);
  
Note: imshow(image,[low, high]) also performs a histogram stretch: all values less than or equal to low will be shown as black;all values greater than or equal to high will be shown as white; the rest will be stretched from 0 to 255.
imshow(pout, [75, 160]) index_pout_eq

2.2 Histogram Matching

See Chapter 5.6 in your textbook.

Sometimes, histogram equalization does not produce the contrast or results that we expect. Consider the following example, taken from Digital Image Processing, Using MATLAB[1]. The original image is courtesy of NASA.

  Image Histogram
Original index_moon_orig index_moon_hist
With Histogram Equalization index_moon_eq index_moon_eq_hist

The above images were created using these steps:

  1. Download Fig0310(a)(MoonPhobos).tiff
  2. Find it in the Current Directory window.
    • Mac OS X Hint: if you don't see the file, you can navigate using Unix commands. It is likely in either the Desktop or Downloads folder which are subfolders of your home folder. To go to your home folder use cd ~.
  3. Double click on Fig0310(a)(MoonPhobos).tiff—it will now appear in your Workspace window with the name Fig03100x28a0x290x28MoonPhobos0x29.
  4. Use the following code in the Command Window.
    • Hint: press tab to autocomplete long names.
        moon=Fig03100x28a0x290x28MoonPhobos0x29;
        imshow(moon);
        figure, imhist(moon);
        moon_eq=histeq(moon);
        figure, imshow(moon_eq);
        figure, imhist(moon_eq);

You will notice that there are very few low intensity values in the histogram equalized image above. To fix this, we will use histogram matching (or histogram specification).

To provide some background, histogram equalization tries to create an equal probability of each intensity occurring. This equal probability is represented by a horizontal line. If you wanted to emulate the results above, you could pass a horizontal line as a second argument to histeq. The following code can be used to yield results similar to the above histogram equalization

    p(1:256)=1  
    g=histeq(moon,p);  
    figure, imshow(g);
    figure, imhist(g);

With histogram matching, we provide a multimodal (multi-peaked) Gaussian function as an argument to histeq. In contrast to the horizontal line of histogram equalization, the Gaussian function specifies that certain intensity ranges will have more probability of occurring than others. To speak intuitively, by passing a multimodal Gaussian function as a second argument to histeq, the modified image histogram will have peaks that approximately match the peaks of the Gaussian function. Note this approximate match in the the moon example below:

index_bimodal_guass

index_moon_match_hist

index_moon_match

In order to create the above results, two M-functions were used[1]

With the above two functions in the "Current Directory", the following calls were made to produce the results in the above table:

    p=manualhist;
    %the following arguments were typed at the prompt  
    % 0.15 0.05 0.75 0.05 1 0.1 0.002  - input is looped don't forget to end with x 
    g=histeq(moon, p);  
    figure, imshow(g);  
    figure, imhist(g);

 


3. The Return of the M-Files

At the end of this section, we are going to take a look at a function, written in MATLAB, that is able to apply all of the intensity transformations mentioned above. Before doing that, the following provides some background information about M-Files and functions.

3.1 Scripts

In Lab 1, we briefly looked at M-Files. Specifically, we saw an m-file script. Remember:

Scripts do not take input arguments or return anything as output. However, anything that is created when the script is run remains in the workspace.

3.2 Functions

By contrast, functions can accept input arguments and return output. Functions are distinguished from other M-files by the keyword function in the first line of code. Let's take a look at a sample function acosd. The code below can be seen in MATLAB using the command:
type acosd

function y = acosd(x)
%ACOSD  Inverse cosine, result in degrees.
%   ACOSD(X) is the inverse cosine, expressed in degrees,
%   of the elements of X.
%
%   Class support for input X:
%      float: double, single
%
%   See also COSD, ACOS.

%   Copyright 1984-2004 The MathWorks, Inc. 
%   $Revision: 1.1.6.3 $  $Date: 2004/06/25 18:51:38 $

if ~isreal(x)
    error('MATLAB:acosd:ComplexInput', 'Argument should be real.');
end

y = 180/pi*acos(x);

A couple of things to note about this code are the following:

3.3 Function Inputs and Outputs

MATLAB functions are very flexible compared to those of many other languages. They can have multiple and varying numbers of inputs—which is pretty normal. They can also have multiple and varying numbers of outputs—which may be new to you. To help you manage this complexity, MATLAB has some special functions. The following table summarizes these functions:

Function Description
nargin returns the number of arguments input into the function
nargout returns the number of arguments output from a function
nargchk(low, high, number)

returns "Not enough input parameters" if number is less than low and returns "Too many input parameters" if number is greater than high. If the value of number is between (or equal to) low and high, then an empty matrix is returned. Often, nargchk is used in conjunction with error. The error function stops execution and prints the message if number is either too high or too low. A sample usage of nargchk inside a function is:
error(nargchk(2, 3, nargin));
%This checks if the input arguments (nargin) are between 2
%and 3 (inclusive). If they aren't, then execution stops and
%an error message (with the string returned from nargchk) is
%printed

varargin accepts a variable number of arguments
varargout returns a variable number of arguments

3.3.1 varargin

To write a function with a variable number of arguments, you can write:
function m=testFunction(varargin)
This means that you could have zero to as many arguments as you want for testFunction. The following code is meant to test this.

function m=testFunction(varargin)
%TESTFUNCTION tests variable number of arguments
%   The following are valid calls to testFunction
%   a=testFunction
%   b=testFunction(val)
%   c=testFunction(val1,val2)
%   d=testFunction(val1,val2,...moreval)


number=nargin;

switch number
    case 0
        m=0;
    case 1
        m=1;
    case 2
        m=2;
    otherwise
        m='more than two';
end

Notice that there are no return statements, instead, m (the return variable) is assigned a value. For instance, if you have zero arguments, then m will be set to 0; if you have one argument, then m will be set to 1; if you have three or more arguments, then m will be assigned a message of: 'more than two'.

If you want to access any one of the arguments sent as a variable argument, you can use the notation varargin{#}. For instance, varargin{1} will return the first argument. The following function constructs m as an array that contains all of the input arguments.

function m=testFunction2(varargin)
%TESTFUNCTION2 displays all of arguments sent to the function
%   The following are valid calls to testFunction2
%   a=testFunction
%   b=testFunction(val)
%   c=testFunction(val1,val2)
%   d=testFunction(val1,val2,moreval)


number=nargin;
% m=[];  %this goes with the commented line below


if number>0
    for i=1:number
       m(i)=varargin{i};  %add one more argument to m
       % m=[m varargin{i}];  %this way works as well
    end
    
else
    m='no arguments';
end

This function uses if/else and for statements. Notice that end is used to end these statements. The for makes use of the colon operator; in each loop, i will have an incremented value from 1 to the number of arguments.

You can use varargin with other arguments, but varargin must come last. For example, if you require one argument and then a variable number of arguments, you would write this:
function m=testFunction3(x, varargin)
In this case, you require at least one argument and it will be stored in x.

3.3.2 Varargout

The idea behind varargout is not far from varargin. However, first we will describe something that you might not be familiar with in other languages—MATLAB can return multiple arguments. To have a function with multiple returns, you can write:
function [m, n] = testFunction4(x)
This indicates that there are two returns to be stored in m and n.

You can also specify a variable number of arguments to be returned. The following code returns a variable number of arguments depending on the number of arguments specified when the function was called:

function [varargout] = testFunction5(x)
%TESTFUNCTION5 tests how variable arguments work
%  testFunction5(x) does not return anything
%  m=testFunction5(x) returns one argument (1+x)
%  [m,n]=testFunction5(x) returns two arguments (1+x and 2+x)
%  [m,n,o]=testFunction5(x) returns three arguments (1+x, 2+x, and 3+x)
%  so on.

number=nargout;

if number>0
    for i=1:number
        varargout{i}=i+x;
    end
end

For instance, if you make a call like:
[m,n,o]=testFunction5(1)
The following will be returned:

m =
     2

n =
     3

o =

     4

5. Exercise

Part 1: Write a Function

Create a MATLAB function called returnSubImage. Given an image and starting x and y locations, return a square subimage. If the size of the subimage is not specified, then the default size is 100 x 100. This function will take, in this order, three arguments and an optional fourth argument:

Include error checking into your function. If the square is outside of the bounds of the original image, then produce an error message.You can use the size function to find out the size of the image with which you are dealing.

To help you create the M-File you can use the MATLAB editor. Access it through File | New | Blank M–File or File | New | Function M-File

File | New | Blank M-File cmd+N

When you have completed your code, try the following calls to returnSubImage; the results should be consistent to what is described below:

  1. % read in the image to get things started
    lifting = imread('liftingbody.png');
    imshow(lifting);

    Original Lifting Body image
  2. % get the subimage of the closer aircraft
    lifting_sub1=returnSubImage(lifting,159,33,330);
    figure, imshow(lifting_sub1)

    Lifting Body sub image
  3. % get the subimage of the second aircraft using the default 100x100 square
    lifting_sub2=returnSubImage(lifting,110,360)
    figure, imshow(lifting_sub2)

    Lifting Body sub image
  4. % test error checking
    lifting_sub=returnSubImage(lifting,159,33,500)

Tips

Part 2: Histograms

  1. Download the following image "flood.jpg" and store it in MATLAB's "Current Directory".
    flood.jpg
  2. Load the image data.
  3. Display the histogram of the image.
  4. Use histogram equalization to create a new image with more contrast.
  5. Display the histogram of the image created in step 4.
  6. Use histogram matching (using M-Files from above notes) to produce an image that looks something like the following:
    food_match.jpg
  7. Display the histogram of the image created in step 6.

Tips

Deliverables:

 


5. References

  1. Digital Image Processing, Using MATLAB, by Rafael C. Gonzalez, Richard E. Woods, and Steven L. Eddins, pp. 86-88
  2. Digital Image Processing, Using MATLAB, by Rafael C. Gonzalez, Richard E. Woods, and Steven L. Eddins, pp. 73-74
  3. Digital Image Processing, Using MATLAB, by Rafael C. Gonzalez, Richard E. Woods, and Steven L. Eddins, pp. 562-563
  4. Getting Started with MATLAB 7 (MATLAB's documentation)--available through MATLAB's help menu or online at:
    http://www.mathworks.com/access/helpdesk/help/techdoc/matlab.html
  5. Image Processing Toolbox, User's Guide, Version 5 (MATLAB's documentation)--available through MATLAB's help menu or online at:
    http://www.mathworks.com/access/helpdesk/help/toolbox/images/