CS425 Lab: Image Processing Toolbox and Histogram


1. Image Processing Toolbox

1.1 What is the Image Processing Toolbox?

1.2 How to read and display an image

1.3 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). Be careful here because you are used to thinking in x and y. The following diagram shows the indices to all four corners of a small image:

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 zero and 255. Zero is full black and 255 is white. Values between zero and 255 are shades of gray.

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

1.4 How to display information about the image

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

1.5 Different types of images

There can be different types of data stored in the image array (or matrix). The following table, modified from 2-3 of the Image Processing Toolbox User's Guide, is meant to summarize the different image types and values.

Image Type Intensity Values/Ranges
Binary
  • logical array of 0's and 1's
  • where 0 is black and 1 is white
Indexed (pseudocolor)
  • single or double arrays [1, p]
  • logical, uint8, or uint16 [0, p-1]
  • where p is the length of the colormap
Grayscale (intensity)
  • single or double arrays [0, 1]
  • uint8 [0, 255]
  • uint16 [0, 65 535]
  • int16 [-32 768, 32 767]
Truecolor (RGB)
  • single or double arrays [0, 1]
  • uint8 [0, 255]
  • uint16 [0, 65 535]

 

1.6 More about Truecolor

Now that you know a little about grayscale and that there are different image types, let's introduce truecolor. With truecolor, 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 red text in the diagram below. For instance, there is a red channel, a green channel and a blue channel.

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?

1.7 convert between different formats


The following table shows how to convert between the different formats given above. All these commands require the Image processing tool box!


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

Operation:

Matlab command:

Convert between intensity/indexed/RGB format to binary format.

dither()

Convert between intensity format to indexed format.

gray2ind()

Convert between indexed format to intensity format.

ind2gray()

Convert between indexed format to RGB format.

ind2rgb()

Convert a regular matrix to intensity format by scaling.

mat2gray()

Convert between RGB format to intensity format.

rgb2gray()

Convert between RGB format to indexed format. rgb2ind()

Convert image to binary image.

im2bw()

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 within 0 and 255 (if you use the uint8 class) or 0 and 1 (if you use the double class).

The following images show the results of some basic format conversion on the RGB image "Bambi"( from the disney website).

Original image: img=imread('bambi.jpg')
graydeer=rgb2gray(img)
deer deer gray
binarydeer=im2bw(img) transposeddeer=graydeer'
deer binary transposed gray deer

 


2. Histograms

2.1. What is an Image Histogram?

Histograms are a way of visualizing the predominate 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 how the x-axis is the intensity value from 0 to 256 (the images are uint8). The y-axis varies depending on the number of pixels in the image.

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

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.2. Histogram Equalization

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 predominately low intensities, histograms will improve the contrast in the image.

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

The following provide some examples of using Histogram Equalization and the results on the images that have predominately low and mid range intensities:

2.2.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.2.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.jpg'); 
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 has the functioning of histogram stretch: all value that no larger than low will be shown as black; 

all pixel value that no less than high will be shown as white; the rest will be stretched to 0~255.
imshow(pout, [75, 160]) index_pout_eq

2.3. Histogram Matching

Sometimes, histogram equalization does not produce the contrast or results that we expect. Consider the following example, taken from page 86 of Digital Image Processing, Using MATLAB, by Rafael C. Gonzalez, Richard E. Woods, and Steven L. Eddins. 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 the following MATLAB calls:

%double click on downloaded Fig0310(a)(MoonPhobos).tif   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 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 from pages 86 to 88 of Digital Image Processing, Using MATLAB, by Rafael C. Gonzalez, Richard E. Woods, and Steven L. Eddins were used. The two functions are:

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  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, any variables that are created when the script is run remain in the workspace.

3.2 Functions

By contrast, functions can accept input arguments and return output. Functions are clearly identified with 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 More about Input and Output Arguments to Functions

MATLAB has some special functions that can be used with regards to the input and output arguments of a function. 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.4 Varargin

We will discuss vargin and varargout a little more here. 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.

It is not necessary, to use varargin only on its own. If you require one argument and then a variable number of arguments, you can write syntax like:
function m=testFunction3(x, varargin)
In this case, you will require at least one argument and it will be stored in x.

3.5 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

3.6 The intrans and changeclass Functions

The following code (taken with permission from pages 73 and 74 of Digital Image Processing, Using MATLAB, by Rafael C. Gonzalez, Richard E. Woods, and Steven L. Eddins) is capable of doing all of the intensity transformations mentioned above.

