[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: proposal for new m-file function
From: |
Olaf Till |
Subject: |
Re: proposal for new m-file function |
Date: |
Thu, 22 Mar 2012 10:06:52 +0100 |
User-agent: |
Mutt/1.5.20 (2009-06-14) |
On Wed, Mar 21, 2012 at 08:33:01PM -0400, Ben Abbott wrote:
> A user's question about fixed points piecewise-linear fitting led to the
> development of a function looks to be a good fit for Octave's core.
>
>
> https://mailman.cae.wisc.edu/pipermail/help-octave/2012-March/050900.html
>
> The idea was to fit a piece-wise polynomial to a set of data. After some
> discussion between myself and Martin Helm, the attached ppfit.m was produced.
> The name was chosen to match the existing ppval().
>
> The function provides a least-squares fit of a 1D interpolation with
> specified break positions to a set of data.
>
> Demos and tests are included.
>
> Any concern about adding this to Octave's core ?
>
> Ben
The loop could be vectorized:
--- original_ppfit.m 2012-03-22 09:18:22.000000000 +0100
+++ ppfit.m 2012-03-22 09:32:04.000000000 +0100
@@ -113,12 +113,8 @@
"Input data value pairs must have the same size.")
endif
## Use linear least squares for the initial guess
- a = zeros (numel (x), numel (xi));
- for ii = 1:numel (xi)
- yi = zeros (numel (xi), 1);
- yi(ii) = 1;
- a(:, ii) = interp1 (xi, yi, x, method, "extrap");
- endfor
+ yi = eye (numel (xi));
+ a = interp1 (xi, yi, x(:), method, "extrap");
## y = (q*r) * yi + errors
ws = warning ("off", "all");
[q, r, p] = qr (a .* w, 0);
but for a bug in interp1 which transposes the returned matrix if both
"spline" and "extrap" are given.
Olaf
--
public key id EAFE0591, e.g. on x-hkp://pool.sks-keyservers.net
signature.asc
Description: Digital signature