home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 9 / 09.iso / l / l048 / 1.ddi / SPLINE.HGH < prev    next >
Encoding:
Text File  |  1986-03-20  |  3.6 KB  |  133 lines
  1. (***********************************************************)
  2. (*                                                         *)
  3. (*                TURBO GRAPHIX version 1.06A              *)
  4. (*                                                         *)
  5. (*                   Cubic spline module                   *)
  6. (*                  Module version  1.06A                  *)
  7. (*                                                         *)
  8. (*                  Copyright (C) 1985 by                  *)
  9. (*                  BORLAND International                  *)
  10. (*                                                         *)
  11. (***********************************************************)
  12.  
  13. procedure Spline(A : PlotArray; N : integer; X1, XM : real;
  14.                  var B : PlotArray; M : integer);
  15. const
  16.   MaxSpline = 50;
  17. type
  18.   Vector = array[1..MaxSpline] of real;
  19. var
  20.   X, Y, Z : Vector;
  21.   I : integer;
  22.   DeltaX : real;
  23.  
  24. procedure Stg(Vector1, Vector2, Vector3 : Vector;
  25.               var Vector4 : Vector; NPts : integer);
  26. var
  27.   I : integer;
  28.   Factor : real;
  29.  
  30. begin
  31.   for I := 2 to NPts do
  32.   begin
  33.     Factor := Vector1[I - 1] / Vector2[I - 1];
  34.     Vector2[I] := Vector2[I] - Factor * Vector3[I - 1];
  35.     Vector4[I] := Vector4[I] - Factor * Vector4[I - 1];
  36.   end;
  37.   Vector4[NPts] := Vector4[NPts] / Vector2[NPts];
  38.   for I := 1 to NPts - 1 do
  39.   Vector4[NPts - I] := (Vector4[NPts - I] - Vector3[NPts - I] *
  40.                         Vector4[NPts - I + 1]) / Vector2[NPts - I];
  41. end; { Stg }
  42.  
  43. procedure Sc(X, Y : Vector; var Z : Vector; NPts : integer);
  44. var
  45.   I : integer;
  46.   D, C : Vector;
  47.  
  48. begin
  49.   D[1] := 1.0;
  50.   C[1] := 0.5;
  51.   Z[1] := 0.5;
  52.   for I := 2 to NPts - 1 do
  53.   begin
  54.     D[I] := 2.0 * (X[I + 1] - X[I - 1]);
  55.     C[I] := X[I + 1] - X[I];
  56.     Z[I] := 6.0 * ((Y[I + 1] - Y[I]) / (X[I + 1] - X[I]) -
  57.                    (Y[I] - Y[I - 1]) / (X[I] - X[I - 1]));
  58.   end;
  59.   D[NPts] := 1.0;
  60.   C[NPts - 1] := 0.0;
  61.   C[NPts] := 0.0;
  62.   Z[NPts] := 0.0;
  63.   Stg(C, D, C, Z, NPts);
  64. end; { Sc }
  65.  
  66. function Si(V : real; X, Y, Z : Vector; NPts : integer) : real;
  67. var
  68.   I, J : integer;
  69.   Dummy, Ai, Hi : real;
  70.  
  71. begin
  72.   if (V > X[1]) and (V < X[NPts]) then
  73.     begin
  74.       J := 1;
  75.       repeat
  76.         J := J + 1;
  77.         I := NPts - J;
  78.         Dummy := V - X[I];
  79.       until (Dummy >= 0.0) or (I = 2);
  80.       Hi := X[I + 1] - X[I];
  81.       Ai := Dummy * (Z[I + 1] - Z[I]) / (6.0 * Hi) + 0.5 * Z[I];
  82.       Ai := Dummy * Ai + (Y[I + 1] - Y[I]) / Hi - Hi *
  83.                          (2.0 * Z[I] + Z[I + 1]) / 6.0;
  84.       Si := Dummy * Ai + Y[I];
  85.     end
  86.   else
  87.     if V = X[1] then
  88.       Si := Y[1]
  89.     else
  90.       Si := Y[NPts];
  91. end; { Si }
  92.  
  93. procedure Sia(X, Y : Vector; NPts : integer; XInt : Vector;
  94.               var YInt : Vector; N : integer);
  95. var
  96.   I : integer;
  97.   V3 : Vector;
  98.  
  99. begin
  100.   Sc(X, Y, V3, NPts);
  101.   for I := 1 to N do
  102.     YInt[I] := Si(XInt[I], X, Y, V3, NPts);
  103. end; { Sia }
  104.  
  105. begin { Spline }
  106.   if (abs(N) >= 2) and (abs(M) >= 2) then
  107.     begin
  108.       if ((X1 >= A[1, 1]) and (XM <= A[N, 1])) and (M >= 2) then
  109.         begin
  110.           DeltaX := (XM - X1) / (M - 1);
  111.           for I := 1 to N do
  112.           begin
  113.             X[I] := A[I, 1];
  114.             Y[I] := A[I, 2];
  115.           end;
  116.           for I := 2 to M - 1 do
  117.             Z[I] := X1 + (I - 1) * DeltaX;
  118.           Z[1] := X1;
  119.           Z[M] := XM;
  120.           Sia(X, Y, N, Z, Y, M);
  121.           for I := 1 to M do
  122.           begin
  123.             B[I, 1] := Z[I];
  124.             B[I, 2] := Y[I];
  125.           end;
  126.         end
  127.       else
  128.         Error(20, 7);
  129.     end
  130.   else
  131.     Error(20, 4);
  132. end; { Spline }
  133.