Implementation of image edge detection by using the canny operator (detailed process with source code)

As the first-order differential of Gaussian function, the well balance between noise suppression and edge detection can be achieved.

Environment: windows xp + matlab 2010b

Time:

The main steps to detect the edge of the canny algorithm are as follows:

1) filter with a 3x3 Gaussian filter to eliminate noise;

2) Calculate the differential approximation between the horizontal and vertical directions for each pixel to obtain the gradient size and direction of the pixel;

3) perform "non-maximum suppression" on the gradient (the non-local maximum value is set to 0 );

4) obtain the gradient threshold twice;

5) connect the edge;

The steps are described in detail below:

0) read the image:

Clear;
Clc;
I =imread('light.jpg ');
K = rgb2y (I); % obtain the h component, that is, the brightness component.

According to the definition of edge,Edge DetectionThe purpose is to mark the point where the brightness of a digital image changes significantly. (Refer to Wikipedia http://zh.wikipedia.org/wiki/%E8%BE%B9%E7%BC%98%E6%A3%80%E6%B5%8B );

 

According to the formula Brightness = 0.3 * R + 0.6 * G + 0.1 * B, the Brightness component y is calculated;

Function k = rgb2y (z)
% I must be an rgb three-dimensional matrix
[M, n, p] = size (z );
K = zeros (m, n );
Z = double (z );
For I = 1: m
For j = 1: n
K (I, j) = 0.3 * z (I, j, 1) + 0.6 * z (I, j, 2) + 0.1 * z (I, j, 3 );
End
End

 

1) filter with a 3x3 Gaussian filter to eliminate noise;

Function j = gaosi (I );
% I must be a two-dimensional double Matrix
J = I;
[H, w] = size (I );
For m = 2; H-1
For n = 2: W-1
J (m, n) = (I (m, n-1) + 2 * I (m, n) + I (m, n + 1)/4; % transverse Gaussian filter
End
End

The preceding custom gaosi function is used for two-dimensional 3x3 image filtering,

K1 = gaosi (k); % Horizontal Filtering
K1 = k1 '; % transpose the image to prepare for the next vertical filter (vertical filter = horizontal filter after transpose)
K1 = gaosi (k1 );
K1 = k1'; % restore

 

2) Calculate the differential approximation between the horizontal and vertical directions for each pixel to obtain the gradient size and direction of the pixel;

The magnitude and direction of the gradient are calculated based on the horizontal and vertical differential approximation obtained by the formula above, yi zhi p and q respectively.

% Calculate the gradient size and direction
[H, w] = size (k );
For m = 2: H-1
For n = 2: W-1
Zz1 = k1 (m, n-1) + k1 (m + 1, n-1 );
Zz2 = k1 (m, n) + k1 (m + 1, n );
Zz3 = k1 (MB, n-1) + k1 (m, n );
Zz4 = k1 (m + 1, n-1) + k1 (m + 1, n );
Kp (m, n) = 0.5*(zz2-zz1 );
Kq (m, n) = 0.5*(zz3-zz4 );
Kfu (m, n) = sqrt (kp (m, n) ^ 2) + (kq (m, n) ^ 2); % gradient size
Kjiao (m, n) = atan (kq (m, n)/(kp (m, n) + 0.001); % Gradient Direction, 0.001 to prevent the denominator from being 0
End
End

3) perform "non-maximum suppression" on the gradient (the non-local maximum value is set to 0 );

1. classify the gradient direction into four main directions: left and right, top and bottom, left oblique, and right oblique.

% Non-maximum suppression
% First divide the gradient direction into four directions: 90,135, (and their reverse extension line)
For m = 2: H-1
For n = 2: W-1
If kjiao (m, n) >=3/8 * pi
Kjiao (m, n) = 2;
Else if kjiao (m, n) >=1/8 * pi
Kjiao (m, n) = 1;
Else if kjiao (m, n) >=- 1/8 * pi
Kjiao (m, n) = 0;
Else if kjiao (m, n) >=- 3/8 * pi
Kjiao (m, n) = 3;
Else
Kjiao (m, n) = 2;
End
End
End
End
End
End

