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