C++實現多項式曲線擬合--polyfit
阿新 • • 發佈:2018-12-30
基本原理:冪函式可逼近任意函式。
上式中,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; }