home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 9 / 09.iso / l / l042 / 1.ddi / CHAP3.ARC / CUBE_FRE.PAS < prev    next >
Encoding:
Pascal/Delphi Source File  |  1987-12-30  |  10.1 KB  |  288 lines
  1. program CubicSplineFree_Prog;
  2.  
  3. {------------------------------------------------------------------------}
  4. {-                                                                      -}
  5. {-     Turbo Pascal Numerical Methods Toolbox                           -}
  6. {-     Copyright (c) 1986, 87 by Borland International, Inc.            -}
  7. {-                                                                      -}
  8. {-       Purpose: This program demonstrates interpolation with          -}
  9. {-                a free cubic spline.                                  -}
  10. {-                                                                      -}
  11. {-       Unit   : Interp    procedure CubicSplineFree                   -}
  12. {-                                                                      -}
  13. {------------------------------------------------------------------------}
  14.  
  15. {$I-}                    { Disable I/O error trapping }
  16. {$R+}                    { Enable range checking }
  17.  
  18. uses
  19.   Interp, Dos, Crt, Common;
  20.  
  21. var
  22.   XData, YData : TNvector;                   { Data points (X,Y) }
  23.   NumPoints : integer;                       { Number of data points }
  24.   Coef0, Coef1, Coef2, Coef3 : TNvector;     { Coefficients of the spline }
  25.   NumInter : integer;                        { Number interpolating points }
  26.   XInter, YInter : TNvector;                 { Interpolating points }
  27.   Error : byte;                              { Flags an error }
  28.  
  29. procedure Initialize(var Coef0     : TNvector;
  30.                      var Coef1     : TNvector;
  31.                      var Coef2     : TNvector;
  32.                      var Coef3     : TNvector;
  33.                      var XData     : TNvector;
  34.                      var YData     : TNvector;
  35.                      var NumPoints : integer;
  36.                      var NumInter  : integer;
  37.                      var XInter    : TNvector;
  38.                      var YInter    : TNvector;
  39.                      var Error     : byte);
  40.  
  41. {----------------------------------------------------------------}
  42. {- Output: Coef0, Coef1, Coef2, Coef3, XData, YData, NumPoints, -}
  43. {-         NumInter, XInter, YInter, Error                      -}
  44. {-                                                              -}
  45. {- This procedure initializes the above variables to zero.      -}
  46. {----------------------------------------------------------------}
  47.  
  48. begin
  49.   FillChar(Coef0, SizeOf(Coef0), 0);
  50.   FillChar(Coef1, SizeOf(Coef1), 0);
  51.   FillChar(Coef2, SizeOf(Coef2), 0);
  52.   FillChar(Coef3, SizeOf(Coef3), 0);
  53.   FillChar(XData, SizeOf(XData), 0);
  54.   FillChar(YData, SizeOf(YData), 0);
  55.   FillChar(XInter, SizeOf(XInter), 0);
  56.   FillChar(YInter, SizeOf(YInter), 0);
  57.   NumPoints := 0;
  58.   NumInter := 0;
  59.   Error := 0;
  60. end; { procedure Initialize }
  61.  
  62. procedure GetData(var NumPoints : integer;
  63.                   var XData     : TNvector;
  64.                   var YData     : TNvector;
  65.                   var NumInter  : integer;
  66.                   var XInter    : TNvector);
  67.  
  68. {--------------------------------------------------------------}
  69. {- Output: NumPoints, NumInter, XData, YData, XInter          -}
  70. {-                                                            -}
  71. {- This procedure reads in data from either the keyboard      -}
  72. {- or a data file.  The number of data points (NumPoints),    -}
  73. {- the data points (XData, YData), the number of interpolated -}
  74. {- points (NumInter) and the X values at which to interpolate -}
  75. {- (XInter) are all read in here.                             -}
  76. {--------------------------------------------------------------}
  77.  
  78. var
  79.   Ch : char;
  80.  
  81. procedure GetTwoVectorsFromFile(var NumPoints : integer;
  82.                                 var XData     : TNvector;
  83.                                 var YData     : TNvector);
  84.  
  85. {-------------------------------------------------------------}
  86. {- Output: NumPoints, XData, YData                           -}
  87. {-                                                           -}
  88. {- This procedure reads in the data points from a data file. -}
  89. {-------------------------------------------------------------}
  90.  
  91. var
  92.   Filename : string[255];
  93.   InFile : text;
  94.  
  95. begin
  96.   Writeln;
  97.   repeat
  98.     Write('File name? ');
  99.     Readln(Filename);
  100.     Assign(InFile,Filename);
  101.     Reset(InFile);
  102.     IOCheck;
  103.   until not IOerr;
  104.   NumPoints := 0;
  105.   while not(EOF(InFile)) do
  106.   begin
  107.     NumPoints := Succ(NumPoints);
  108.     Readln(InFile, XData[NumPoints], YData[NumPoints]);
  109.     IOCheck;
  110.   end;
  111.   Close(InFile);
  112. end; { procedure GetTwoVectorsFromFile }
  113.  
  114. procedure GetTwoVectorsFromKeyboard(var NumPoints : integer;
  115.                                     var XData     : TNvector;
  116.                                     var YData     : TNvector);
  117.  
  118. {--------------------------------------------------------------}
  119. {- Output: NumPoints, XData, YData                            -}
  120. {-                                                            -}
  121. {- This procedure reads in the data points from the keyboard. -}
  122. {--------------------------------------------------------------}
  123.  
  124. var
  125.   Term : integer;
  126.  
  127. begin
  128.   NumPoints := 0;
  129.   Writeln;
  130.   repeat
  131.     Write('Number of points (0-', TNArraySize, ')? ');
  132.     Readln(NumPoints);
  133.     IOCheck;
  134.   until ((NumPoints >= 0) and (NumPoints <= TNArraySize) and not IOerr);
  135.   Writeln;
  136.   Write('Enter the X ');
  137.   Writeln('and Y values, separated by a space (not a comma):');
  138.   for Term := 1 to NumPoints do
  139.   repeat
  140.     Write('X[', Term, '], Y[', Term, ']:');
  141.     Read(XData[Term], YData[Term]);
  142.     Writeln;
  143.     IOCheck;
  144.   until not IOerr;
  145. end; { procedure GetTwoVectorsFromKeyboard }
  146.  
  147. procedure GetOneVectorFromFile(var NumInter : integer;
  148.                                var XInter   : TNvector);
  149.  
  150.  
  151. {------------------------------------------}
  152. {- Output: NumInter, XInter               -}
  153. {-                                        -}
  154. {- This procedure reads in the points at  -}
  155. {- which to interpolate from a data file. -}
  156. {------------------------------------------}
  157.  
  158. var
  159.   Filename : string[255];
  160.   InFile : text;
  161. begin
  162.   Writeln;
  163.   repeat
  164.     Write('File name? ');
  165.     Readln(Filename);
  166.     Assign(InFile, Filename);
  167.     Reset(InFile);
  168.     IOCheck;
  169.   until not IOerr;
  170.   NumInter := 0;
  171.   while not(EOF(InFile)) do
  172.   begin
  173.     NumInter := Succ(NumInter);
  174.     Readln(InFile, XInter[NumInter]);
  175.     IOCheck;
  176.   end;
  177.   Close(InFile);
  178. end; { procedure GetOneVectorFromFile }
  179.  
  180. procedure GetOneVectorFromKeyboard(var NumInter : integer;
  181.                                    var XInter   : TNvector);
  182.  
  183. {-------------------------------------------}
  184. {- Output: NumInter, XInter                -}
  185. {-                                         -}
  186. {- This procedure reads in the points at   -}
  187. {- which to interpolate from the keyboard. -}
  188. {-------------------------------------------}
  189.  
  190. var
  191.   Term : integer;
  192.  
  193. begin
  194.   NumInter := 0;
  195.   Writeln;
  196.   repeat
  197.     Write('Number of points (0-',TNArraySize,')?');
  198.     Readln(NumInter);
  199.     IOCheck;
  200.   until ((NumInter >= 0) and (NumInter <= TNArraySize) and not IOerr);
  201.   Writeln;
  202.   for Term := 1 to NumInter do
  203.   repeat
  204.     Write('Point ', Term, ':');
  205.     Readln(XInter[Term]);
  206.     IOCheck;
  207.   until not IOerr;
  208. end; { procedure GetOneVectorFromKeyboard }
  209.  
  210. begin { procedure GetData }
  211.   case InputChannel('Input Data Points From') of
  212.     'K' : GetTwoVectorsFromKeyboard(NumPoints, XData, YData);
  213.     'F' : GetTwoVectorsFromFile(NumPoints, XData, YData);
  214.   end;
  215.   Writeln;
  216.   case InputChannel('Input Interpolated Points From') of
  217.     'K' : GetOneVectorFromKeyboard(NumInter, XInter);
  218.     'F' : GetOneVectorFromFile(NumInter, XInter);
  219.   end;
  220.   GetOutputFile(OutFile);
  221. end; { procedure GetData }
  222.  
  223. procedure Results(NumPoints : integer;
  224.               var XData     : TNvector;
  225.               var YData     : TNvector;
  226.               var Coef0     : TNvector;
  227.               var Coef1     : TNvector;
  228.               var Coef2     : TNvector;
  229.               var Coef3     : TNvector;
  230.                   NumInter  : integer;
  231.               var XInter    : TNvector;
  232.               var YInter    : TNvector;
  233.                   Error     : byte);
  234.  
  235. {------------------------------------------------------------}
  236. {- this procedure outputs the results to the device OutFile -}
  237. {------------------------------------------------------------}
  238.  
  239. var
  240.   Index : integer;
  241.  
  242. begin
  243.   Writeln(OutFile);
  244.   Writeln(OutFile);
  245.   Writeln(OutFile, 'Data :           X                   Y');
  246.   for Index := 1 to NumPoints do
  247.     Writeln(OutFile, Index:3, ':    ', XData[Index] : 15 : 10,
  248.                      '      ', YData[Index] : 15 : 10);
  249.   Writeln(OutFile);
  250.   if Error >= 1 then
  251.     DisplayError;
  252.   case Error of
  253.     0 : begin
  254.           Writeln(OutFile, 'Splines:', ' ':6, 'Coef0', ' ':12, 'Coef1',
  255.                            ' ':13, 'Coef2', ' ':14, 'Coef3');
  256.           for Index := 1 to NumPoints-1 do
  257.             Writeln(OutFile, '  ', Index : 3, ':', '  ', Coef0[Index]:15:10,
  258.                              ' ':3, Coef1[Index]:15:10, ' ':3,
  259.                              Coef2[Index]:15:10, ' ':3, Coef3[Index]:15:10);
  260.           Writeln(OutFile);
  261.           Writeln(OutFile, 'Interpolated Points:    X                   Y');
  262.           for Index := 1 to NumInter do
  263.             Writeln(OutFile, Index:10, ':    ',  XInter[Index] : 15 : 10,
  264.                              '      ', YInter[Index] : 15 : 10);
  265.         end;
  266.  
  267.     1 : Writeln(OutFile, 'The X points must be unique.');
  268.  
  269.     2 : Writeln(OutFile,
  270.                        'The X points must be in increasing sequential order.');
  271.  
  272.     3 : Writeln(OutFile, 'There must be at least two data points.');
  273.  
  274.   end; { case }
  275. end; { procedure Results }
  276.  
  277. begin { program CubicSplineFree }
  278.   ClrScr;
  279.   Initialize(Coef0, Coef1, Coef2, Coef3, XData, YData, NumPoints,
  280.              NumInter, XInter, YInter, Error);
  281.   GetData(NumPoints, XData, YData, NumInter, XInter);
  282.   CubicSplineFree(NumPoints, XData, YData, NumInter, XInter,
  283.                   Coef0, Coef1, Coef2, Coef3, YInter, Error);
  284.   Results(NumPoints, XData, YData, Coef0, Coef1, Coef2, Coef3,
  285.           NumInter, XInter, YInter, Error);
  286.   Close(OutFile);
  287. end. { program CubicSplineFree }
  288.