home *** CD-ROM | disk | FTP | other *** search
/ Simtel MSDOS 1992 June / SIMTEL_0692.cdr / msdos / turbopas / pas_sci.arc / CFIT2.PAS < prev    next >
Encoding:
Pascal/Delphi Source File  |  1985-07-28  |  1.6 KB  |  94 lines
  1. program cfit2;        { -> 142+147 }
  2. { plot service included }
  3. { Pascal program to perform a linear least-squares fit }
  4.  
  5. const    max    = 20;
  6.  
  7. type    index    = 1..max;
  8.     ary    = array[index] of real;
  9.  
  10. var    x,y,y_calc    : ary;
  11.     n        : integer;
  12.     first,done    : boolean;
  13.     seed,a,b    : real;
  14.  
  15. function random(dummy: integer): real;
  16. { random number 0-1 }
  17. { define seed=4.0 as global }
  18.  
  19. const    pi    = 3.14159;
  20.  
  21. var    x    : real;
  22.     i    : integer;
  23.  
  24. begin    { RANDOM }
  25.   x:=seed+pi;
  26.   x:=exp(5.0*ln(x));
  27.   seed:=x-trunc(x);
  28.   random:=seed
  29. end;    { RANDOM }
  30.  
  31.  
  32.  
  33. procedure get_data(var x,y: ary;
  34.            var n: integer);
  35. { get values for n and arrays x,y }
  36. { y is randomly scattered about a straight line }
  37.  
  38. const    a = 2.0;
  39.     b = 5.0;
  40.  
  41. var    i,j    : integer;
  42.     fudge    : real;
  43.  
  44. begin
  45.   write('Fudge? (<0 to terminate) ');
  46.   readln(fudge);
  47.   if fudge<0.0 then done:=true
  48.   else
  49.     begin
  50.       repeat
  51.     write('How many points? ');
  52.     readln(n)
  53.       until (n>2) and (n<=max);
  54.       if first then first:=false else ClrScr;
  55.       for i:=1 to n do
  56.     begin
  57.       j:=n+1-i;
  58.       x[i]:=j;
  59.       y[i]:=(a+b*j)*(1.0+(2.0*random(0)-1.0)*fudge)
  60.       end    { for-loop }
  61.     end        { if }
  62. end;        { procedure get_data }
  63.  
  64.  
  65. procedure write_data;
  66. { print out the answers }
  67. var    i    : integer;
  68.  
  69. begin
  70.   writeln;
  71.   writeln('  I      X      Y');
  72.   for i:=1 to n do
  73.     writeln(i:3,x[i]:8:1,y[i]:9:2);
  74.   writeln
  75. end;        { write_data }
  76.  
  77. {$I PLOT.LIB }
  78.  
  79. begin    { MAIN program }
  80.   first:=true;
  81.   ClrScr;
  82.   seed:=4.0;
  83.   done:=false;
  84.   repeat
  85.     get_data(x,y,n);
  86.     if not done then
  87.       begin
  88.     write_data;
  89.     plot(x,y,y,-n);
  90.     { ***** --->  more lines to be added here ********* }
  91.     end
  92.   until done
  93. end.
  94.