斑点检测
2015-08-15 08:43
232 查看
参考:http://blog.csdn.net/abcd1992719g/article/details/27071273
http://www.cnblogs.com/ronny/p/3895883.html#3036526
这两位博主都实现了这个算法,里面有很多细节,还是很有收获的。
View Code
http://www.cnblogs.com/ronny/p/3895883.html#3036526
这两位博主都实现了这个算法,里面有很多细节,还是很有收获的。
#include "opencv2/core/core.hpp" #include "opencv2/features2d/features2d.hpp" #include "opencv2/highgui/highgui.hpp" #include "opencv2/calib3d/calib3d.hpp" #include "opencv2/nonfree/nonfree.hpp" #include<opencv2/highgui/highgui.hpp> #include<opencv2/imgproc/imgproc.hpp> #include <iostream> using namespace cv; using namespace std; /******************************Blob特征实现**********************************************/ /*****生成log sigma******/ vector<double> create_sigma(double start, double step, double end) { vector<double> sigma; while (start <= end + 1e-8) { double s = exp(start); sigma.push_back(s); start += step; } return sigma; } /********生成log算子模板************/ Mat create_LOG(int size, double sigma) { Mat H(size, size, CV_64F); int cx = (size - 1) / 2; int cy = (size - 1) / 2; double sum = 0; for (int i = 0; i < H.rows; i++) { for (int j = 0; j < H.cols; j++) { int nx = i - cx; int ny = j - cy; double s = -1 / (3.1415926 * sigma*sigma*sigma*sigma); s = s* (1 - (nx*nx + ny*ny) / (2 * sigma*sigma)); s = s*exp(-(nx*nx + ny*ny) / (2 * sigma*sigma)); sum += s; H.at<double>(i, j) = s; } } double mean = sum / (size*size); for (int i = 0; i < H.rows; i++) { for (int j = 0; j < H.cols; j++) H.at<double>(i, j) -= mean; } return H; } struct blob { int x; int y; double sigma;//sigma用于计算半径r=1.5*sigma int intensity;//intensity用于当两个 blob重叠度高时取intensity大的(梯度大) }; /*****生成LOG卷积后的scale_space******/ vector<Mat> creat_scale_space(Mat& im_in, vector<double>& sigma) { vector<Mat> scale_space; for (int i = 0; i < sigma.size(); i++) { int n = ceil(sigma[i] * 3) * 2 + 1; Mat LOG = create_LOG(n, sigma[i]); Mat dst; //为了使得response一致,我们对做卷积的图片要先乘上sigma^2,这一步叫做scale normalize filter2D(im_in.mul(sigma[i] * sigma[i]), dst, -1, LOG, Point(-1, -1)); scale_space.push_back(dst); } return scale_space; } /*******检测所有可能的blob**********/ vector<struct blob> detect_blobs(Mat& im_in, double thresh, vector<double> et, int padding = 10) { vector<struct blob> blobs; Mat addborder, norm; /*The function copies the source image into the middle of the destination image. The areas to the left, to the right, above and below the copied source image will be filled with extrapolated pixels*/ copyMakeBorder(im_in, addborder, padding, padding, padding, padding, BORDER_CONSTANT, Scalar(255)); normalize(addborder, norm, 1, 0, NORM_MINMAX, CV_64F); vector<Mat> all_ims = creat_scale_space(norm, et); for (int i = 0; i < et.size(); i++) { Mat grayscale, mx; normalize(all_ims[i], grayscale, 0, 255, NORM_MINMAX, CV_8U, Mat()); Mat structedElement(3, 3, CV_8U, Scalar(1)); dilate(grayscale, mx, structedElement);//默认就是3*3,第三个参数可省略 grayscale = (grayscale == mx)&(all_ims[i]>thresh);//取大于threshold并且是周围最大 for (int j = 0; j < norm.rows; j++) { for (int k = 0; k < norm.cols; k++) { if (grayscale.at<double>(j, k)>0)//是局部最大值点 { struct blob b; b.x = j - padding; b.y = k - padding; b.sigma = et[i]; b.intensity = all_ims[i].at<double>(j, k); blobs.push_back(b); } } } } return blobs; } /*****************删除重叠度大的blob********************/ vector<struct blob> prune_blobs(vector<struct blob>blobs_in) { vector<bool> keep(blobs_in.size(), true); for (int i = 0; i < blobs_in.size(); i++) { for (int j = 0; j < blobs_in.size(); j++) { if ((i == j) || (keep[i] == false) || (keep[j] == false))//?????????? continue; struct blob blobi = blobs_in[i]; struct blob blobj = blobs_in[j]; int xi = blobi.x; int yi = blobi.y; double ri = blobi.sigma; int xj = blobj.x; int yj = blobj.y; double rj = blobj.sigma; //算重叠度 double d = sqrt((xj - xi)*(xj - xi) + (yj - yi)*(yj - yi)); double rirj = ri + rj; double overlap_percent = rirj / d; if (overlap_percent > 2.0) { if (blobi.intensity > blobj.intensity) { keep[i] = true; keep[j] = false; } else { keep[i] = false; keep[j] = true; } } } } vector<struct blob> blobs_out; for (int i = 0; i < blobs_in.size(); i++) if (keep[i]) blobs_out.push_back(blobs_in[i]); return blobs_out; } int main() { Mat src, gray; src = imread("1.jpg"); cvtColor(src, gray, CV_RGB2GRAY); vector<double> sigma = create_sigma(1.0, 0.2, 3.0); vector<struct blob> blobs = detect_blobs(gray, 0.2, sigma); vector<struct blob> blobs_trimmed = prune_blobs(blobs); for (int i = 0; i < blobs_trimmed.size(); i++) { struct blob b = blobs_trimmed[i]; circle(src, Point(b.y, b.x), 1.5*b.sigma, Scalar(0), 2, 8, 0); } namedWindow("result", 1); imshow("result", src); waitKey(0); return 0; }
View Code
相关文章推荐
- spring mvc 控制器方法传递一些经验对象的数组
- POJ 3268 Silver Cow Party(最短路dijkstra)
- LeetCode Power of Two
- LeetCode Power of Two
- hdu5358
- hdu5358
- 玩转Node.js - 03. 第一个I/O!
- PHP的Undefined variable错误怎么解决?
- block详解
- UVA 624 CD (DP)
- UVA 10954
- HDU 1014 Uniform Generator
- ie自带打印
- PHP面试题基础问题
- OC语言-06-OC语言-block与protocol
- iOS 浅赋值、深复制、完全复制的知识点梳理验证(附加归档解档)
- 【iOS】随机三角瓷砖布局算法
- BZOJ1029
- js的隐含参数(arguments,callee,caller)使用方法
- 最短路算法之dijkstra