help-octave
[Top][All Lists]

## faster polyfit

 From: Dr. G. Buerger Subject: faster polyfit Date: Mon, 9 Jun 1997 16:17:16 +0100

```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

```