OpenCV找圆

#include "opencv2/imgproc/imgproc.hpp"
#include "opencv2/highgui/highgui.hpp"
#include <stdlib.h>
#include <stdio.h>
#include <iostream>
using namespace cv;
using namespace std;
void LeastSquaresCircleFitting(vector<Point2d> Points, Point2d &Centroid, double& dRadius);

int main()

{
  Mat src_gray, src_binary, src_blur, detected_edges;
  Mat src = imread("月亮.bmp", 1);
  cvtColor(src, src_gray, CV_BGR2GRAY);
  medianBlur(src_gray, src_blur, 21);
  threshold(src_blur, src_binary, 50, 255, 0);
  imwrite("src_binary.jpg", src_binary);
  Canny(src_binary, detected_edges, 50, 150, 3);
  imwrite("筋片detected_edges.jpg", detected_edges);
  //----------------------------------------------------------------//
  vector<Point2d> m_Points; 
  //扫描取点，从左到右，逐行扫描 
  for (int i = 100; i <1900; ++i)
  {
    uchar* data = detected_edges.ptr<uchar>(i);
    for (int j = 1570; j < 3000; ++j)
    {
      if (data[j] == 255)
      {
        m_Points.push_back(Point(j, i));
        break;
      }
    }
  }

  //扫描取点，从右到左，逐行扫描 
  for (int i = 100; i < 1900; ++i)
  {
    uchar* data = detected_edges.ptr<uchar>(i);
    for (int j = 1570; j > 50; --j)
    {
      if (data[j] == 255)
      {
        m_Points.push_back(Point(j, i));
        break;
      }
    }
  }

  //扫描后的点保存在图像scan_dst中
  Mat scan_dst = Mat::zeros(cv::Size(src_gray.cols, src_gray.rows), src_gray.type());
  vector<Point2d>::iterator iter1;
  vector<Point2d>::iterator end1 = m_Points.end();
  for (iter1 = m_Points.begin(); iter1 != end1; ++iter1)
  {
    scan_dst.at<uchar>((*iter1).x, (*iter1).y) = 255;
  }
  //显示扫描后的图
  //namedWindow("Result1", 0); imshow("Result1", scan_dst);
  imwrite("scan_dst.jpg", scan_dst);


  //  //求半径
  double m_dRadius_DianPian = 0;
  Point2d m_Centroid_DianPian;
  LeastSquaresCircleFitting(m_Points, m_Centroid_DianPian, m_dRadius_DianPian);
  /最小二乘拟合圆
  cout << "最小二乘拟合结果 " << endl << " " << endl << endl;
  cout << "垫片直径 = " << endl << " " << m_dRadius_DianPian * 2 * 2.46 << endl << endl;
  cout << "m_Centroid.x = " << endl << " " << m_Centroid_DianPian.x << endl << endl;
  cout << "m_Centroid.y = " << endl << " " << m_Centroid_DianPian.y << endl << endl;

  //**********在原图画拟合圆

  //***垫片图上画圆
  circle(src, Point2d(m_Centroid_DianPian.x, m_Centroid_DianPian.y), m_dRadius_DianPian, Scalar(255, 0, 0), 10, 8, 0);
  //在原图画一个圆圈点
  cv::Point2d point1;//特征点，用以画在图像中  
  point1.x = m_Centroid_DianPian.x;//特征点在图像中横坐标  
  point1.y = m_Centroid_DianPian.y;//特征点在图像中纵坐标  
  cv::circle(src, point1, 4, cv::Scalar(255, 0, 0), 10, 8, 0);
  //*****//
  //图片的缩小与放大
  Mat scale_dst = Mat::zeros(1032, 1544, CV_8UC3); //转化为1544*1032大小的
  resize(src, scale_dst, scale_dst.size());
  namedWindow("src", 0); imshow("src", scale_dst);
  waitKey(0);
  return 0;
}

void LeastSquaresCircleFitting(vector<Point2d> Points, Point2d &Centroid, double& dRadius)
{
  int iNum = (int)Points.size();
  double X1 = 0.0;
  double Y1 = 0.0;
  double X2 = 0.0;
  double Y2 = 0.0;
  double X3 = 0.0;
  double Y3 = 0.0;
  double X1Y1 = 0.0;
  double X1Y2 = 0.0;
  double X2Y1 = 0.0;
  vector<Point2d>::iterator iter;
  vector<Point2d>::iterator end = Points.end();
  for (iter = Points.begin(); iter != end; ++iter)
  {
    X1 = X1 + (*iter).x;
    Y1 = Y1 + (*iter).y;
    X2 = X2 + (*iter).x * (*iter).x;
    Y2 = Y2 + (*iter).y * (*iter).y;
    X3 = X3 + (*iter).x * (*iter).x * (*iter).x;
    Y3 = Y3 + (*iter).y * (*iter).y * (*iter).y;
    X1Y1 = X1Y1 + (*iter).x * (*iter).y;
    X1Y2 = X1Y2 + (*iter).x * (*iter).y * (*iter).y;
    X2Y1 = X2Y1 + (*iter).x * (*iter).x * (*iter).y;
  }
  double C = 0.0;
  double D = 0.0;
  double E = 0.0;
  double G = 0.0;
  double H = 0.0;
  double a = 0.0;
  double b = 0.0;
  double c = 0.0;
  C = iNum * X2 - X1 * X1;
  D = iNum * X1Y1 - X1 * Y1;
  E = iNum * X3 + iNum * X1Y2 - (X2 + Y2) * X1;
  G = iNum * Y2 - Y1 * Y1;
  H = iNum * X2Y1 + iNum * Y3 - (X2 + Y2) * Y1;
  a = (H * D - E * G) / (C * G - D * D);
  b = (H * C - E * D) / (D * D - G * C);
  c = -(a * X1 + b * Y1 + X2 + Y2) / iNum;
  double A = 0.0;
  double B = 0.0;
  double R = 0.0;
  A = a / (-2);
  B = b / (-2);
  R = double(sqrt(a * a + b * b - 4 * c) / 2);
  Centroid.x = A;
  Centroid.y = B;
  dRadius = R;
}