C++实现多项式曲线拟合--polyfit

18 篇文章 8 订阅
2 篇文章 0 订阅

基本原理:幂函数可逼近任意函数。

上式中,N表示多项式阶数,实际应用中一般取3或5;

假设N=5,则:

共有6个未知数,仅需6个点即可求解;

可表示为矩阵方程:

Y的维数为[R*1],U的维数[R * 6],K的维数[6 * 1]。

R> 6时,超定方程求解:


下面是使用C++实现的多项式拟合的程序,程序中使用opencv进行矩阵运算和图像显示。程序分别运行了N=3,5,7,9时的情况,结果如下:



#include <opencv2\opencv.hpp>
#include <iostream>
#include <vector>
using namespace cv;
using namespace std;

Mat polyfit(vector<Point>& in_point, int n);

int main()
{
	//数据输入
	Point in[19] = { Point(50,120),Point(74,110),Point(98,100),Point(122,100),Point(144,80)
		,Point(168,80),Point(192,70),Point(214,50),Point(236,40),Point(262,20)
		,Point(282,20),Point(306,30),Point(328,40),Point(356,50),Point(376,50)
		,Point(400,50),Point(424,50),Point(446,40),Point(468,30) };

	vector<Point> in_point(begin(in),end(in));
	
	//n:多项式阶次
	int n = 9;
	Mat mat_k = polyfit(in_point, n);


	//计算结果可视化
	Mat out(150, 500, CV_8UC3,Scalar::all(0));

	//画出拟合曲线
	for (int i = in[0].x; i < in[size(in)-1].x; ++i)
	{
		Point2d ipt;
		ipt.x = i;
		ipt.y = 0;
		for (int j = 0; j < n + 1; ++j)
		{
			ipt.y += mat_k.at<double>(j, 0)*pow(i,j);
		}
		circle(out, ipt, 1, Scalar(255, 255, 255), CV_FILLED, CV_AA);
	}

	//画出原始散点
	for (int i = 0; i < size(in); ++i)
	{
		Point ipt = in[i];
		circle(out, ipt, 3, Scalar(0, 0, 255), CV_FILLED, CV_AA);
	}

	imshow("9次拟合", out);
	waitKey(0);
	
	return 0;
}

Mat polyfit(vector<Point>& in_point, int n)
{
	int size = in_point.size();
	//所求未知数个数
	int x_num = n + 1;
	//构造矩阵U和Y
	Mat mat_u(size, x_num, CV_64F);
	Mat mat_y(size, 1, CV_64F);

	for (int i = 0; i < mat_u.rows; ++i)
		for (int j = 0; j < mat_u.cols; ++j)
		{
			mat_u.at<double>(i, j) = pow(in_point[i].x, j);
		}

	for (int i = 0; i < mat_y.rows; ++i)
	{
		mat_y.at<double>(i, 0) = in_point[i].y;
	}

	//矩阵运算,获得系数矩阵K
	Mat mat_k(x_num, 1, CV_64F);
	mat_k = (mat_u.t()*mat_u).inv()*mat_u.t()*mat_y;
	cout << mat_k << endl;
	return mat_k;
}

  • 43
    点赞
  • 314
    收藏
    觉得还不错? 一键收藏
  • 19
    评论
