интерполяция кубическим сплайном
От: D.S.Khan  
Дата: 04.03.08 08:43
Оценка:
Имеем постоянно обновляющиеся значения, которые необходимо быстро интерполировать.
К примеру по последним 4-м точкам. Т.е. постоянно первые значения удаляются и добавляются новые.
Для обычного сплайна по фиксированным точкам составил следующее, нужно дописать функцию добавления и удаления точек в массив с пересчетом сплайна:
//---------------------------------------------------------------------------

#include <vcl.h>
#pragma hdrstop

#include "Unit1.h"
#include "math.h"
//---------------------------------------------------------------------------
#pragma package(smart_init)
#pragma resource "*.dfm"
TForm1 *Form1;
сlass Spline2D
{
public:
//нулевой коэффициет не используем
static const int num = 100;
float mx[num], my[num]; // 1..10 -контрольные точки сплайна
float a[num],b[num],c[num],d[num]; // 1..10
float alpha[num],beta[num]; // 1..10

static const int n=num-1; // 1..10

float t; //--------то что мы выводим
//---------------------------------------------------------------------------
float delt(int i)
{
return (my[i+1]-my[i]);
}
//---------------------------------------------------------------------------
float delt2(int i)
{
return (my[i+2]-2*my[i+1]+my[i]);
}
//---------------------------------------------------------------------------
void progonka(void)
{
beta[1]=-0.25;
alpha[1]=(3*delt2(1))/4;//(3*delt2(1))/(4*h*h);

for(int k=3;k<=n;k++) // 3..10
{
beta[k-1] =(-1)/(4+beta[k-2]); //beta[10]=0 ???
alpha[k-1]=( ((3*delt2(k-2))/pow((mx[k-1]-mx[k-2]),2)-alpha[k-2]) )/(4+beta[k-2]);
}
}
//---------------------------------------------------------------------------
void calc(void)
{
progonka();
c[n]=0; //10
for(int k=n; k>=2 ;k--) //10..2
{
c[k-1]=beta[k-1]*c[k]+alpha[k-1]; //9..1
}
c[0]=0;
for(int k=2;k<=n;k++) //1..10
{
b[k]=delt(k-1)/(mx[k]-mx[k-1])+((mx[k]-mx[k-1])*(2*c[k]+c[k-1]))/3;
d[k]=(c[k]-c[k-1])/(3*(mx[k]-mx[k-1]));
}
for(int k=0;k<=n; k++)
a[k]=my[k];
}
//---------------------------------------------------------------------------
};

Spline2D spl;
//---------------------------------------------------------------------------
__fastcall TForm1::TForm1(TComponent* Owner)
: TForm(Owner)
{
// инициализация графика
graph.CreateChart(Form1,10,10,2);
graph.ChangeLeftAxis(-200,200);
graph.SetSeriesColor(0,clGreen);
graph.SetSeriesColor(1,clRed);
graph.LinePoints(true);
graph.SlideAll(false);
graph.AllInScreen(true);
// заполнение массивов
spl.mx[0]=0;
spl.my[0]=0;
for(int k=0;k<spl.num;k++) // 0..9
{
spl.mx[k+1]=k+k; //массив аргументов
spl.my[k+1]=random(90); //массив значений
graph.AddXYPoint(0,spl.mx[k+1],spl.my[k+1]); //добавление точек на график
}
}
//---------------------------------------------------------------------------
void __fastcall TForm1::Button1Click(TObject *Sender)
{
graph.ClearAllGraphs(); // очистка графика
// генерация новых точек сплайна
spl.mx[0]=0;
spl.my[0]=0;
for(int k=0;k<spl.num;k++) // 0..9
{
spl.mx[k+1]=k;
spl.my[k+1]=random(160);
graph.AddXYPoint(0,spl.mx[k+1],spl.my[k+1]);// добавление точек на график
}
}
//---------------------------------------------------------------------------
void __fastcall TForm1::Button2Click(TObject *Sender)
{
spl.calc(); // расчет коэффициентов сплайна
float step=0.1f; // с каким шагом интерполируем
for(int k=2;k<spl.num;k++)
{
for(spl.t=spl.mx[k-1];spl.t<(spl.mx[k]);spl.t+=step)
{
// добавление необходимых точек на график
graph.AddXYPoint(1,spl.t,(spl.a[k]
+ spl.b[k]*(spl.t-spl.mx[k])
+ spl.c[k]*pow((spl.t-spl.mx[k]),2)
+ spl.d[k] * pow((spl.t-spl.mx[k]),3)));
}
}
}
//---------------------------------------------------------------------------
 
Подождите ...
Wait...
Пока на собственное сообщение не было ответов, его можно удалить.