实现Harris角点检测器，按照标准材料，最终检测到边缘

问题描述 投票：0回答：1

嗨，感谢您的关注，我正在尝试实现哈里斯角点检测算法，并且一直遵循大多数来源给出的标准算法，其中之一是Wiki页面我遵循该算法的所有过程，并且我使用sobel 滤波器用于计算梯度 gx 和 gy，但结果最终是一个边缘检测器（最终比我实现的精明边缘检测器好得多：O）我想知道我做错了什么。

核心.h

#ifndef CORE
#define CORE

#include <opencv2/core.hpp>
#include <opencv2/highgui.hpp>
#include <opencv2/imgcodecs.hpp>
#include <opencv2/imgproc.hpp>
#include <opencv2/opencv.hpp>
using namespace std;
using namespace cv;

#endif

Utils.h

#ifndef UTILS
#define UTILS

#include "Core.h"

static void displayImage(Mat &img, unsigned int time = 0,
                         string title = "frame") {
  imshow(title, img);
  waitKey(time);
}

static bool isValidPoint(Mat &img, int x, int y) {
  int rows = img.rows;
  int cols = img.cols;
  return (x >= 0 and x < rows and y >= 0 and y < cols);
}

static void performDiff(Mat img1, Mat img2) {
  int rows = img1.rows;
  int cols = img2.cols;
  Mat diff(rows, cols, CV_8U);
  bool zero = true;
  for (int i = 0; i < rows; i++) {
    for (int j = 0; j < cols; j++) {
      int d = (static_cast<int>(img1.at<u_char>(i, j)) -
               static_cast<int>(img2.at<u_char>(i, j)));
      diff.at<u_char>(i, j) = d;
      if (d != 0)
        zero = false;
    }
  }
  displayImage(diff, 0, "diff image");
}

static int getValIntOf(Mat &img, int x, int y) {
  return static_cast<int>(img.at<u_char>(x, y));
}

#endif

高斯变换.h

#ifndef GUASSIAN_TRANSFORM
#define GUASSIAN_TRANSFORM


static vector<vector<double>> getGuassianKernal(int size, double sigma) {
  vector<vector<double>> Kernal(size, vector<double>(size, 0.0));
  double sum = 0.0;
  int center = size / 2;
  for (int i = 0; i < size; i++) {
    for (int j = 0; j < size; j++) {
      int x = i - center;
      int y = j - center;
      Kernal[i][j] = exp(-(x * x + y * y) / (2 * sigma * sigma));
      sum += Kernal[i][j];
    }
  }
  if (sum != 0) {
    for (int i = 0; i < size; i++) {
      for (int j = 0; j < size; j++) {
        Kernal[i][j] /= sum;
      }
    }
  }
  return Kernal;
}


static Mat guassianFilterTransform(Mat &img, int FiltSize = 3,
                                   double Sigma = 1.5) {
  vector<vector<double>> filter = getGuassianKernal(FiltSize, Sigma);
  int trows = img.rows - FiltSize + 1;
  int tcols = img.cols - FiltSize + 1;
  Mat transformed(trows, tcols, CV_8U);
  for (int i = 1; i <= trows - 2; i++) {
    for (int j = 1; j <= tcols - 2; j++) {
      double tval = 0;
      for (int x = -1; x <= 1; x++) {
        for (int y = -1; y <= 1; y++) {
          tval = tval + (filter[x + 1][y + 1] *
                         static_cast<double>(img.at<u_char>(i + x, j + y)));
        }
      }
      transformed.at<u_char>(i, j) = static_cast<u_char>(tval);
    }
  }
  return (transformed);
}

#endif

哈里斯.cpp

#ifndef HARRIS_CORNER
#define HARRIS_CORNER

#include "../../BasicTransforms/GaussianTransform.h"
#include "../../Utils/Core/Core.h"
#include "../../Utils/Core/Utils.h"

double computeCornerResponseFunc(vector<vector<double>> &StructMat) {
  double det =
      StructMat[0][0] * StructMat[1][1] - StructMat[0][1] * StructMat[1][0];
  double trace = StructMat[0][0] + StructMat[1][1];
  return (det + (0.04) * (trace * trace));
}

