X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=octave_packages%2Fm%2Fpolynomial%2Fpolyaffine.m;fp=octave_packages%2Fm%2Fpolynomial%2Fpolyaffine.m;h=973268ea9799df94c8d27470b06a89db4d0b58d6;hb=1c0469ada9531828709108a4882a751d2816994a;hp=0000000000000000000000000000000000000000;hpb=63de9f36673d49121015e3695f2c336ea92bc278;p=CreaPhase.git diff --git a/octave_packages/m/polynomial/polyaffine.m b/octave_packages/m/polynomial/polyaffine.m new file mode 100644 index 0000000..973268e --- /dev/null +++ b/octave_packages/m/polynomial/polyaffine.m @@ -0,0 +1,88 @@ +## Copyright (C) 2009-2012 Tony Richardson, Jaroslav Hajek +## +## This file is part of Octave. +## +## Octave is free software; you can redistribute it and/or modify it +## under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 3 of the License, or (at +## your option) any later version. +## +## Octave is distributed in the hope that it will be useful, but +## WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +## General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with Octave; see the file COPYING. If not, see +## . + +## -*- texinfo -*- +## @deftypefn {Function File} {} polyaffine (@var{f}, @var{mu}) +## Return the coefficients of the polynomial vector @var{f} after an affine +## transformation. If @var{f} is the vector representing the polynomial f(x), +## then @code{@var{g} = polyaffine (@var{f}, @var{mu})} is the vector +## representing: +## +## @example +## g(x) = f( (x - @var{mu}(1)) / @var{mu}(2) ) +## @end example +## +## @seealso{polyval, polyfit} +## @end deftypefn + + +function g = polyaffine (f, mu) + + if (nargin != 2) + print_usage (); + endif + + if (! isvector (f)) + error ("polyaffine: F must be a vector"); + endif + + if (! isvector (mu) || length (mu) != 2) + error ("polyaffine: MU must be a two-element vector"); + endif + + lf = length (f); + + ## Ensure that f is a row vector + if (rows (f) > 1) + f = f.'; + endif + + g = f; + + ## Scale. + if (mu(2) != 1) + g = g ./ (mu(2) .^ (lf-1:-1:0)); + endif + + ## Translate. + if (mu(1) != 0) + w = (-mu(1)) .^ (0:lf-1); + ii = lf:-1:1; + g = g(ii) * (toeplitz (w) .* pascal (lf, -1)); + g = g(ii); + endif + +endfunction + + +%!demo +%! f = [1/5 4/5 -7/5 -2]; +%! g = polyaffine (f, [1, 1.2]); +%! x = linspace (-4, 4, 100); +%! plot(x, polyval (f, x), x, polyval (g, x)); +%! legend ("original", "affine"); +%! axis ([-4 4 -3 5]); +%! grid ("on"); + +%!test +%! f = [1/5 4/5 -7/5 -2]; +%! mu = [1, 1.2]; +%! g = polyaffine (f, mu); +%! x = linspace (-4, 4, 100); +%! assert (polyval (f, x, [], mu), polyval (g, x), 1e-10); +