This file is indexed.

/usr/share/octave/packages/signal-1.3.2/zp2sos.m is in octave-signal 1.3.2-5.

This file is owned by root:root, with mode 0o644.

The actual contents of the file can be viewed below.

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
## Copyright (C) 2005 Julius O. Smith III <jos@ccrma.stanford.edu>
##
## This program 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.
##
## This program 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
## this program; if not, see <http://www.gnu.org/licenses/>.

## -*- texinfo -*-
## @deftypefn  {Function File} {[@var{sos}, @var{g}] =} zp2sos (@var{z})
## @deftypefnx {Function File} {[@var{sos}, @var{g}] =} zp2sos (@var{z}, @var{p})
## @deftypefnx {Function File} {[@var{sos}, @var{g}] =} zp2sos (@var{z}, @var{p}, @var{k})
## @deftypefnx {Function File} {@var{sos} =} zp2sos (@dots{})
## Convert filter poles and zeros to second-order sections.
##
## INPUTS:
## @itemize
## @item
## @var{z} = column-vector containing the filter zeros
## @item
## @var{p} = column-vector containing the filter poles
## @item
## @var{k} = overall filter gain factor
## If not given the gain is assumed to be 1.
## @end itemize
##
## RETURNED:
## @itemize
## @item
## @var{sos} = matrix of series second-order sections, one per row:
## @example
## @var{sos} = [@var{B1}.' @var{A1}.'; ...; @var{BN}.' @var{AN}.']
## @end example
## where
## @code{@var{B1}.' = [b0 b1 b2] and @var{A1}.' = [1 a1 a2]} for
## section 1, etc.  The b0 entry must be nonzero for each section.
## See @code{filter} for documentation of the second-order direct-form filter
## coefficients @var{B}i and %@var{A}i, i=1:N.
##
## @item
## @var{g} is the overall gain factor that effectively scales
## any one of the @var{B}i vectors.
## @end itemize
##
## If called with only one output argument, the overall filter gain is
## applied to the first second-order section in the matrix @var{sos}.
##
## EXAMPLE:
## @example
##   [z, p, k] = tf2zp ([1 0 0 0 0 1], [1 0 0 0 0 .9]);
##   [sos, g] = zp2sos (z, p, k)
##
## sos =
##    1.0000    0.6180    1.0000    1.0000    0.6051    0.9587
##    1.0000   -1.6180    1.0000    1.0000   -1.5843    0.9587
##    1.0000    1.0000         0    1.0000    0.9791         0
##
## g =
##     1
## @end example
##
## @seealso{sos2pz, sos2tf, tf2sos, zp2tf, tf2zp}
## @end deftypefn

function [sos,g] = zp2sos(z,p,k)

  if nargin<3, k=1; endif
  if nargin<2, p=[]; endif

  [zc,zr] = cplxreal(z(:));
  [pc,pr] = cplxreal(p(:));

  ## zc,zr,pc,pr

  nzc=length(zc);
  npc=length(pc);

  nzr=length(zr);
  npr=length(pr);

  ## Pair up real zeros:
  if nzr
    if mod(nzr,2)==1, zr=[zr;0]; nzr=nzr+1; endif
    nzrsec = nzr/2;
    zrms = -zr(1:2:nzr-1)-zr(2:2:nzr);
    zrp = zr(1:2:nzr-1).*zr(2:2:nzr);
  else
    nzrsec = 0;
  endif

  ## Pair up real poles:
  if npr
    if mod(npr,2)==1, pr=[pr;0]; npr=npr+1; endif
    nprsec = npr/2;
    prms = -pr(1:2:npr-1)-pr(2:2:npr);
    prp = pr(1:2:npr-1).*pr(2:2:npr);
  else
    nprsec = 0;
  endif

  nsecs = max(nzc+nzrsec,npc+nprsec);

  ## Convert complex zeros and poles to real 2nd-order section form:
  zcm2r = -2*real(zc);
  zca2 = abs(zc).^2;
  pcm2r = -2*real(pc);
  pca2 = abs(pc).^2;

  sos = zeros(nsecs,6);
  sos(:,1) = ones(nsecs,1); # all 2nd-order polynomials are monic
  sos(:,4) = ones(nsecs,1);

  nzrl=nzc+nzrsec; # index of last real zero section
  nprl=npc+nprsec; # index of last real pole section

  for i=1:nsecs

    if i<=nzc # lay down a complex zero pair:
      sos(i,2:3) = [zcm2r(i) zca2(i)];
    elseif i<=nzrl # lay down a pair of real zeros:
      sos(i,2:3) = [zrms(i-nzc) zrp(i-nzc)];
    endif

    if i<=npc # lay down a complex pole pair:
      sos(i,5:6) = [pcm2r(i) pca2(i)];
    elseif i<=nprl # lay down a pair of real poles:
      sos(i,5:6) = [prms(i-npc) prp(i-npc)];
    endif
  endfor

  ## If no output argument for the overall gain, combine it into the
  ## first section.
  if (nargout < 2)
    sos(1,1:3) *= k;
  else
    g = k;
  endif

endfunction

%!test
%! B=[1 0 0 0 0 1]; A=[1 0 0 0 0 .9];
%! [z,p,k] = tf2zp(B,A);
%! [sos,g] = zp2sos(z,p,k);
%! [Bh,Ah] = sos2tf(sos,g);
%! assert({Bh,Ah},{B,A},100*eps);