简单来说,Canny检测的步骤如下:
(1)使用低通滤波器抑制噪声
(2)计算梯度幅值和梯度方向
(3)非极大值抑制
(4)检测边缘
下面给出每一步的代码。
抑制噪声:
Mat imageSource = imread("D:\\image.jpg", 0), imageGaussian; //读取图像 //imshow("Image", imageSource); GaussianBlur(imageSource, imageGaussian, Size(5, 5), 0); //高斯滤波平滑图像 //imshow("Gaussian Image", imageGaussian);
计算X和Y梯度幅值,梯度方向:
/** imageInput - 输入图像 imageX - 输出X方向的梯度幅值 imageY - 输出Y方向的梯度幅值 pointDirection - 每一点的梯度方向 **/ void SobelXY(Mat imageInput, Mat &imageX, Mat &imageY, double *pointDirection) { imageX = imageInput.clone(); imageY = imageInput.clone(); for (int i = 0; i < imageInput.rows*imageInput.cols; i++) { pointDirection[i] = 0; } double dx, dy; int k = 0; for (int i = 1; i < imageInput.rows - 1; i++) { for (int j = 1; j < imageInput.cols - 1; j++) { dx = imageInput.at<uchar>(i + 1, j - 1) + 2 * imageInput.at<uchar>(i + 1, j) + imageInput.at<uchar>(i + 1, j + 1) - (imageInput.at<uchar>(i - 1, j - 1) + 2 * imageInput.at<uchar>(i - 1, j) + imageInput.at<uchar>(i - 1, j + 1)); dy = imageInput.at<uchar>(i - 1, j - 1) + 2 * imageInput.at<uchar>(i, j - 1) + imageInput.at<uchar>(i + 1, j - 1) - (imageInput.at<uchar>(i - 1, j + 1) + 2 * imageInput.at<uchar>(i, j + 1) + imageInput.at<uchar>(i + 1, j + 1)); imageX.at<uchar>(i, j) = saturate_cast<uchar>(dx); imageY.at<uchar>(i, j) = saturate_cast<uchar>(dy); if (dx == 0) dx = 0.00000000000000001; //避免除0 pointDirection[k] = atan(dy / dx) * 57.3; //弧度转角度 pointDirection[k] += 90; k++; } } convertScaleAbs(imageX, imageX); convertScaleAbs(imageY, imageY); }
梯度幅值:
若梯度幅值为G,X方向为Gx,Y方向为Gy。则有G = sqrt(Gx^2 + Gy^2)。也有简单的算法G = Gx + Gy
/** imageX - X方向的梯度幅值 imageY - Y方向的梯度幅值 imageOutput - 输出图像 **/ void SobelXYAdd(Mat imageX, Mat imageY, Mat &imageOutput) { imageOutput = imageX.clone(); for (int i = 0; i < imageX.rows; i++) { uchar *ptrX = imageX.ptr<uchar>(i); uchar *ptrY = imageY.ptr<uchar>(i); uchar *ptr = imageOutput.ptr<uchar>(i); for (int j = 0; j < imageX.cols; j++) { ptr[j] = saturate_cast<uchar>(sqrt(ptrX[j] * ptrX[j] + ptrY[j] * ptrY[j])); } } convertScaleAbs(imageOutput, imageOutput); }
非极大值抑制:沿着该点的梯度方向,比较前后两个点的梯度幅值,若该点大于前后两点的梯度幅值则保留,否则置为0
此处我贴个图用以说明前后两点究竟是哪两点以及其梯度幅值的求法
上图中g1,g2,g3,g4,c都为像素点,其中点c是我们所要判断是否为极大值的点。显然,g1,g2,g3,g4在点c相邻的八个点。手绘的点A,B即是点c的前后两点。设点g1,g2的梯度幅值为val_g1,val_g2
则点A的梯度幅值val_A = m*val_g1+(1-m)val_g2,其中m = tan(θ),显然θ可由该点的梯度方向求得
即此处采用简单的线性插值来求出点A的梯度幅值。类似的B也可这样求得。下面给出代码:
1 /** 2 imageInput - 输入图像 3 imageOutput - 输出图像 4 pointDirection - 各个点的梯度方向 5 **/ 6 void NotMaxValue(Mat imageInput, Mat &imageOutput, double *pointDirection) 7 { 8 imageOutput = imageInput.clone(); 9 int k = 0; 10 for (int i = 1; i < imageInput.rows - 1; i++) 11 { 12 for (int j = 1; j < imageInput.cols - 1; j++) 13 { 14 uchar leftTop = imageInput.at<uchar>(i - 1, j - 1); 15 uchar top = imageInput.at<uchar>(i - 1, j); 16 uchar rightTop = imageInput.at<uchar>(i - 1, j + 1); 17 uchar left = imageInput.at<uchar>(i, j - 1); 18 uchar val = imageInput.at<uchar>(i, j); 19 uchar right = imageInput.at<uchar>(i, j + 1); 20 uchar leftBottom = imageInput.at<uchar>(i + 1, j - 1); 21 uchar bottom = imageInput.at<uchar>(i + 1, j); 22 uchar rightBottom = imageInput.at<uchar>(i + 1, j + 1); 23 int direction = pointDirection[k]; 24 if (direction > 0 && direction <= 45) 25 { 26 if (val <= right + tan(direction)*(rightTop - right) || 27 val <= left + tan(direction)*(leftBottom - left)) 28 { 29 imageInput.at<uchar>(i, j) = 0; 30 } 31 } 32 if (direction > 45 && direction <= 90) 33 { 34 if (val <= top + tan(90 - direction)*(rightTop - top) || 35 val <= bottom + tan(90 - direction)*(leftBottom - bottom)) 36 { 37 imageInput.at<uchar>(i, j) = 0; 38 } 39 } 40 if (direction > 90 && direction <= 135) 41 { 42 if (val <= top + tan(direction - 90)*(leftTop - top) || 43 val <= bottom + tan(direction - 90)*(rightBottom - bottom)) 44 { 45 imageInput.at<uchar>(i, j) = 0; 46 } 47 } 48 if (direction > 135 && direction <= 180) 49 { 50 if (val <= left + tan(180 - direction)*(leftTop - left) || 51 val <= right + tan(180 - direction)*(rightBottom - right)) 52 { 53 imageInput.at<uchar>(i, j) = 0; 54 } 55 } 56 k++; 57 } 58 } 59 }
边缘检测:
这里采用了双阈值,有阈值L,H。其中L < H
若某点的梯度幅值大于H,则认为该点是强边缘,若小于L则认为该点非边缘。若在两个阈值之间,则认为是弱边缘,根据其与强边缘是否连通来判断其是否是边缘。
void DoubleThreshold(Mat &imageInput, int lowThreshold, int highThreshold) { for (int i = 0; i < imageInput.rows; i++) { uchar *ptr = imageInput.ptr<uchar>(i); for (int j = 0; j < imageInput.cols; j++) { if (ptr[j] < lowThreshold) ptr[j] = 0; else if (ptr[j] >= highThreshold) ptr[j] = 255; } } } bool LinkWithMaxValue(Mat image, int i, int j) { uchar val = image.at<uchar>(i, j); uchar leftTop = image.at<uchar>(i - 1, j - 1); uchar top = image.at<uchar>(i - 1, j); uchar rightTop = image.at<uchar>(i - 1, j + 1); uchar left = image.at<uchar>(i, j - 1); uchar right = image.at<uchar>(i, j + 1); uchar leftBottom = image.at<uchar>(i + 1, j - 1); uchar bottom = image.at<uchar>(i + 1, j); uchar rightBottom = image.at<uchar>(i + 1, j + 1); return (leftTop == 255 || top == 255 || rightTop == 255 || left == 255 || right == 255 || leftBottom == 255 || bottom == 255 || rightBottom == 255); } void DoubleThresholdLink(Mat &imageOutput) { int dx[] = { -1, 0, 1 }; int dy[] = { -1, 0, 1 }; bool *visited = new bool[imageOutput.rows*imageOutput.cols]; stack<Point> mst; for (int i = 0; i < imageOutput.rows*imageOutput.cols; i++) visited[i] = 0; for (int i = 1; i < imageOutput.rows - 1; i++) { for (int j = 1; j < imageOutput.cols - 1; j++) { if (visited[i*imageOutput.rows + j]) continue; visited[i*imageOutput.rows + j] = 1; uchar val = imageOutput.at<uchar>(i, j); if (!val || val == 255) continue; if (LinkWithMaxValue(imageOutput, i, j)) { imageOutput.at<uchar>(i, j) = 255; while (!mst.empty()) { Point temp = mst.top(); mst.pop(); if (visited[temp.x*imageOutput.rows + temp.y]) continue; visited[temp.x*imageOutput.rows*temp.x + temp.y] = 1; int data = imageOutput.at<uchar>(temp.x, temp.y); if (!data) continue; imageOutput.at<uchar>(temp.x, temp.y) = 255; for (int r = 0; r < 3; r++) { for (int c = 0; c < 3; c++) { uchar changedX = temp.x + dx[r]; uchar changedY = temp.y + dy[c]; if (visited[changedX*imageOutput.rows + changedY]) continue; mst.push(Point(changedX, changedY)); } } } } else imageOutput.at<uchar>(i, j) = 0; } } }
最后我贴下整个程序的代码:
1 #include "core/core.hpp" 2 #include "highgui/highgui.hpp" 3 #include "imgproc/imgproc.hpp" 4 #include "iostream" 5 #include "math.h" 6 #include "stack" 7 8 using namespace std; 9 using namespace cv; 10 11 void SobelXY(Mat imageInput, Mat &imageX, Mat &imageY, double *pointDirection); 12 13 void SobelXYAdd(Mat imageX, Mat imageY, Mat &imageOutput); 14 15 void NotMaxValue(Mat imageInput, Mat &imageOutput, double *pointDirection); 16 17 bool LinkWithMaxValue(Mat image, int i, int j); 18 19 void DoubleThreshold(Mat &imageInput, int lowThreshold, int highThreshold); 20 21 void DoubleThresholdLink(Mat &imageOutput); 22 23 int main(int argc, char **argv) 24 { 25 Mat imageSource = imread("D:\\image.jpg", 0), imageGaussian; //读取图像 26 //imshow("Image", imageSource); 27 GaussianBlur(imageSource, imageGaussian, Size(5, 5), 0); //高斯滤波平滑图像 28 //imshow("Gaussian Image", imageGaussian); 29 Mat imageSobelX, imageSobelY, imageSobel; 30 double *pointDirection = new double[imageGaussian.rows * imageGaussian.cols]; 31 SobelXY(imageGaussian, imageSobelX, imageSobelY, pointDirection); //计算X和Y方向的梯度和梯度方向 32 //imshow("Sobel X", imageSobelX); 33 //imshow("Sobel Y", imageSobelY); 34 SobelXYAdd(imageSobelX, imageSobelY, imageSobel); //计算梯度幅值 35 //imshow("Sobel", imageSobel); 36 Mat imageMaxValue; 37 NotMaxValue(imageSobel, imageMaxValue, pointDirection); //抑制非极大值 38 //imshow("MaxValue", imageMaxValue); 39 DoubleThreshold(imageMaxValue, 50, 100); 40 DoubleThresholdLink(imageMaxValue); 41 imshow("Canny", imageMaxValue); 42 Canny(imageGaussian, imageGaussian, 50, 100); 43 imshow("OpenCV Canny", imageGaussian); 44 waitKey(); 45 return 0; 46 } 47 48 /** 49 imageInput - 输入图像 50 imageX - 输出X方向的梯度幅值 51 imageY - 输出Y方向的梯度幅值 52 pointDirection - 每一点的梯度方向 53 **/ 54 void SobelXY(Mat imageInput, Mat &imageX, Mat &imageY, double *pointDirection) 55 { 56 imageX = imageInput.clone(); 57 imageY = imageInput.clone(); 58 for (int i = 0; i < imageInput.rows*imageInput.cols; i++) 59 { 60 pointDirection[i] = 0; 61 } 62 double dx, dy; 63 int k = 0; 64 for (int i = 1; i < imageInput.rows - 1; i++) 65 { 66 for (int j = 1; j < imageInput.cols - 1; j++) 67 { 68 dx = imageInput.at<uchar>(i + 1, j - 1) + 2 * imageInput.at<uchar>(i + 1, j) + imageInput.at<uchar>(i + 1, j + 1) - 69 (imageInput.at<uchar>(i - 1, j - 1) + 2 * imageInput.at<uchar>(i - 1, j) + imageInput.at<uchar>(i - 1, j + 1)); 70 71 dy = imageInput.at<uchar>(i - 1, j - 1) + 2 * imageInput.at<uchar>(i, j - 1) + imageInput.at<uchar>(i + 1, j - 1) - 72 (imageInput.at<uchar>(i - 1, j + 1) + 2 * imageInput.at<uchar>(i, j + 1) + imageInput.at<uchar>(i + 1, j + 1)); 73 74 imageX.at<uchar>(i, j) = saturate_cast<uchar>(dx); 75 imageY.at<uchar>(i, j) = saturate_cast<uchar>(dy); 76 if (dx == 0) dx = 0.00000000000000001; //避免除0 77 78 pointDirection[k] = atan(dy / dx) * 57.3; //弧度转角度 79 pointDirection[k] += 90; 80 k++; 81 } 82 } 83 convertScaleAbs(imageX, imageX); 84 convertScaleAbs(imageY, imageY); 85 } 86 /** 87 imageX - X方向的梯度幅值 88 imageY - Y方向的梯度幅值 89 imageOutput - 输出图像 90 **/ 91 void SobelXYAdd(Mat imageX, Mat imageY, Mat &imageOutput) 92 { 93 imageOutput = imageX.clone(); 94 for (int i = 0; i < imageX.rows; i++) 95 { 96 uchar *ptrX = imageX.ptr<uchar>(i); 97 uchar *ptrY = imageY.ptr<uchar>(i); 98 uchar *ptr = imageOutput.ptr<uchar>(i); 99 for (int j = 0; j < imageX.cols; j++) 100 { 101 ptr[j] = saturate_cast<uchar>(sqrt(ptrX[j] * ptrX[j] + ptrY[j] * ptrY[j])); 102 } 103 } 104 convertScaleAbs(imageOutput, imageOutput); 105 } 106 107 /** 108 imageInput - 输入图像 109 imageOutput - 输出图像 110 pointDirection - 各个点的梯度方向 111 **/ 112 void NotMaxValue(Mat imageInput, Mat &imageOutput, double *pointDirection) 113 { 114 imageOutput = imageInput.clone(); 115 int k = 0; 116 for (int i = 1; i < imageInput.rows - 1; i++) 117 { 118 for (int j = 1; j < imageInput.cols - 1; j++) 119 { 120 uchar leftTop = imageInput.at<uchar>(i - 1, j - 1); 121 uchar top = imageInput.at<uchar>(i - 1, j); 122 uchar rightTop = imageInput.at<uchar>(i - 1, j + 1); 123 uchar left = imageInput.at<uchar>(i, j - 1); 124 uchar val = imageInput.at<uchar>(i, j); 125 uchar right = imageInput.at<uchar>(i, j + 1); 126 uchar leftBottom = imageInput.at<uchar>(i + 1, j - 1); 127 uchar bottom = imageInput.at<uchar>(i + 1, j); 128 uchar rightBottom = imageInput.at<uchar>(i + 1, j + 1); 129 int direction = pointDirection[k]; 130 if (direction > 0 && direction <= 45) 131 { 132 if (val <= right + tan(direction)*(rightTop - right) || 133 val <= left + tan(direction)*(leftBottom - left)) 134 { 135 imageInput.at<uchar>(i, j) = 0; 136 } 137 } 138 if (direction > 45 && direction <= 90) 139 { 140 if (val <= top + tan(90 - direction)*(rightTop - top) || 141 val <= bottom + tan(90 - direction)*(leftBottom - bottom)) 142 { 143 imageInput.at<uchar>(i, j) = 0; 144 } 145 } 146 if (direction > 90 && direction <= 135) 147 { 148 if (val <= top + tan(direction - 90)*(leftTop - top) || 149 val <= bottom + tan(direction - 90)*(rightBottom - bottom)) 150 { 151 imageInput.at<uchar>(i, j) = 0; 152 } 153 } 154 if (direction > 135 && direction <= 180) 155 { 156 if (val <= left + tan(180 - direction)*(leftTop - left) || 157 val <= right + tan(180 - direction)*(rightBottom - right)) 158 { 159 imageInput.at<uchar>(i, j) = 0; 160 } 161 } 162 k++; 163 } 164 } 165 } 166 167 void DoubleThreshold(Mat &imageInput, int lowThreshold, int highThreshold) 168 { 169 for (int i = 0; i < imageInput.rows; i++) 170 { 171 uchar *ptr = imageInput.ptr<uchar>(i); 172 for (int j = 0; j < imageInput.cols; j++) 173 { 174 if (ptr[j] < lowThreshold) ptr[j] = 0; 175 else if (ptr[j] >= highThreshold) ptr[j] = 255; 176 } 177 } 178 } 179 180 bool LinkWithMaxValue(Mat image, int i, int j) 181 { 182 uchar val = image.at<uchar>(i, j); 183 uchar leftTop = image.at<uchar>(i - 1, j - 1); 184 uchar top = image.at<uchar>(i - 1, j); 185 uchar rightTop = image.at<uchar>(i - 1, j + 1); 186 uchar left = image.at<uchar>(i, j - 1); 187 uchar right = image.at<uchar>(i, j + 1); 188 uchar leftBottom = image.at<uchar>(i + 1, j - 1); 189 uchar bottom = image.at<uchar>(i + 1, j); 190 uchar rightBottom = image.at<uchar>(i + 1, j + 1); 191 return (leftTop == 255 || top == 255 || rightTop == 255 || left == 255 || 192 right == 255 || leftBottom == 255 || bottom == 255 || rightBottom == 255); 193 } 194 195 void DoubleThresholdLink(Mat &imageOutput) 196 { 197 int dx[] = { -1, 0, 1 }; 198 int dy[] = { -1, 0, 1 }; 199 bool *visited = new bool[imageOutput.rows*imageOutput.cols]; 200 stack<Point> mst; 201 for (int i = 0; i < imageOutput.rows*imageOutput.cols; i++) visited[i] = 0; 202 for (int i = 1; i < imageOutput.rows - 1; i++) 203 { 204 for (int j = 1; j < imageOutput.cols - 1; j++) 205 { 206 if (visited[i*imageOutput.rows + j]) continue; 207 visited[i*imageOutput.rows + j] = 1; 208 uchar val = imageOutput.at<uchar>(i, j); 209 if (!val || val == 255) continue; 210 if (LinkWithMaxValue(imageOutput, i, j)) 211 { 212 imageOutput.at<uchar>(i, j) = 255; 213 while (!mst.empty()) 214 { 215 Point temp = mst.top(); mst.pop(); 216 if (visited[temp.x*imageOutput.rows + temp.y]) continue; 217 visited[temp.x*imageOutput.rows*temp.x + temp.y] = 1; 218 int data = imageOutput.at<uchar>(temp.x, temp.y); 219 if (!data) continue; 220 imageOutput.at<uchar>(temp.x, temp.y) = 255; 221 for (int r = 0; r < 3; r++) 222 { 223 for (int c = 0; c < 3; c++) 224 { 225 uchar changedX = temp.x + dx[r]; 226 uchar changedY = temp.y + dy[c]; 227 if (visited[changedX*imageOutput.rows + changedY]) continue; 228 mst.push(Point(changedX, changedY)); 229 } 230 } 231 } 232 } 233 else imageOutput.at<uchar>(i, j) = 0; 234 } 235 } 236 }View Code
参考博客:
博客1:https://blog.csdn.net/dcrmg/article/details/52344902
博客2:https://www.cnblogs.com/love6tao/p/5152020.html
其中,文中所用图片来自博客2