/usr/share/octave/packages/control-3.0.0/thiran.m is in octave-control 3.0.0-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 | ## Copyright (C) 2013-2015 Thomas Vasileiou
##
## This file is part of LTI Syncope.
##
## LTI Syncope 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.
##
## LTI Syncope 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 LTI Syncope. If not, see <http://www.gnu.org/licenses/>.
## -*- texinfo -*-
## @deftypefn {Function File} {@var{sys} =} thiran (@var{tau}, @var{tsam})
## Approximation of continuous-time delay using a discrete-time
## allpass Thiran filter.
##
## Thiran filters can approximate continuous-time delays that are
## non-integer multiples of the sampling time (fractional delays).
## This approximation gives a better matching of the phase shift
## between the continuous- and the discrete-time system.
## If there is no fractional part in the delay, then the standard
## discrete-time delay representation is used.
##
## @strong{Inputs}
## @table @var
## @item tau
## A continuous-time delay, given in time units (seconds).
## @item tsam
## The sampling time of the resulting Thiran filter.
## @end table
##
## @strong{Outputs}
## @table @var
## @item sys
## Transfer function model of the resulting filter.
## The order of the filter is determined automatically.
## @end table
##
## @strong{Example}
## @example
## @group
## octave:1> sys = thiran (1.33, 0.5)
##
## Transfer function 'sys' from input 'u1' to output ...
##
## 0.003859 z^3 - 0.03947 z^2 + 0.2787 z + 1
## y1: -----------------------------------------
## z^3 + 0.2787 z^2 - 0.03947 z + 0.003859
##
## Sampling time: 0.5 s
## Discrete-time model.
## @end group
## @end example
## @example
## @group
## octave:2> sys = thiran (1, 0.5)
##
## Transfer function 'sys' from input 'u1' to output ...
##
## 1
## y1: ---
## z^2
##
## Sampling time: 0.5 s
## Discrete-time model.
## @end group
## @end example
##
## @seealso{absorbdelay, pade}
## @end deftypefn
## Author: Thomas Vasileiou <thomas-v@wildmail.com>
## Created: January 2013
## Version: 0.1
function sys = thiran (del, tsam)
## check args
if (nargin != 2)
print_usage ();
endif
if (! issample (del, 0))
error ("thiran: the delay parameter 'tau' must be a non-negative scalar.");
endif
if (! issample (tsam))
error ("thiran: the second parameter 'tsam' is not a valid sampling time.");
endif
if (del == 0)
sys = tf (1);
return;
endif
## find fractional and discrete delay
N = floor (del/tsam + eps); # put eps or sometimes it misses
d = del - N*tsam;
## check if we do need a thiran filter
if (d/tsam < eps)
sys = tf (1, [1, zeros(1, N)], tsam);
else
## make filter order ~ del to ensure stability
N = N + 1; # order of the filter
d = del/tsam;
tmp = ((d-N+(0:N).') * ones (1,N)) ./ (d-N + bsxfun (@plus, 1:N, (0:N).'));
a = horzcat (1, (-1).^(1:N) .* bincoeff (N, 1:N) .* prod (tmp));
sys = tf (fliplr (a), a, tsam);
endif
endfunction
%!shared num, den, expc
%! expc = [1, 0.5294, -0.04813, 0.004159];
%! sys = thiran (2.4, 1);
%! [num, den] = tfdata (sys, "vector");
%!assert (den, expc, 1e-4);
%!assert (num, fliplr (expc), 1e-4);
%!shared num, den
%! sys = thiran (0.5, 0.1);
%! [num, den] = tfdata (sys, "vector");
%!assert (den, [1, 0, 0, 0, 0, 0], eps);
%!assert (num, 1, eps);
%!error (thiran (-1, 1));
%!error (thiran (1, -1));
%!error (thiran ([1 2 3], 1));
|