[수학] 다항함수 보간법(Polynomial Interpolation)

다항함수 보간법

•  주어진 점들을 지나는 다항식을 찾는 보간법이다.
•  n개의 데이터들이 존재한다면, 모든 n개의 점들을 통과하는 다항식은 차수가 n-1으로서 유일한 다항식이다.

 

n개의 점을 통과하는 다항식 f(x)

 

 

 

 

보간 다항식 구성

 

n개의 점에 대한 다항식을 풀어쓰면 다음과 같다.

 

n개의 다항식

 

위 다항식을 행렬로 표현하면 다음과 같다.

 

행렬로 표현

 

이때,

 위 식이 성립되므로 각 점(x1~xn)을 통과하는 다항식을 구할 수 있다.

 

 

 

 

OpenCV Code

 n개의 점을 지나는 다항식을 구하는 함수를 OpenCV로 작성해보자.

 

void GetPolynomialPrameters(const vector<cv::Point>& pts, vector<double>& parameters) {
	int nCount = pts.size();
    if (nCount < 2)
        return;

	// fx = ax * params
	cv::Mat fx		= cv::Mat::zeros(cv::Size(1, nCount), CV_64FC1);
	cv::Mat params	= cv::Mat::zeros(cv::Size(1, nCount), CV_64FC1);
	cv::Mat ax		= cv::Mat::zeros(cv::Size(nCount, nCount), CV_64FC1);

	for (int i = 0; i < nCount; ++i) {
		fx.at<double>(0, i) = pts[i].y;

		for (int j = 0; j < nCount; ++j) {
			ax.at<double>(j, i) = std::pow(pts[j].x, i);
		}
	}

	// params = ax^-1 * fx
	cv::Mat inv_ax = cv::Mat::zeros(cv::Size(nCount, nCount), CV_64FC1);
	cv::invert(ax, inv_ax, cv::DECOMP_SVD);

	params = inv_ax * fx;

	parameters.clear();
	for (int i = 0; i < nCount; ++i) {
		parameters.push_back(params.at<double>(0, i));
	}

	return;
}

 

 

먼저 vector로 받은 각 점의 개수를 구하고 fx, params, ax를 각각 사이즈에 맞게 초기화한다.

	int nCount = pts.size();
    if (nCount < 2)
        return;

	// fx = ax * params
	cv::Mat fx		= cv::Mat::zeros(cv::Size(1, nCount), CV_64FC1);
	cv::Mat params	= cv::Mat::zeros(cv::Size(1, nCount), CV_64FC1);
	cv::Mat ax		= cv::Mat::zeros(cv::Size(nCount, nCount), CV_64FC1);

 

fx에는 각 점들의 y값, ax에는 위에서 언급한 다항식 구조에 맞게 x의 값을 넣는다.

	for (int i = 0; i < nCount; ++i) {
		fx.at<double>(0, i) = pts[i].y;

		for (int j = 0; j < nCount; ++j) {
			ax.at<double>(j, i) = std::pow(pts[j].x, i);
		}
	}

 

fx = ax * params 이므로 inv_ax * fx = params이다.

(이때 inv_ax는 ax의 역행렬인데 역행렬이 존재하지 않을 수도 있기 때문에 특이값 분해(SVD, cv::DECOMP_SVD)를 이용하여 근사치를 구하게 된다.)

	// params = ax^-1 * fx
	cv::Mat inv_ax = cv::Mat::zeros(cv::Size(nCount, nCount), CV_64FC1);
	cv::invert(ax, inv_ax, cv::DECOMP_SVD);

	params = inv_ax * fx;

 

 

 

 

 

 

 

 

 

끝.