二维点拟合椭圆

原理

椭圆的一般方程：
$$A{x}^2+Bxy+C{y}^2+Dx+Ey+1=0$$
长轴倾角：
$$\theta =\frac{1}{2}arctan\frac{B}{A-C}$$
椭圆几何中心：
$$\begin{gathered} X_c=\frac{BE-2CD}{4AC-B^2} \\ {Y}_c=\frac{BD-2AE}{4AC-B^2} \end{gathered}$$
长轴短轴分别为：
$$\begin{gathered} a^2=\frac{2(A{X}_c^2+C{Y}_c^2+B{X}_c{Y}_c-1)}{A+C+\sqrt{(A-C{)}^2+B^2}} \\ {b}^2=\frac{2(A{X}_c^2+C{Y}_c^2+B{X}_c{Y}_c-1)}{A+C-\sqrt{(A-C{)}^2+B^2}} \end{gathered}$$
离心率：
$$e=\frac{c}{a}=\frac{\sqrt{a^2-b^2}}{a}$$

假设有n个椭圆上的点需要拟合：
$$(x_i,y_i),i=1,2,...,n$$

则：
$$\left[\begin{array}{ccccc}{x}_i^2& x_i{y}_i& y_i^2& x& y\end{array}\right]\left[\begin{array}{c}A\\ B\\ C\\ D\\ E\end{array}\right]=-1$$

代码

//椭圆的中心坐标
cv::Point2f measure::EllipticCentrePoint(std::vector<cv::Point2f> &points, cv::Point2f &CentrePoints)
{
	//float A, B, C, D, E;
	cv::Mat InputArray1 = cv::Mat(points.size(), 5, CV_32F, cv::Scalar(0));
	cv::Mat OutputArray;

	for (int i = 0; i < points.size(); i++)
	{
		//Ax^2+Bxy+Cy^2+Dx+Ey+1=0椭圆方程
		InputArray1.at<float>(i, 0) = float(pow(points[i].x, 2));//x^2
		InputArray1.at<float>(i, 1) = points[i].x *points[i].y;//xy
		InputArray1.at<float>(i, 2) = pow(points[i].y, 2);//y^2
		InputArray1.at<float>(i, 3) = points[i].x;//x
		InputArray1.at<float>(i, 4) = points[i].y;//y
		//InputArray1.at<float>(i, 4) = 1;//1
	}
	cv::Mat InputArray2 = cv::Mat(points.size(), 1, CV_32F, cv::Scalar(-1));
	cv::solve(InputArray1, InputArray2, OutputArray, CV_SVD);
	float A = OutputArray.at<float>(0, 0);
	float B = OutputArray.at<float>(1, 0);
	float C = OutputArray.at<float>(2, 0);
	float D = OutputArray.at<float>(3, 0);
	float E = OutputArray.at<float>(4, 0);
	float EllipseCenterX = float(B * E - 2 * C * D) / (4 * A * C - B * B);
	float EllipseCenterY = float(B * D - 2 * A * E) / (4 * A * C - B * B);
	CentrePoints.x = EllipseCenterX;
	CentrePoints.y = EllipseCenterY;
	return CentrePoints;
}