#ifndef FUNCTION_H_ #define FUNCTION_H_ #include #include #include "polyfit.h" #include using namespace std; dxs::dxs() { ifstream fin("多项式拟合.txt"); fin>>n; x=new float[n]; y=new float[n]; for(int i=0;i>x[i]; } for(i=0;i>y[i]; } cout<>nn; m=nn+1; u=new float*[m]; for(i=0;i<m;i++) { u[i]=new float[m+1]; }//创建m行,m+1列数组 } void dxs::dfine() { for(int i=0;i<m;i++) { for(int j=0;j<m+1;j++) { u[i][j]=0; } } for(i=0;i<m;i++) { for(int j=0;j<m;j++) { for(int k=0;k<n;k++) { u[i][j]=u[i][j]+pow(x[k],j+i); } } } for(i=0;i<m;i++) { for(int k=0;k<n;k++) { u[i][m]=u[i][m]+pow(x[k],i)*y[k]; } } } void dxs::show() { for(int i=0;i<m;i++) { for(int j=0;j<m+1;j++) { cout<<u[i][j]<<" ";//<<endl; } cout<<endl; } ////显示具有m行m+1列u数组的各元素值 } void dxs::select_main(int k,float **p,int m) { double d; d=*(*(p+k)+k); //cout<<d; int l=k; int i=k+1; for(;i fabs(d)) { d=*(*(p+i)+k); l=i; } else continue; } if(d==0) cout<<"错误"; else { if(k!=l) { for(int j=k;j<m+1;j++) { double t; t=*(*(p+l)+j); *(*(p+l)+j)=*(*(p+k)+j); *(*(p+k)+j)=t; } } } } void dxs::gaosi() { for(int k=0;k<m;k++) { select_main(k,u,m);//调用列主元函数 for(int i=1+k;i<m;i++) { // *(*(p+i)+k)=(float) *(*(p+i)+k) / *(*(p+k)+k); u[i][k]=(float) u[i][k] / u[k][k]; } for(i=k+1;i<m;i++) { for(int j=k+1;j=0;i--) { float a=0; for(int j=i+1;j<m;j++) { //a=a + (*(*(p+i)+j) * *(*(p+j)+m)); a=a+u[i][j] * u[j][m]; } //*(*(p+i)+n-1)= (*(*(p+i)+n-1) - a) / *(*(p+i)+i); u[i][m]= (u[i][m] -a) / u[i][i]; } cout<<"方程组的解为:"<<endl; for(i=0;i<m;i++) { cout<<"a"<<i+1<<"="; cout<<u[i][m]<<endl; // l[i]=*(*(p+i)+n-1); } cout<<"y="<<u[0][m]; for(i=1;i<m;i++) { cout<<showpos<<u[i][m]<<"x"; if(i!=1)cout<<"^"<<noshowpos<<i; } cout<<endl; } dxs::~dxs() { delete[]x,y; delete []*u; } #endif
### 回答1: 我可以为你提供最小二乘法求解二次多项式拟合方程的算法。这是一种常见的方法,可以用来拟合常见的曲线,如抛物线、圆弧等。 首先,需要准备输入数据,包括至少两个点对 (x, y)。这些点对用来拟合二次多项式方程,其中 x 为自变量,y 为因变量。 然后,假设你要拟合的二次多项式方程的形式为: y = ax^2 + bx + c 求解过程如下: 1. 定义三个变量 a、b、c,分别用来存储二次多项式方程的系数。 2. 定义三个变量 sum_xx、sum_xy、sum_yy,分别用来存储以下三个式子的值: - sum_xx = Σ(x^2) - sum_xy = Σ(x*y) - sum_yy = Σ(y^2) 3. 计算 a、b、c 的值。 - a = (sum_xx * sum_y - sum_xy * sum_x) / (n * sum_xx - sum_x^2) - b = (n * sum_xy - sum_x * sum_y) / (n * sum_xx - sum_x^2) - c = (sum_y - b * sum_x - a * sum_xx) / n 其中,n 为点对 (x, y) 的个数,sum_x、sum_y 分别为所有 x 的和以及所有 y 的和。 最后,使用求得的 a、b、c 值,就可以得到拟合的二次多项式方程了。 以下是使用 C 语言实现的示例代 ### 回答2: 要用C语言实现二次多项式拟合方程的算法,首先需要明确的是二次多项式可以表示为y = ax^2 + bx + c的形式。具体实现的算法步骤如下: 1. 定义三个变量a、b、c来表示二次多项式的系数,初始值都设置为0。 2. 使用最小二乘法来进行拟合,即找到最接近给定数据点的二次多项式,并确定最优的系数a、b、c。 3. 读取并保存要进行拟合的数据点,包括x和y的值,可以使用数组或链表来保存。 4. 计算数据点的个数n。 5. 根据最小二乘法的原理,计算公式如下: - 计算各个变量的和: - 计算x和y的和:sum_x、sum_y - 计算x的平方和:sum_xx - 计算x*y的和:sum_xy - 根据公式计算a、b、c的值: - 计算a = (n*sum_xy - sum_x*sum_y) / (n*sum_xx - sum_x^2) - 计算b = (sum_y - a*sum_x) / n - 计算c = (sum_xx*sum_y - sum_xy*sum_x) / (n*sum_xx - sum_x^2) 6. 输出最终的拟合方程为y = ax^2 + bx + c,即输出a、b、c的值。 以上就是使用C语言实现二次多项式拟合方程的算法。在实际实现中,可以根据具体的需求来进行扩展和优化,例如添加输入、输出函数,增加异常处理等。 ### 回答3: 二次多项式拟合方程是指通过一组离散的数据点,在最小二乘准则下,找到一个二次多项式来拟合这些数据。 实现这个算法可以按照以下步骤进行: 1. 定义需要拟合的数据点坐标数组,包括x和y的值。 2. 计算x和y的平均值,分别记为x_mean和y_mean。 3. 计算x和y的偏差数组,即x的偏差数组为(x_i - x_mean),y的偏差数组为(y_i - y_mean),其中x_i和y_i分别为每个数据点的x和y坐标。 4. 计算x的偏差平方和数组,即(x_i - x_mean)²,记为x_deviation_square_sum。计算x的偏差乘以y的偏差数组,即(x_i - x_mean) * (y_i - y_mean),记为xy_deviation_sum。 5. 根据最小二乘准则,计算二次多项式的系数a、b和c。其中a = xy_deviation_sum / x_deviation_square_sum,b = y_mean - a * x_mean,c = y_mean - a * x_mean² - b * x_mean。 6. 得到二次多项式方程为f(x) = ax² + bx + c。 在C语言中,可以使用数组和循环来实现上述算法。具体代码如下: ```c #include <stdio.h> void quadraticFit(double x[], double y[], int n, double *a, double *b, double *c) { double x_mean = 0, y_mean = 0; double x_deviation_square_sum = 0, xy_deviation_sum = 0; // 计算x和y的平均值 for (int i = 0; i < n; i++) { x_mean += x[i]; y_mean += y[i]; } x_mean /= n; y_mean /= n; // 计算偏差平方和和偏差乘积和 for (int i = 0; i < n; i++) { double x_deviation = x[i] - x_mean; double y_deviation = y[i] - y_mean; x_deviation_square_sum += x_deviation * x_deviation; xy_deviation_sum += x_deviation * y_deviation; } // 计算二次多项式的系数 *a = xy_deviation_sum / x_deviation_square_sum; *b = y_mean - (*a) * x_mean; *c = y_mean - (*a) * x_mean * x_mean - (*b) * x_mean; } int main() { double x[] = {1, 2, 3, 4, 5}; double y[] = {2, 4, 6, 8, 10}; int n = sizeof(x) / sizeof(x[0]); double a, b, c; quadraticFit(x, y, n, &a, &b, &c); printf("拟合方程为:f(x) = %fx² + %fx + %f\n", a, b, c); return 0; } ``` 以上代码可以拟合一组简单的数据点,输出二次多项式方程。实际应用中,可以根据实际数据点的情况进行修改和优化。
评论 19
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值