faster polyfit

[Dr. G. Buerger]

Hi octave,

2 things. First, I would like to know if it is a general problem of killing long
loops with ^C? On my machine, the loop only stops after the regular (very long)
time.

Second, I've slightly modified the code polyfit.m. Somehow, it consumes an
enormous amount of time calculating the following:

______
...

  [Q, R] = qr (C);
  p = R \ (Q' * y);
______


Here is my alternative polfit.m which seems to work much faster.

regards,

        Gerd

________
## usage:  polfit (x, y, n)
##
## Returns the coefficients of a polynomial p(x) of degree n that
## minimizes sumsq (p(x(i)) - y(i)), i.e., that best fits the data
## in the least squares sense.

function p = polfit (x, y, n)


  if (nargin != 3)
    usage ("polfit (x, y, n)");
  endif
  if (! (is_vector (x) && is_vector (y) && size (x) == size (y)))
    error ("polfit: x and y must be vectors of the same size");
  endif
  if (! (is_scalar (n) && n >= 0 && ! isinf (n) && n == round (n)))
    error ("polfit: n must be a nonnegative integer");
  endif

  l = length (x);
  x = reshape (x, l, 1);
  y = reshape (y, l, 1);

  C = (x * ones (1, n+1)) .^ (ones (l, 1) * (0 : n));
  b = - C' * y;

  ## Compute polynomial coeffients, making returned value compatible
  ## with Matlab.

  A = C'*C;
  L = chol(A);
  y = inv(L') * b;
  p = - inv(L) * y;

  if (! prefer_column_vectors)
    p = p';
  endif

endfunction