BASIC IMAGE PROCESSING TUTORIAL

Tutorials on standard use of MATLAB can be found at here while the specific use for image processing is here.
 

Fourier Transform

To compute the Fourier Transform of an image do the following

     % Prepare image
     f = zeros(30,30);
     f(5:24,13:17) = 1;
     imshow(f,'notruesize')

     % Compute Fourier Transform
     F = fft2(f,256,256);
     F = fftshift(F); % Center FFT

     % Measure the minimum and maximum value of the transform amplitude
     min(min(abs(F))) %   0
     max(max(abs(F))) % 100
     imshow(abs(F),[0,100]); colormap(jet); colorbar
     imshow(log(1+abs(F)),[0,3]); colormap(jet); colorbar
            * What is the main difference between representing the amplitude and its logarithm?

     % Look at the phases
     imshow(angle(F),[-pi,pi]); colormap(jet); colorbar

            * Try with other images
     f = imread('saturn.tif');
     f = ind2gray(f,gray(256));

            * Or more patterned from Yamitani (at Caltech Univ.) homepage.
            * Checkered images
     function I=checkered(N,s)
     % This function returns a checkered image of size N
     % and with a check size of s

     I=zeros(N,N);
     for i=1:N,
         for j=1:N,
             if (mod(floor(i/s)+floor(j/s),2)==0) I(i,j)=1; end;
         end;
     end;

     f=checkered(256,32); imshow(f);
 

Sampling rate (Downsampling)

Downsampling can be done in two fashions: by interpolation and by sampling theory.

By interpolation (explore different possibilities of interpolation schemes used by imresize)

     f=imread('blood1.tif');
     f=ind2gray(f,gray(256));
     f=f(1:256,1:256);
     f1=imresize(f,0.2);      % Downsample to a 1/5th of the size
     f2=imresize(f1,5);       % Go back to original size
     figure; imshow(f);
     figure; imshow(f1);
     figure; imshow(f2);

Via sampling theory

     f=imread('blood1.tif');
     f=ind2gray(f,gray(256));
     f=f(1:256,1:256);

     function Idown=downsampling(I,m)
     % Downsample the square image I by a factor of m

     [N,M]=size(I);

     % Apply ideal filter
     w=1/m;
     F=fftshift(fft2(I));
     for i=1:N
         for j=1:N
             r2=(i-round(N/2))^2+(j-round(N/2))^2;
             if (r2>round((N/2*w)^2)) F(i,j)=0; end;
         end;
     end;
     Idown=real(ifft2(fftshift(F)));

     % Now downsample
     Idown=imresize(Idown,[N/m,N/m],'nearest');

     % ====================================================================
     function Iup=upsampling(I,m)
     % Upsample the square image I by a factor of m

     [N,M]=size(I);
     Iup=zeros(m*N,m*N);

     % Expand input image
     for i=1:N
         for j=1:N
            Iup(m*(i-1)+1,m*(j-1)+1)=I(i,j);
         end;
     end;

     % Ideal filter
     [N,M]=size(Iup);
     w=1/m;
     F=fftshift(fft2(Iup));
     for i=1:N
         for j=1:N
             r2=(i-round(N/2))^2+(j-round(N/2))^2;
             if (r2>round((N/2*w)^2)) F(i,j)=0; end;
         end;
     end;
     Iup=(m*m)*abs(ifft2(fftshift(F)));

     f1=downsampling(f,2);
     f2=upsampling(f1,2);
     figure; imshow(f);
     figure; imshow(f1);
     figure; imshow(f2);
 

Image Quantization

     f=load('saturn.tif');
     imshow(f); colormap(gray);
     imshow(f); colormap(gray(32));
     imshow(f); colormap(gray(4));
     imshow(f); colormap(jet(16)); % False color
 

Fourier space filtering

     function Ifilter=ideal_lowpass(I,w)
     % ideally filter image I up to a frequency w

     [N,M]=size(I);
     F=fftshift(fft2(I));
     for i=1:N
         for j=1:N
             r2=(i-round(N/2))^2+(j-round(N/2))^2;
             if (r2>round((N*w)^2)) F(i,j)=0; end;
         end;
     end;
     Ifilter=real(ifft2(fftshift(F)));

     f=imread('blood1.tif');
     f=ind2gray(f,gray(256));
     f=f(1:256,1:256);
     figure; imshow(ideal_lowpass(f,0.25));
     figure; imshow(ideal_lowpass(f,0.10));

Real space filtering

Next example shows you how to perform a low pass filter via convolution with a small kernel

     f=imread('blood1.tif');
     f=ind2gray(f,gray(256));
     H=0.25*[0.25 0.5 0.25; 0.5 1 0.5; 0.25 0.5 0.25];
     f1=filter2(H,f);
     figure; imshow(f);
     figure; imshow(f1);

The frequency respose of the filter can be shown with

     freqz2(H);
 

Filter design

Histogram operations

     imhist(f);

while imadjust and histeq changes the histogram shape or limits

     f1=imadjust(f,[0.4 0.6],[0 1]);
     f2=histeq(f);
     figure; imhist(f);
     figure; imhist(f1);
     figure; imhist(f2);
     figure; imshow(f);
     figure; imshow(f1);
     figure; imshow(f2);
 

Image enhancement

• Thresholding


f=imread('blood1.tif');
f=ind2gray(f,gray(256));
imhist(f);
f2=max(f1,0.4); % remove all values below 0.4
imshow(f);
figure; imshow(f2);
 

• Median filtering for salt & pepper noise


f1=imnoise(f,'salt & pepper',0.05);
f2=medfilt2(f1,[3 3]);
figure; imshow(f);
figure; imshow(f1);
figure; imshow(f2);
 

• Adaptative filtering for speckle noise

 

 

f1=imnoise(f,'speckle',0.04);
f2=wiener2(f1,[7 7]);
figure; imshow(f);
figure; imshow(f1);
figure; imshow(f2);
 

• Edge reinforcement

 

 

f=imread('saturn.tif');
f=in2gray(f,gray(256));
f1=filter2(fspecial('unsharp'),f);
figure; imshow(f);
figure; imshow(f1);

Edge detection

Segmentation

1. Thresholding

All values above some threshold are set to 1 and those below the threshold are set to 0 or viceversa.

Morphological operations

Image processing practices by other people