octave_packages/communications-1.1.1/@galois/roots.m

   1 ## Copyright (C) 2002 David Bateman
   2 ##
   3 ## This program is free software; you can redistribute it and/or modify it under
   4 ## the terms of the GNU General Public License as published by the Free Software
   5 ## Foundation; either version 3 of the License, or (at your option) any later
   6 ## version.
   7 ##
   8 ## This program is distributed in the hope that it will be useful, but WITHOUT
   9 ## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  10 ## FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
  11 ## details.
  12 ##
  13 ## You should have received a copy of the GNU General Public License along with
  14 ## this program; if not, see <http://www.gnu.org/licenses/>.
  15
  16 ## -*- texinfo -*-
  17 ## @deftypefn {Function File} {} roots (@var{v})
  18 ##
  19 ## For a vector @var{v} with @math{N} components, return
  20 ## the roots of the polynomial over a Galois Field
  21 ## @iftex
  22 ## @tex
  23 ## $$
  24 ## v_1 z^{N-1} + \cdots + v_{N-1} z + v_N.
  25 ## $$
  26 ## @end tex
  27 ## @end iftex
  28 ## @ifinfo
  29 ##
  30 ## @example
  31 ## v(1) * z^(N-1) + ... + v(N-1) * z + v(N).
  32 ## @end example
  33 ## @end ifinfo
  34 ##
  35 ## The number of roots returned and their value will be determined
  36 ## by the order and primitive polynomial of the Galios Field
  37 ## @end deftypefn
  38
  39 function r = roots (v)
  40
  41   if (nargin != 1)
  42     error("usage: r = roots(v)");
  43   endif
  44
  45   if (!isgalois(v))
  46     error("roots: argument must be a galois variable");
  47   endif
  48
  49   if (min (size (v)) > 1 || nargin != 1)
  50     usage ("roots (v), where v is a galois vector");
  51   endif
  52
  53   v = reshape (v, 1, length(v));
  54   m = v.m;
  55   prim_poly = v.prim_poly;
  56   n = 2^m - 1;
  57   poly = v;
  58   nr = 0;
  59   t = 0;
  60   r = [];
  61
  62   while ((t <= n)  && (length(poly) > 1))
  63     [npoly, nrem] = deconv(poly,gf([1,t],m,prim_poly));
  64     if (any(nrem))
  65       t = t + 1;
  66     else
  67       nr = nr + 1;
  68       r(nr) = t;
  69       poly = npoly;
  70     endif
  71   end
  72
  73   r = gf(r,m,prim_poly);
  74
  75 endfunction