Mat detectHarrisCorners(Mat &img) {
  int filtSize = 3;
  vector<vector<int>> filtery, filterx;
  filtery = {{-1, -2, -1}, {0, 0, 0}, {1, 2, 1}};
  filterx = {{-1, 0, 1}, {-2, 0, 2}, {-1, 0, 1}};
  int trows = img.rows - filtSize + 1;
  int tcols = img.cols - filtSize + 1;
  Mat gradMag(trows, tcols, CV_8U);
  Mat gradDir(trows, tcols, CV_64F);
  Mat result(trows, tcols, CV_8U);
  double Thrs = 200;
  vector<vector<double>> gaussianKernal = getGuassianKernal(3, 1.5);
  for (int i = 1; i <= trows - 2; i++) {
    for (int j = 1; j <= tcols - 2; j++) {
      double gradX = 0;
      double gradY = 0;
      for (int x = -1; x <= 1; x++) {
        for (int y = -1; y <= 1; y++) {
          gradX =
              gradX + gaussianKernal[x + 1][y + 1] *
                          (filterx[x + 1][y + 1] *
                           static_cast<double>(img.at<u_char>(i + x, j + y)));
          gradY =
              gradY + gaussianKernal[x + 1][y + 1] *
                          (filtery[x + 1][y + 1] *
                           static_cast<double>(img.at<u_char>(i + x, j + y)));
        }
      }
      vector<vector<double>> strtMtrx(2, vector<double>(2, 0.0));
      strtMtrx[0][0] = gradX * gradX;
      strtMtrx[1][1] = gradY * gradY;
      strtMtrx[0][1] = strtMtrx[1][0] = gradX * gradY;
      double CornerResponseFunc = computeCornerResponseFunc(strtMtrx);
      if (CornerResponseFunc >= Thrs) {
        result.at<u_char>(i, j) = 255;
      }
      double tval = sqrt(gradX * gradX + gradY * gradY);
      if (gradX == 0.00) {
        gradDir.at<double>(i, j) = (gradY >= 0 ? M_PI / 2.0 : -M_PI / 2.0);
      } else {
        gradDir.at<double>(i, j) = atan(gradY / gradX);
      }
      gradDir.at<double>(i, j) *= (180.0 / M_PI);
      gradDir.at<double>(i, j) = abs(gradDir.at<double>(i, j));
      gradMag.at<u_char>(i, j) = static_cast<u_char>(tval);
    }
  }
  for (int i = 1; i < trows - 1; i++) {
    for (int j = 1; j < tcols - 1; j++) {
      double angle = gradDir.at<double>(i, j);
      double neigh1, neigh2;
      if ((angle >= 0 and angle < 22.5) or (angle >= 157.5 && angle < 180)) {
        neigh1 = gradMag.at<u_char>(i, j + 1);
        neigh2 = gradMag.at<u_char>(i, j - 1);
      } else if ((angle >= 22.5 and angle < 67.5)) {
        neigh1 = gradMag.at<u_char>(i - 1, j + 1);
        neigh2 = gradMag.at<u_char>(i + 1, j - 1);
      } else if ((angle >= 67.5 and angle < 112.5)) {
        neigh1 = gradMag.at<u_char>(i - 1, j);
        neigh2 = gradMag.at<u_char>(i + 1, j);
      } else if ((angle >= 112.5 and angle < 157.5)) {
        neigh1 = gradMag.at<u_char>(i - 1, j - 1);
        neigh2 = gradMag.at<u_char>(i + 1, j + 1);
      } else {
        neigh1 = neigh2 = 0;
      }
      if (gradMag.at<u_char>(i, j) >= max(neigh1, neigh2))
        result.at<u_char>(i, j) = result.at<u_char>(i, j);
      else
        result.at<u_char>(i, j) = 0;
    }
  }
  return result;
}

int main() {
  string imPath =
      "/home/panirpal/workspace/Projects/ComputerVision/data/chess2.jpg";
  Mat img = imread(imPath, IMREAD_GRAYSCALE);
  if (!img.empty()) {
    Mat gaussian = guassianFilterTransform(img);
    Mat out = detectHarrisCorners(gaussian);
    displayImage(out);
  } else
    cerr << "image not found! exiting...";
}
#endif

结果：

预期输出：

尝试在 

computeCornerResponseFunc
中尝试使用诸如阈值和 k 值 0.04 之类的神奇数字，但最终检测到边缘却没有用。

c++ opencv image-processing edge-detection corner-detection
1个回答
1
投票

您的代码存在一些问题。最重要的是，您没有对平方梯度图像应用平滑。

          gradX =
              gradX + gaussianKernal[x + 1][y + 1] *
                          (filterx[x + 1][y + 1] *
                           static_cast<double>(img.at<u_char>(i + x, j + y)));

在这里，您正在计算图像与由 Sobel 核和高斯核的元素相乘形成的核的卷积。这只是意味着您没有使用 Sobel 内核，而是完全不同的东西，我无法立即分析。

逐元素乘法和卷积是两种不同的东西，并且它们不能交换。结构张量 M （如您提供的链接中所述）的第一个元素是 G*(IxIx)。 Ix 是 Sobel 滤波器的结果。 IxIx 是它的平方（结果的元素级乘法，就像你所做的那样）。但G*部分是与高斯的卷积。这种平滑必须在元素乘法之后应用。这与 (G*Ix)(G*Ix) 不同！

还有一个问题：

static_cast<u_char>(tval)

这会将双精度 

tval
转换为 8 位无符号整数。 
tval
保证为非负数，但它很可能大于 255。如果是这样，转换将换行，值 256 将被写为 0。您必须先限制该值：

static_cast<u_char>(std::min(tval,255.0))


performDiff
中的这个更糟糕：

      int d = (static_cast<int>(img1.at<u_char>(i, j)) -
               static_cast<int>(img2.at<u_char>(i, j)));
      diff.at<u_char>(i, j) = d

减法的结果不能用8位无符号整数表示。当 

img2
较大时，
diff
将具有一些正值，例如 -1 将表示为 255。您需要取绝对差，并进行与我上面所示相同的钳位。

调用 

getGuassianKernal(3, 1.5)
不会返回高斯分布。用 sigma 1.5 到 3 个像素截断高斯，留下的东西比高斯更类似于具有均匀权重的过滤器。我在here更详细地解释了这一点。对于 1.5 的西格玛，我们通常使用大小为 
2*ceil(3*1.5)+1
的内核，即 11。至少您需要 
2*ceil(2*1.5)+1 = 7
（不太精确，但速度更快）。

最新问题
© www.soinside.com 2019 - 2024. All rights reserved.