function g = intrans(f, varargin)
%INTRANS Performs intensity (gray-level) transformations.
%  G = INTRANS(F, 'neg') computes the negative of input image F
%
%  G = INTRANS(F, 'log', C, CLASS) computes C*log(1 + F) and
%  multiplies the result by (positive) constant C.  If the last two
%  parameters are omitted, C defaults to 1.  Because the log is used 
%  frequently to display Fourier spectra, parameter CLASS offers the 
%  option to specify the class of the output as 'uint8' or 
%  'uint16'. If parameter CLASS is omitted, the output is of the 
%  same class as the input.
%
%  G = INTRANS(F, 'gamma', GAM) performs a gamma transformation on
%  the input image using parameter GAM (a required input).
%
%  G = INTRANS(F, 'stretch', M, E) computes a contrast-stretching
%  transformation using the expression 1./(1 + (M./(F + eps)).^E).
%  Parameter M must be in the range [0, 1].  The default value for 
%  M is mean2(im2double(F)), and the default value for E is 4.
%
%  For the 'neg', 'gamma', and 'stretch' transformations, double
%  input images whose maximum value is greater than 1 are scaled 
%  first using MAT2GRAY. Other images are converted to double first
%  using IM2DOUBLE.  For the 'log' transformation, double images are
%  transformed without being scaled; other images are converted to 
%  double first using IM2DOUBLE.
%
%  The output is of the same class as the input, except if a 
%  different class is specified for the 'log' option.

%  Verify the correct number of inputs.
error(nargchk(2, 4, nargin))

%  Store the class of the input for use later.
classin = class(f);

%  If the input is of class double, and it is outside the range
%  [0, 1], and the specified transformation is not 'log', convert the 
%  input to the range [0, 1].
if strcmp(class(f), 'double') & max(f(:)) > 1 ...
      ~strcmp(varargin{1}, 'log')
   f = mat2gray(f);
else % Convert to double, regardless of class(f).
   f = im2double(f);
end

%  Determine the type of transformation specified.
method = varargin{1};

%  Perform the intensity transformation specified.
switch method
case 'neg'
   g = imcomplement(f);
case 'log'
   if length(varargin) == 1
      c = 1;
   elseif length(varargin) == 2
      c = varargin{2};
   elseif length(varargin) == 3
      c = varargin{2};
      classin = varargin{3};
   else
      error('Incorrect number of inputs for the log option.')
   end
   g = c*(log(1 + double(f)));
case 'gamma'
   if length(varargin) < 2
      error('Not enough inputs for the gamma option.')
   end
   gam = varargin{2};
   g = imadjust(f, [ ], [ ], gam);
case 'stretch'
   if length(varargin) == 1
      % Use defaults.
      m = mean2(f);
      E = 4.0;
   elseif length(varargin) == 3
      m = varargin{2};
      E = varargin{3};
   else error('Incorrect number of inputs for the stretch option.')
   end
   g = 1./(1 + (m./(f + eps)).^E);
otherwise
   error('Unknown enhancement method.')
end

%  Convert to the class of the input image.
g = changeclass(classin, g);

The above code relies on a function called changeclass, which is on pages 562 and 563 of the Digital Image Processing, Using MATLAB, by Rafael C. Gonzalez, Richard E. Woods, and Steven L. Eddins:

function image = changeclass(class, varargin)
%CHANGECLASS changes the storage class of an image.
%   I2 = CHANGECLASS(CLASS, I);
%   RGB2 = CHANGECLASS(CLASS, RGB);
%   BW2 = CHANGECLASS(CLASS, BW);
%   X2 = CHANAGECLASS(CLASS, X, 'indexed');
%   Copyright 1993-2002 The MathWorks, Inc
%   $Revision: 1.2 $ $Date:2003/02/19 22:09:58 $

switch class
case 'uint8'
   image = im2uint8(varargin{:});
case 'uint16'
   image = im2uint16(varargin{:});
case 'double'
   image = im2double(varargin{:});
otherwise
   error('Unsupported IPT data class.');
end

The first part of the intrans function describes how to use it. The following table provides some examples of using intrans to correspond to the four Intensity Transformation Functions. Assume that I=imread('tire.tif');

Transformation Intensity Transformation Function Corresponding intrans Call
photographic negative neg=imcomplement(I); neg=intrans(I,'neg');
logarithmic I2=im2double(I);
log=5*log(1+I2);
log=intrans(I,'log',5);
gamma gamma=imadjust(I,[],[],0.4); gamma=intrans(I,'gamma',0.4);
contrast-stretching I2=im2double(I);
contrast=1./(1+(0.2./(I2+eps)).^5);
contrast=intrans(I,'stretch',0.2,5);

 


4. References


5. Exercise

Part 1

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 -> M-File

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.jpg');
    imshow(lifting);

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

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

  4. % test the error checking
    lifting_sub=returnSubImage(lifting,159,33,500)

    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.

      Deliverables:

      • Histogram of original image
      • Histogram of image produced from histogram equalization
      • Image produced from histogram matching
      • Histogram of image produced from histogram-matching