octave_packages/signal-1.1.3/sos2tf.m

   1 %% Copyright (C) 2005 Julius O. Smith III <jos@ccrma.stanford.edu>
   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} {[@var{B}, @var{A}] =} sos2tf (@var{sos}, @var{Bscale})
  18 %% Convert series second-order sections to direct form @math{H(z) = B(z)/A(z)}.
  19 %%
  20 %% INPUTS:
  21 %% @itemize
  22 %%
  23 %% @item
  24 %% @var{sos} = matrix of series second-order sections, one per row:@*
  25 %% @var{sos} = [@var{B1}.' @var{A1}.'; ...; @var{BN}.' @var{AN}.'], where@*
  26 %% @code{@var{B1}.'==[b0 b1 b2] and @var{A1}.'==[1 a1 a2]} for
  27 %% section 1, etc.@*
  28 %% b0 must be nonzero for each section.@*
  29 %% See @code{filter()} for documentation of the
  30 %% second-order direct-form filter coefficients @var{B}i and @var{A}i.
  31 %%
  32 %% @item
  33 %% @var{Bscale} is an overall gain factor that effectively scales
  34 %% the output @var{B} vector (or any one of the input @var{B}i vectors).
  35 %% @end itemize
  36 %%
  37 %% RETURNED:
  38 %% @var{B} and @var{A} are vectors specifying the digital filter @math{H(z) = B(z)/A(z)}.
  39 %% See @code{filter()} for further details.
  40 %%
  41 %% @seealso{tf2sos zp2sos sos2pz zp2tf tf2zp}
  42 %% @end deftypefn
  43
  44 function [B,A] = sos2tf(sos, Bscale = 1)
  45
  46   if (nargin < 1 || nargin > 2)
  47     print_usage;
  48   endif
  49
  50   [N,M] = size(sos);
  51
  52   if M~=6
  53     error('sos2tf: sos matrix should be N by 6');
  54   end
  55
  56   A = 1;
  57   B = 1;
  58
  59   for i=1:N
  60     B = conv(B, sos(i,1:3));
  61     A = conv(A, sos(i,4:6));
  62   end
  63
  64   nB = length(B);
  65   while nB && B(nB)==0
  66     B=B(1:nB-1);
  67     nB=length(B);
  68   end
  69
  70   nA = length(A);
  71   while nA && A(nA)==0
  72     A=A(1:nA-1);
  73     nA=length(A);
  74   end
  75   B = B * Bscale;
  76
  77 endfunction
  78
  79 %!test
  80 %! B=[1   1];
  81 %! A=[1 0.5];
  82 %! [sos,g] = tf2sos(B,A);
  83 %! [Bh,Ah] = sos2tf(sos,g);
  84 %! assert({Bh,Ah},{B,A},10*eps);
  85
  86 %!test
  87 %! B=[1 0 0 0 0   1];
  88 %! A=[1 0 0 0 0 0.9];
  89 %! [sos,g] = tf2sos(B,A);
  90 %! [Bh,Ah] = sos2tf(sos,g);
  91 %! assert({Bh,Ah},{B,A},100*eps);