modules/features2d/src/blobdetector.cpp

   1 /*M///////////////////////////////////////////////////////////////////////////////////////
   2 //
   3 //  IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
   4 //
   5 //  By downloading, copying, installing or using the software you agree to this license.
   6 //  If you do not agree to this license, do not download, install,
   7 //  copy or use the software.
   8 //
   9 //
  10 //                           License Agreement
  11 //                For Open Source Computer Vision Library
  12 //
  13 // Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
  14 // Copyright (C) 2009, Willow Garage Inc., all rights reserved.
  15 // Third party copyrights are property of their respective owners.
  16 //
  17 // Redistribution and use in source and binary forms, with or without modification,
  18 // are permitted provided that the following conditions are met:
  19 //
  20 //   * Redistribution's of source code must retain the above copyright notice,
  21 //     this list of conditions and the following disclaimer.
  22 //
  23 //   * Redistribution's in binary form must reproduce the above copyright notice,
  24 //     this list of conditions and the following disclaimer in the documentation
  25 //     and/or other materials provided with the distribution.
  26 //
  27 //   * The name of the copyright holders may not be used to endorse or promote products
  28 //     derived from this software without specific prior written permission.
  29 //
  30 // This software is provided by the copyright holders and contributors "as is" and
  31 // any express or implied warranties, including, but not limited to, the implied
  32 // warranties of merchantability and fitness for a particular purpose are disclaimed.
  33 // In no event shall the Intel Corporation or contributors be liable for any direct,
  34 // indirect, incidental, special, exemplary, or consequential damages
  35 // (including, but not limited to, procurement of substitute goods or services;
  36 // loss of use, data, or profits; or business interruption) however caused
  37 // and on any theory of liability, whether in contract, strict liability,
  38 // or tort (including negligence or otherwise) arising in any way out of
  39 // the use of this software, even if advised of the possibility of such damage.
  40 //
  41 //M*/
  42
  43 #include "precomp.hpp"
  44 #include <iterator>
  45
  46 //#define DEBUG_BLOB_DETECTOR
  47
  48 #ifdef DEBUG_BLOB_DETECTOR
  49 #  include "opencv2/opencv_modules.hpp"
  50 #  ifdef HAVE_OPENCV_HIGHGUI
  51 #    include "opencv2/highgui/highgui.hpp"
  52 #  else
  53 #    undef DEBUG_BLOB_DETECTOR
  54 #  endif
  55 #endif
  56
  57 using namespace cv;
  58
  59 /*
  60 *  SimpleBlobDetector
  61 */
  62 SimpleBlobDetector::Params::Params()
  63 {
  64     thresholdStep = 10;
  65     minThreshold = 50;
  66     maxThreshold = 220;
  67     minRepeatability = 2;
  68     minDistBetweenBlobs = 10;
  69
  70     filterByColor = true;
  71     blobColor = 0;
  72
  73     filterByArea = true;
  74     minArea = 25;
  75     maxArea = 5000;
  76
  77     filterByCircularity = false;
  78     minCircularity = 0.8f;
  79     maxCircularity = std::numeric_limits<float>::max();
  80
  81     filterByInertia = true;
  82     //minInertiaRatio = 0.6;
  83     minInertiaRatio = 0.1f;
  84     maxInertiaRatio = std::numeric_limits<float>::max();
  85
  86     filterByConvexity = true;
  87     //minConvexity = 0.8;
  88     minConvexity = 0.95f;
  89     maxConvexity = std::numeric_limits<float>::max();
  90 }
  91
  92 void SimpleBlobDetector::Params::read(const cv::FileNode& fn )
  93 {
  94     thresholdStep = fn["thresholdStep"];
  95     minThreshold = fn["minThreshold"];
  96     maxThreshold = fn["maxThreshold"];
  97
  98     minRepeatability = (size_t)(int)fn["minRepeatability"];
  99     minDistBetweenBlobs = fn["minDistBetweenBlobs"];
 100
 101     filterByColor = (int)fn["filterByColor"] != 0 ? true : false;
 102     blobColor = (uchar)(int)fn["blobColor"];
 103
 104     filterByArea = (int)fn["filterByArea"] != 0 ? true : false;
 105     minArea = fn["minArea"];
 106     maxArea = fn["maxArea"];
 107
 108     filterByCircularity = (int)fn["filterByCircularity"] != 0 ? true : false;
 109     minCircularity = fn["minCircularity"];
 110     maxCircularity = fn["maxCircularity"];
 111
 112     filterByInertia = (int)fn["filterByInertia"] != 0 ? true : false;
 113     minInertiaRatio = fn["minInertiaRatio"];
 114     maxInertiaRatio = fn["maxInertiaRatio"];
 115
 116     filterByConvexity = (int)fn["filterByConvexity"] != 0 ? true : false;
 117     minConvexity = fn["minConvexity"];
 118     maxConvexity = fn["maxConvexity"];
 119 }
 120
 121 void SimpleBlobDetector::Params::write(cv::FileStorage& fs) const
 122 {
 123     fs << "thresholdStep" << thresholdStep;
 124     fs << "minThreshold" << minThreshold;
 125     fs << "maxThreshold" << maxThreshold;
 126
 127     fs << "minRepeatability" << (int)minRepeatability;
 128     fs << "minDistBetweenBlobs" << minDistBetweenBlobs;
 129
 130     fs << "filterByColor" << (int)filterByColor;
 131     fs << "blobColor" << (int)blobColor;
 132
 133     fs << "filterByArea" << (int)filterByArea;
 134     fs << "minArea" << minArea;
 135     fs << "maxArea" << maxArea;
 136
 137     fs << "filterByCircularity" << (int)filterByCircularity;
 138     fs << "minCircularity" << minCircularity;
 139     fs << "maxCircularity" << maxCircularity;
 140
 141     fs << "filterByInertia" << (int)filterByInertia;
 142     fs << "minInertiaRatio" << minInertiaRatio;
 143     fs << "maxInertiaRatio" << maxInertiaRatio;
 144
 145     fs << "filterByConvexity" << (int)filterByConvexity;
 146     fs << "minConvexity" << minConvexity;
 147     fs << "maxConvexity" << maxConvexity;
 148 }
 149
 150 SimpleBlobDetector::SimpleBlobDetector(const SimpleBlobDetector::Params &parameters) :
 151 params(parameters)
 152 {
 153 }
 154
 155 void SimpleBlobDetector::read( const cv::FileNode& fn )
 156 {
 157     params.read(fn);
 158 }
 159
 160 void SimpleBlobDetector::write( cv::FileStorage& fs ) const
 161 {
 162     params.write(fs);
 163 }
 164
 165 void SimpleBlobDetector::findBlobs(const cv::Mat &image, const cv::Mat &binaryImage, vector<Center> &centers) const
 166 {
 167     (void)image;
 168     centers.clear();
 169
 170     vector < vector<Point> > contours;
 171     Mat tmpBinaryImage = binaryImage.clone();
 172     findContours(tmpBinaryImage, contours, CV_RETR_LIST, CV_CHAIN_APPROX_NONE);
 173
 174 #ifdef DEBUG_BLOB_DETECTOR
 175     //  Mat keypointsImage;
 176     //  cvtColor( binaryImage, keypointsImage, CV_GRAY2RGB );
 177     //
 178     //  Mat contoursImage;
 179     //  cvtColor( binaryImage, contoursImage, CV_GRAY2RGB );
 180     //  drawContours( contoursImage, contours, -1, Scalar(0,255,0) );
 181     //  imshow("contours", contoursImage );
 182 #endif
 183
 184     for (size_t contourIdx = 0; contourIdx < contours.size(); contourIdx++)
 185     {
 186         Center center;
 187         center.confidence = 1;
 188         Moments moms = moments(Mat(contours[contourIdx]));
 189         if (params.filterByArea)
 190         {
 191             double area = moms.m00;
 192             if (area < params.minArea || area >= params.maxArea)
 193                 continue;
 194         }
 195
 196         if (params.filterByCircularity)
 197         {
 198             double area = moms.m00;
 199             double perimeter = arcLength(Mat(contours[contourIdx]), true);
 200             double ratio = 4 * CV_PI * area / (perimeter * perimeter);
 201             if (ratio < params.minCircularity || ratio >= params.maxCircularity)
 202                 continue;
 203         }
 204
 205         if (params.filterByInertia)
 206         {
 207             double denominator = sqrt(pow(2 * moms.mu11, 2) + pow(moms.mu20 - moms.mu02, 2));
 208             const double eps = 1e-2;
 209             double ratio;
 210             if (denominator > eps)
 211             {
 212                 double cosmin = (moms.mu20 - moms.mu02) / denominator;
 213                 double sinmin = 2 * moms.mu11 / denominator;
 214                 double cosmax = -cosmin;
 215                 double sinmax = -sinmin;
 216
 217                 double imin = 0.5 * (moms.mu20 + moms.mu02) - 0.5 * (moms.mu20 - moms.mu02) * cosmin - moms.mu11 * sinmin;
 218                 double imax = 0.5 * (moms.mu20 + moms.mu02) - 0.5 * (moms.mu20 - moms.mu02) * cosmax - moms.mu11 * sinmax;
 219                 ratio = imin / imax;
 220             }
 221             else
 222             {
 223                 ratio = 1;
 224             }
 225
 226             if (ratio < params.minInertiaRatio || ratio >= params.maxInertiaRatio)
 227                 continue;
 228
 229             center.confidence = ratio * ratio;
 230         }
 231
 232         if (params.filterByConvexity)
 233         {
 234             vector < Point > hull;
 235             convexHull(Mat(contours[contourIdx]), hull);
 236             double area = contourArea(Mat(contours[contourIdx]));
 237             double hullArea = contourArea(Mat(hull));
 238             double ratio = area / hullArea;
 239             if (ratio < params.minConvexity || ratio >= params.maxConvexity)
 240                 continue;
 241         }
 242
 243         center.location = Point2d(moms.m10 / moms.m00, moms.m01 / moms.m00);
 244
 245         if (params.filterByColor)
 246         {
 247             if (binaryImage.at<uchar> (cvRound(center.location.y), cvRound(center.location.x)) != params.blobColor)
 248                 continue;
 249         }
 250
 251         //compute blob radius
 252         {
 253             vector<double> dists;
 254             for (size_t pointIdx = 0; pointIdx < contours[contourIdx].size(); pointIdx++)
 255             {
 256                 Point2d pt = contours[contourIdx][pointIdx];
 257                 dists.push_back(norm(center.location - pt));
 258             }
 259             std::sort(dists.begin(), dists.end());
 260             center.radius = (dists[(dists.size() - 1) / 2] + dists[dists.size() / 2]) / 2.;
 261         }
 262
 263         centers.push_back(center);
 264
 265 #ifdef DEBUG_BLOB_DETECTOR
 266         //    circle( keypointsImage, center.location, 1, Scalar(0,0,255), 1 );
 267 #endif
 268     }
 269 #ifdef DEBUG_BLOB_DETECTOR
 270     //  imshow("bk", keypointsImage );
 271     //  waitKey();
 272 #endif
 273 }
 274
 275 void SimpleBlobDetector::detectImpl(const cv::Mat& image, std::vector<cv::KeyPoint>& keypoints, const cv::Mat&) const
 276 {
 277     //TODO: support mask
 278     keypoints.clear();
 279     Mat grayscaleImage;
 280     if (image.channels() == 3)
 281         cvtColor(image, grayscaleImage, CV_BGR2GRAY);
 282     else
 283         grayscaleImage = image;
 284
 285     vector < vector<Center> > centers;
 286     for (double thresh = params.minThreshold; thresh < params.maxThreshold; thresh += params.thresholdStep)
 287     {
 288         Mat binarizedImage;
 289         threshold(grayscaleImage, binarizedImage, thresh, 255, THRESH_BINARY);
 290
 291 #ifdef DEBUG_BLOB_DETECTOR
 292         //    Mat keypointsImage;
 293         //    cvtColor( binarizedImage, keypointsImage, CV_GRAY2RGB );
 294 #endif
 295
 296         vector < Center > curCenters;
 297         findBlobs(grayscaleImage, binarizedImage, curCenters);
 298         vector < vector<Center> > newCenters;
 299         for (size_t i = 0; i < curCenters.size(); i++)
 300         {
 301 #ifdef DEBUG_BLOB_DETECTOR
 302             //      circle(keypointsImage, curCenters[i].location, curCenters[i].radius, Scalar(0,0,255),-1);
 303 #endif
 304
 305             bool isNew = true;
 306             for (size_t j = 0; j < centers.size(); j++)
 307             {
 308                 double dist = norm(centers[j][ centers[j].size() / 2 ].location - curCenters[i].location);
 309                 isNew = dist >= params.minDistBetweenBlobs && dist >= centers[j][ centers[j].size() / 2 ].radius && dist >= curCenters[i].radius;
 310                 if (!isNew)
 311                 {
 312                     centers[j].push_back(curCenters[i]);
 313
 314                     size_t k = centers[j].size() - 1;
 315                     while( k > 0 && centers[j][k].radius < centers[j][k-1].radius )
 316                     {
 317                         centers[j][k] = centers[j][k-1];
 318                         k--;
 319                     }
 320                     centers[j][k] = curCenters[i];
 321
 322                     break;
 323                 }
 324             }
 325             if (isNew)
 326             {
 327                 newCenters.push_back(vector<Center> (1, curCenters[i]));
 328                 //centers.push_back(vector<Center> (1, curCenters[i]));
 329             }
 330         }
 331         std::copy(newCenters.begin(), newCenters.end(), std::back_inserter(centers));
 332
 333 #ifdef DEBUG_BLOB_DETECTOR
 334         //    imshow("binarized", keypointsImage );
 335         //waitKey();
 336 #endif
 337     }
 338
 339     for (size_t i = 0; i < centers.size(); i++)
 340     {
 341         if (centers[i].size() < params.minRepeatability)
 342             continue;
 343         Point2d sumPoint(0, 0);
 344         double normalizer = 0;
 345         for (size_t j = 0; j < centers[i].size(); j++)
 346         {
 347             sumPoint += centers[i][j].confidence * centers[i][j].location;
 348             normalizer += centers[i][j].confidence;
 349         }
 350         sumPoint *= (1. / normalizer);
 351         KeyPoint kpt(sumPoint, (float)(centers[i][centers[i].size() / 2].radius));
 352         keypoints.push_back(kpt);
 353     }
 354
 355 #ifdef DEBUG_BLOB_DETECTOR
 356     namedWindow("keypoints", CV_WINDOW_NORMAL);
 357     Mat outImg = image.clone();
 358     for(size_t i=0; i<keypoints.size(); i++)
 359     {
 360         circle(outImg, keypoints[i].pt, keypoints[i].size, Scalar(255, 0, 255), -1);
 361     }
 362     //drawKeypoints(image, keypoints, outImg);
 363     imshow("keypoints", outImg);
 364     waitKey();
 365 #endif
 366 }