由点坐标拟合圆和计算圆度

由点坐标拟合圆和计算圆度

1.根据坐标点拟合圆
最小二乘法拟合圆

   void FitCircle(const vector<Point> pts,double& circleX,double& circleY,double& radius)
   {
	   circleX=0;
	   circleY=0;
	   radius=0;
	   if(pts.size()<3)
	   {
		   return;
	   }

	   int i=0;
	   double X1=0,Y1=0,X2=0,Y2=0,X3=0,Y3=0;

	   double X1Y1=0,X1Y2=0,X2Y1=0;

	   for (i=0;i<pts.size();i++)
	   {
		   X1+=pts.at(i).x;
		   Y1+=pts.at(i).y;

		   X2+=std::pow(pts.at(i).x*1.0,2);
		   Y2+=std::pow(pts.at(i).y*1.0,2);

		   X3+=std::pow(pts.at(i).x*1.0,3);
		   Y3+=std::pow(pts.at(i).y*1.0,3);

		   X1Y1+=pts.at(i).x*pts.at(i).y;
		   X1Y2+=pts.at(i).x*std::pow(pts.at(i).y*1.0,2);
		   X2Y1+=std::pow(pts.at(i).x*1.0,2)*pts.at(i).y;
	   }

	   double C,D,E,G,H,N;
	   double a,b,c;
	   N = pts.size();
	   C = N*X2 - X1*X1;
	   D= N*X1Y1 - X1*Y1;
	   E = N*X3+N*X1Y2-(X2+Y2)*X1;
	   G = N*Y2 - Y1*Y1;
	   H = N*X2Y1 + N*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)/N;

		circleX = a/(-2);
		circleY = b/(-2);
		radius = std::sqrt(a*a+b*b-4*c)/2;


   }

2.圆度计算方式1

计算区域的面积和近似最小外接圆的面积的比值

  double CalCircularity(const vector<Point> contours,double area)
  {

	  double dCircular = 0;
	  int nCount = contours.size();
	  double dSumX = 0, dSumY = 0;

	  for (size_t i = 0; i < nCount; i++)
	  {
		  dSumX += contours.at(i).x;
		  dSumY += contours.at(i).y;
	  }
	  double dDistSum = 0;
	  double dtemp = 0;
	  double meanX = dSumX * 1.0 / nCount;
	  double meanY = dSumY * 1.0 / nCount;
	  double dMaxDist=0;
	  for (size_t i = 0; i < nCount; i++)
	  {
		  dtemp = pow(contours.at(i).x - meanX, 2) + pow(contours.at(i).y - meanY, 2);
		  dtemp=sqrt(dtemp);
		  if (dMaxDist<dtemp)
		  {
			  dMaxDist = dtemp;
		  }
		
	  }

	  double dRatio = area/(dMaxDist*dMaxDist*CV_PI);

	  dCircular = STDMIN(1,dRatio);

	  return dCircular;
  }

3.圆度计算方式2
计算区域轮廓上各点到区域中心的平均距离（Distance）与这些距离的标准差（Sigma）之间的关系

 double CalRoundess(const vector<Point> contours)
 {
	 double dCircular = 0;
	 int nCount = contours.size();
	 double dSumX = 0, dSumY = 0;
	 for (size_t i = 0; i < nCount; i++)
	 {
		 dSumX += contours.at(i).x;
		 dSumY += contours.at(i).y;
	 }
	 double dDistSum = 0;
	 double dtemp = 0;
	 double meanX = dSumX * 1.0 / nCount;
	 double meanY = dSumY * 1.0 / nCount;

	 for (size_t i = 0; i < nCount; i++)
	 {
		 dtemp = pow(contours.at(i).x - meanX, 2) + pow(contours.at(i).y - meanY, 2);

		 dDistSum += sqrt(dtemp);
	 }

	 double dMeanDist = dDistSum * 1.0 / nCount;

	 double dSigmaSum = 0;
	 for (size_t i = 0; i < nCount; i++)
	 {
		 dtemp = pow(contours.at(i).x - meanX, 2) + pow(contours.at(i).y - meanY, 2);

		 dSigmaSum +=pow( sqrt(dtemp)- dMeanDist,2);
	 }

	 double dMeanSigma = sqrt(dSigmaSum * 1.0 / nCount);

	 dCircular = 1 - dMeanSigma / dMeanDist;

	 return dCircular;
 }