三次樣條插值算法C++實現

三次樣條插值算法

1 總體說明

三次樣條插值算法是一種計算量和效果都比較理想的插值算法。關於三次樣條插值算法的原理這里不做過多的解釋，下面的代碼是我在網上收集了兩種C++實現版本的基礎上自己整合的一個版本。由於本人剛接觸C++不久，水平有限。沒有使用模板機制將代碼做的更通用。關於算法實現有下面幾點說明。

1. 所有有關的類都被包含到SplineSpace命名空間中。
2. SplineSpace中一個有三個類分別是異常類(SplineFailure),接口類(SplineInterface)和實現類(Spline)。有一個枚舉類型說明邊界條件(BoundaryCondition),取值為：GivenFirstOrder和GivenSecondOrder。分別對應I型邊界條件和II型邊界條件。
3. 接口類定義了Spline在實現的過程中必須要有的三個方法：單點插值、多點插值和自動生成插值序列。
4. 異常類是可能被實現類拋出的類，如果在實現類的運行過程中出現了已知數據過少構造失敗、使用了外插值、設定輸出點數過少等行為會拋出該類。因此應該將插值的過程用try...catch(SplineFailure sf)包裹起來。如：

double x0[2]={1,2};
double y0[2]={3,4};

try
{
  SplineInterface* sp = new Spline(x0,y0,2);
  //...
}
catch(SplineFailure sf)
{
  cout<<sf.GetMessage()<<endl;
}

上面代碼就會拋出異常並顯示“構造失敗,已知點數過少”。

2 插值方法調用

2.1單點插值

調用方法如下：

#include <iostream>
#include "Spline.h"

using namespace std;
using namespace SplineSpace;

int main(void)
{
  //單點插值測試

	double x0[5]={1,2,4,5,6};		//已知的數據點
	double y0[5]={1,3,4,2,5};
	try
	{
		//Spline sp(x0,y0,5,GivenSecondOrder,0,0);
		SplineInterface* sp = new Spline(x0,y0,5);	//使用接口，且使用默認邊界條件
		double x=4.5;
		double y;
		sp->SinglePointInterp(x,y);	//求x的插值結果y
		cout<<"x="<<x<<"時的插值結果為:"<<y<<endl;
	}
	catch(SplineFailure sf)
	{
		cout<<sf.GetMessage()<<endl;
	}
	getchar();	//程序暫停
}

此時屏幕會輸出"x=4.5時的插值結果為2.71107"。

1. 可以直接構造Spline對象進行使用，也可以將它轉化成對應的接口使用。
2. 默認邊界條件為II型邊界條件(已知邊界的二階導數)，默認左右邊界二階導數為都0，即自然邊界條件。
3. 已知的數據點數要至少為3（3個點對應一條曲線）。

2.2多點插值

調用方法如下（省略了和上面重復的代碼部分）：

//多點插值測試

double x0[5]={1,2,4,5,6};		//已知的數據點
double y0[5]={1,3,4,2,5};

double x[4] = {1.5,2.5,3.5,4.5};	//插值點
double y[4];
double leftBound=0,RightBound=0;	//邊界導數

try
{
  Spline sp(x0,y0,5,GivenSecondOrder,leftBound,RightBound);
  sp.MultiPointInterp(x,4,y);			//求x的插值結果y
  for(int i = 0;i < 4;i++)
  {
    cout<<"x="<<x[i]<<"時的插值結果為:"<<y[i]<<endl;
  }
}
catch(SplineFailure sf)
{
  cout<<sf.GetMessage()<<endl;
}
getchar();	//程序暫停

顯示結果如下：

x=1.5時的插值結果為:2.01383
x=2.5時的插值結果為:3.8978
x=3.5時的插值結果為:4.62372
x=4.5時的插值結果為:2.71107

2.3自動生成插值序列

自動生成插值序列方法會根據指定的插值點數，在插值自變量區間等間距的生成插值點。這個方法在繪圖的時候非常好用。指定的插值點數越多，繪制出來的圖形質量越好。調用方法如下：

//自動插值測試

double x0[5]={1,2,4,5,6};		//已知的數據點
double y0[5]={1,3,4,2,5};

double x[10];	//插值點
double y[10];

try
{
  SplineInterface* sp = new Spline(x0,y0,5);	//使用接口，且使用默認邊界條件
  sp->AutoInterp(10,x,y);			//求x的插值結果y

  for(int i = 0;i < 10;i++)
    cout<<x[i]<<",";
  cout<<endl;
  for(int i = 0;i < 10;i++)
    cout<<y[i]<<",";
}
catch(SplineFailure sf)
{
  cout<<sf.GetMessage()<<endl;
}
getchar();	//程序暫停

顯示結果如下：

1,1.55556,2.11111,2.66667,3.22222,3.77778,4.33333,4.88889,5.44444,6,
1,2.12528,3.21225,4.14572,4.639,4.38218,3.1524,2.02937,2.79183,5,

將上面結果復制粘貼到matlab中繪圖，並且和matlab中的三次樣條擬合結果作比較。matlab代碼如下:

clear
clc

x0=[1,2,4,5,6];         %已知的自變量
y0=[1,3,4,2,5];         %已知的因變量

pp=csape(x0,y0,'second',[0,0]); %二階邊界條件,且導數都為0
x = linspace(1,6,200);
y = ppval(pp,x);

%C++算法得出了的結果
xCpp = [1,1.55556,2.11111,2.66667,3.22222,3.77778,4.33333,4.88889,5.44444,6];
yCpp = [1,2.12528,3.21225,4.14572,4.639,4.38218,3.1524,2.02937,2.79183,5];

h=figure;
plot(xCpp,yCpp,'*',x,y);
legend('c++插值結果','matlab插值結果')
title('插值結果')
h.Name = '插值結果';
h.NumberTitle = 'off';

結果如下：

可以看到C++中的算法得到的結果和matlab中得到的結果相同。
C++源程序代碼：
有CSDN下載積分的童鞋可以使用右邊的下載鏈接支持一下：三次樣條插值C++代碼
沒有積分的童鞋也可以直接在本站下載源碼：三次樣條插值C++代碼