Determine whether the vertex is the local maximum value (Gradient Direction) in the eight-neighbor region based on the four directions after the division. For example, if the gradient direction is a point in the left and right directions, determine whether it is greater than the value of the Left or Right two points. If not, set the value of this point to 0.

% Determine by direction
K2 = k1;
For m = 2: H-1
For N = 2: W-1
If kjiao (m, n) = 0
If K1 (m, n)> K1 (M, n-1) & K1 (m, n)> K1 (m, n + 1 );
Else K2 (m, n) = 0;
End
End
If kjiao (m, n) = 1
If K1 (m, n)> K1 (m + 1, n-1) & K1 (m, n)> K1 (m-1, n + 1 );
Else K2 (m, n) = 0;
End
End
If kjiao (m, n) = 2
If K1 (m, n)> K1 (m-1, n) & K1 (m, n)> K1 (m + 1, n );
Else K2 (m, n) = 0;
End
End
If kjiao (m, n) = 3
If K1 (m, n)> K1 (S-1, n-1) & K1 (m, n)> K1 (m + 1, n + 1 );
Else K2 (m, n) = 0;
End
End
End
End

4) obtain the gradient threshold twice;

Using two thresholds t1 and t2 (T2> T1, generally t2 = 2 * T1), we set the gray scale of the pixel whose gradient value is smaller than T1 to 0 to get image 1, then, we set the gray scale of the pixel whose gradient value is less than t2 to 0 to get Image 2. Because Image 2 has a high threshold and less noise (but it also loses useful edge information, while image 1 has a low threshold and retains more information, therefore, we can connect the edge of an Image Based on Image 2 and supplemented with image 1.

% Two Threshold Segmentation
K3 = k2; % matrix separated by T1 as the threshold
K4 = k2; % matrix after T2 is used as the threshold
T1 = 50;
T2 = 2 * T1;
For m = 2: H-1
For N = 2: W-1
If KFU (m, n) <t1
K3 (m, n) = 0;
End
If KFU (m, n) <t2
K4 (m, n) = 0;
End
End
End

5) connect the edge;

A. Scan Image 2. When we encounter a non-zero pixel P (tracking the contour line starting with P until the end of the contour line Q;
B. in image 1, we examine the eight neighborhoods of point P corresponding to point P in Image 2. If there is a non-zero pixel Q' in the eight neighborhoods of point P, include it in Image 2 as point R, starting from R (Repeat step a until we can't continue in Image 1 and Image 2;
C. we have ended the connection to the contour line containing P, marked it as accessed, and returned to step a to find the next contour line, repeating step (a) (B) (c) step 2 until the new contour line is no longer found in Image 2.

Findline. M:

Function [ff, flag1] = findline (K3, K4, flag, m, n)

Flag1 = flag;
M1 = m + 1; n1 = n + 1;
While (M ~ = M1 | n ~ = N1) % if neither m nor n changes, line has reached the end point
Flagg = 0;
For I = 1:3
If (Flagg = 1) break;
End
For J = 1:3
If K3 (m-2 + I, N-2 + J )~ = 0
K4 (m-2 + I, N-2 + J) = 255;
M1 = m-2 + I; n1 = n-2 + J; % new [M, N] Point
Flag1 (m, n) = 1; % mark Detected
Flagg = 1; break;
End
End
End
M = m1; n = n1;
End
Ff = k4;

Write in the main function:

Figure;
Subplot (221); imshow (I); title ('original image ');
Subplot (222); imshow (k3, []); title ('image segmentation with a threshold of 50 ');
Subplot (223); imshow (k4, []); title ('image segmentation with a threshold of 100 ');
Flag = zeros (h, w); % indicates whether the vertex has been detected. 1 indicates that the vertex has been detected.
For m = 2: H-1
For n = 2: W-1
If k4 (m, n )~ = 0 & flag (m, n) = 0
[K4, flag] = findline (k3, k4, flag, m, n );
End
End
End
Subplot (224); imshow (k4, []); title ('corrected segmented image ');

So far, the program has been completed.

The effect is as follows:

 

You can define the threshold value in the program.

I hope our predecessors can learn image processing experience. Thank you.

 

Reprinted Please note: http://www.cnblogs.com/blue-lg/archive/2011/12/25/2300899.html