Aqua Phoenix
     >>  Lectures >>  Matlab 10  


10.2 Signal Processing - Digital Image processing

Digital Image processing refers to the manipulation of image data, which is motivated by conversion between the spatial and frequency domains. This section focuses on examples in smoothing, sharpening, and edge detection. Most of these signal processing techniques require the application of filters to an image. In the spatial domain, a filter is a matrix, typically 3x3 or 5x5, which is applied to input pixels, one at a time, resulting in an output pixel in the filtered image. Any operation in the spatial domain corresponds to an operation that can be described in the frequency domain, and hence filters tend to be explained in both domains.

Images can be represented in several color spaces - the most common are grayscale and RGB (red, green, blue), where RGB can be thought of as 3-layers of intensity representing red, green, and blue components.

Matlab offers a fairly simple function that can read images of many popular file formats: matrix = imread(filename); Depending on file type and color space, the returned matrix is either a 2D matrix of intensity values (grayscale images), or index values (indexed images, like GIF), or a 3D matrix of RGB values. Red, green, and blue intensities, when mixed together, result in the full color spectrum of all representable colors.

pomegranate = imread('pomegranateSeeds.jpg');
subplot(2,2,1), image(pomegranate), title('Tasty Pomegranate Seeds (RGB)');
subplot(2,2,2), imshow(pomegranate(:,:,1)), title('Intensity Image: Red Layer');
subplot(2,2,3), imshow(pomegranate(:,:,2)), title('Intensity Image: Green Layer');
subplot(2,2,4), imshow(pomegranate(:,:,3)), title('Intensity Image: Blue Layer');
An example of a color photograph of a bowl with pomegranate seeds. The matrix representing this image is a 3D structure, holding data for rows, columns, and color component. Shown in the figure are the intensity values for all 3 color layers.

Function image interprets a 3D matrix as an RGB color image and displays its content.

Function imshow displays a 2D matrix as an intensity image (grayscale or indexed images).

Figure 10.12
Click image to enlarge, or click here to open
Images are represented with a specific kind of data type, which allows for non-negative, full integer intensity values between 0 and 255. The corresponding data type is named uint8 (unsigned integer with 8 bit representation).

Matrices loaded from images using imread are automatically typed as uint8. When manipulating or building image data matrices in Matlab, data type double is typically used. The use of different data types requires a facility for conversion, which are available with the following type conversion functions:

  • uint8(x) converts matrix x to uint8 data type (0-255)
  • double(x) converts matrix x to double data type
  • logical(x) converts matrix x to logical data type (0s and 1s)

Before displaying an image data matrix, the data should be converted to uint8:

image(uint8(x)) displays double-type matrix x.

10.2.1 Filters

Common 2D filters can be built in Matlab by using built-in function fspecial (special filters).

fspecial(filtername, paramaters, ...) = matrix of values representing the filter.

  • fspecial defines the following common filters:
  • average : averaging filter
  • disk : circular averaging filter
  • gaussian : Gaussian lowpass filter
  • laplacian : filter approximating the 2-D Laplacian operator
  • log : Laplacian of Gaussian filter
  • motion : motion filter
  • prewitt : Prewitt horizontal edge-emphasizing filter
  • sobel : Sobel horizontal edge-emphasizing filter
  • unsharp : unsharp contrast enhancement filter

Filters are applied to multi-dimensional images (RGB images = 3D matrices) using function imfilter , and are applied to intensity images (2D matrices) using function filter2 .

10.2.2 Smoothing / Blurring

The process of smoothing or blurring and image supresses noise and small fluctuations. In the frequency domain, this process refers to the supression of high frequencies.

A smoothing filter can be built in Matlab by using function fspecial (special filters):

gaussianFilter = fspecial('gaussian', [7, 7], 5)
builds a gaussian filter matrix of 7 rows and 7 columns, with standard deviation of 5.

Figure 10.13
Click image to enlarge, or click here to open
An example of the application of a Gaussian filter to the multi-dimensional image of pomegranate seeds:

gaussianPomegranate = imfilter(pomegranate, gaussianFilter, 'symmetric', 'conv');
subplot(1,2,1), image(pomegranate);
subplot(1,2,2), image(gaussianPomegranate), title('Blurred Pomegranate, blur matrix size 7');
Figure 10.14
Click image to enlarge, or click here to open
An example of the application of the same Gaussian filter to an intensity image (we take the red layer of the previous image):

gaussianPomegranate = filter2(gaussianFilter, pomegranate(:,:,1));
subplot(1,2,1), imshow(pomegranate(:,:,1));
subplot(1,2,2), imshow(uint8(gaussianPomegranate)), title('Blurred Pomegranate, blur matrix size 7');
Figure 10.15
Click image to enlarge, or click here to open

10.2.3 Edge Detection

The process of edge detection attenuates high fluctuations in color, i.e. dramatic change in intensity. In the frequency domain, this process refers to the attenuation of high frequencies.

Matlab includes the built-in function edge designed for edge detection. It supports the following types of edge detectors:

  • sobel
  • prewitt
  • roberts
  • log (laplacian)
  • canny
  • zerocross

For demonstrative purposes, we scale down the image of pomegranate seeds:

smallPomegranate = pomegranate(1:4:size(pomegranate,1), 1:4:size(pomegranate,2), :);
subplot(2,2,1), image(smallPomegranate);
Figure 10.16
Click image to enlarge, or click here to open
We now extract edges for each color layer (red, green, and blue), and place the edge intensity images in a matrix which will resemble a multi-dimensional color (RGB) image:

