unit USpline;

interface

uses
  Windows, Messages, SysUtils, Classes, Graphics, Controls, Forms, Dialogs,
  Menus, ComCtrls, ToolWin, Buttons, StdCtrls;

type
  TKoords2 = array[0..1] of array of extended;
  TKoord = record x,y: extended; end;

  TFormSpline = class(TForm)
    procedure FormCreate(Sender: TObject);
  private
    { Private-Deklarationen }
  public
    { Public-Deklarationen }
  end;

  TSplineXY = class(TObject)
  private
    a,b,c,d,p,St: TKoords2;
    AnzP: integer;
    GP: extended;
  public
    constructor Create;
    procedure Init(TheSt: TKoords2; TheGP: extended);
    procedure InitFelder(i: integer);
    function Calc(TheXt:extended):TKoord;
  end;

var
  FormSpline: TFormSpline;
  SplineXY: TSplineXY;

implementation

{$R *.DFM}


constructor TSplineXY.Create;
begin
  GP:=0;
end;

procedure TSplineXY.Init(TheSt: TKoords2; TheGP: extended);
var i: integer;
begin
  GP:=TheGP;
  AnzP:=Length(TheSt[0]);
  for i:=0 to 1 do begin
    St[i]:=TheSt[i];
    SetLength(a[i],AnzP); SetLength(b[i],AnzP);
    SetLength(c[i],AnzP); SetLength(d[i],AnzP);
    SetLength(p[i],AnzP); p[i,0]:=0.0; p[i,AnzP-1]:=0.0;
    InitFelder(i);
  end;
end;

procedure TSplineXY.InitFelder(i: integer);
var k: integer; z,h,m2,GP2,GPP2,m1,GP1,GPP1: extended;
begin
  GPP1:=0; GP1:=0; m1:=0;
  c[i,0]:=0; d[i,0]:=0;
  for k:=0 to AnzP-2 do begin
    m2:=St[i,k+1]-St[i,k];
    GP2:=GP+2;
    GPP2:=3.0/(GP*(GP+3.0)+3.0);
    a[i,k]:=GPP2/6.0;
    b[i,k]:=a[i,k];
    if k > 0 then begin
      z:=1.0/(GPP1*GP1+GPP2*GP2-GPP1*c[i,k-1]);
      c[i,k]:=z*GPP2;
      h:=6.0*(m2-m1);
      if k = 1 then h:=h-GPP1*p[i,0];
      if k = AnzP-2 then h:=h-GPP2*p[i,AnzP-1];
      d[i,k]:=z*(h-GPP1*d[i,k-1]);
    end;
    GP1:=GP2; GPP1:=GPP2; m1:=m2;
  end;
  p[i,AnzP-2]:=d[i,AnzP-2];
  if AnzP > 3 then for k:=AnzP-3 downto 1 do p[i,k]:=d[i,k]-c[i,k]*p[i,k+1];
  for k:=0 to AnzP-2 do begin
    c[i,k]:=a[i,k]*p[i,k];    d[i,k]:=b[i,k]*p[i,k+1];
    a[i,k]:=St[i,k]-c[i,k];   b[i,k]:=St[i,k+1]-d[i,k];
  end;
end;

function TSplineXY.Calc(TheXt:extended):TKoord;
var k: integer; t,u: extended;
begin
  k:=trunc(TheXt);
  t:=TheXt-k;
  u:=1.0-t;
  Result.x:=a[0,k]*u+b[0,k]*t+c[0,k]*u*u*u/(GP*t+1)+d[0,k]*t*t*t/(GP*u+1);
  Result.y:=a[1,k]*u+b[1,k]*t+c[1,k]*u*u*u/(GP*t+1)+d[1,k]*t*t*t/(GP*u+1);
end;
                 

procedure TFormSpline.FormCreate(Sender: TObject);
begin
  SplineXY:=TSplineXY.Create;
end;

end.