edges = [];
edges(:,:,1) = edge(smallPomegranate(:,:,1), 'sobel');
edges(:,:,2) = edge(smallPomegranate(:,:,2), 'sobel');
edges(:,:,3) = edge(smallPomegranate(:,:,3), 'sobel');
subplot(2,2,2), imshow(edges(:,:,1)), title('Pomegranate Edges Channel R');
subplot(2,2,3), imshow(edges(:,:,2)), title('Pomegranate Edges Channel G');
subplot(2,2,4), imshow(edges(:,:,3)), title('Pomegranate Edges Channel B');
Figure 10.17
Click image to enlarge, or click here to open
Figure 10.18
Click image to enlarge, or click here to open
We can also show the composite (color) edge images, which combine information from all three channels.

subplot(1,2,1), image(edges), title('Pomegranate Edges RGB');
edgesAll = edges(:,:,1) + edges(:,:,2) + edges(:,:,3);
subplot(1,2,2), imshow(edgesAll), title('Combined edges in one intensity image');
The figure shows a color edge image and a combined intensity image.

Figure 10.19
Click image to enlarge, or click here to open

10.2.4 Sharpening

The process of sharpening is related to edge detection - changes in color are attenuated to create an effect of sharper edges.

Using a fspecial , we create a filter for sharpening an image. The special filter is ironically named 'unsharp':

sharpFilter = fspecial('unsharp');
subplot(2,2,1), image(pomegranate), title('Original Pomegranate Seeds');
sharp = imfilter(pomegranate, sharpFilter, 'replicate');
subplot(2,2,2), image(sharp), title('Sharpened Pomegranate');
sharpMore = imfilter(sharp, sharpFilter, 'replicate');
subplot(2,2,3), image(sharpMore), title('Excessive sharpening attenuates noise');
Figure 10.20
Click image to enlarge, or click here to open
Figure 10.21
Click image to enlarge, or click here to open

10.2.5 Filter and Sharpening in more detail

The process of sharpening an image is quite interesting and ironic: In order to sharpen an image, it is first blurred, edges are detected in the blurred version, and finally added to the blurred image to create a sharper image.

The original image:

subplot(2,2,1), image(pomegranate), title('Original Pomegranate');
To demonstrate filters in more detail, we will create a blur (Gaussian) filter manually. The following values are used.

gaussianFilter = [1,4,7,4,1;4,20,33,20,4;7,33,55,33,7;4,20,33,20,4;1,4,7,4,1]
Values in a Gaussian filter are used as weights to mix a given input pixel and its neighboring pixels to create an output pixel which has been "smudged" with its neighborhood. Accordingly, the center value is the largest, corresponding to an input pixel. Values decrease towards the edge of the filter matrix, corresponding to neighboring pixels.

The matrix must be normalized such that the sum of values is equal to one. If the sum of weights is smaller or greater than one, the brightness of the image is decreased or increased.

gaussianFilter = gaussianFilter / sum(sum(gaussianFilter));
Convolve image using the filter:

gaussianPomegranate = imfilter(pomegranate, gaussianFilter);
subplot(2,2,2), image(gaussianPomegranate), title('Gaussian Pomegranate');
Create a filter for edge detection. We will use a Laplacian filter. Note that the sum of values is equal to zero, which means that overall brightess is not preserved. In fact, the resulting image is mostly black with only a few lines denoting edges.

edgeFilter = [-1,-1,-1;-1,8,-1;-1,-1,-1]
edgePomegranate = imfilter(gaussianPomegranate, edgeFilter);
subplot(2,2,3), image(edgePomegranate), title('Edges from Blurred Pomegranate');
Finally, we add the edge image to the blurred image to create the final sharp version:

sharpPomegranate = gaussianPomegranate + edgePomegranate;
subplot(2,2,4), image(sharpPomegranate), title('Sharp Pomegranate (sharper than original)');
Figure 10.22
Click image to enlarge, or click here to open
Figure 10.23
Click image to enlarge, or click here to open

10.2.6 Line detection

Line detection is a special type of edge detection. In edge detection, a pixel is attenuated, if there is a dramatic change in color in any direction. For line detection, the direction in which a color change is considered is restricted.

The following tables outline a Laplacian for edge detection, and their line detection counterparts.

Edge detection Horizonal line detection Vertical line detection 45 degree (upward slope) 45 degree (downward slope)
Table 10.1: Laplacian edge and line detectors
horizontalKernel = [-1,-1,-1;2,2,2;-1,-1,-1];
verticalKernel = [-1,2,-1;-1,2,-1;-1,2,-1];
diagUpKernel = [-1,-1,2;-1,2,-1;2,-1,-1];
diagDownKernel = [2,-1,-1;-1,2,-1;-1,-1,2];
image(building), title('Building');
Figure 10.24
Click image to enlarge, or click here to open
horizontalBuilding = imfilter(building, horizontalKernel);
subplot(2,2,1), image(horizontalBuilding), title('Horizontal lines');
verticalBuilding = imfilter(building, verticalKernel);
subplot(2,2,2), image(verticalBuilding), title('Vertical lines');
diagDownBuilding = imfilter(building, diagDownKernel);
subplot(2,2,3), image(diagDownBuilding), title('Diagonals down slope');
diagUpBuilding = imfilter(building, diagUpKernel);
subplot(2,2,4), image(diagUpBuilding), title('Diagonals up slope');
Applying these filters to an image containing differently slanted edges reveals differently emphasized edge types:

Figure 10.25
Click image to enlarge, or click here to open