Mercurial > hg > ltpda

comparison m-toolbox/classes/+utils/@math/pfallpsymz2.m @ 0:f0afece42f48

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Import.
author Daniele Nicolodi <nicolodi@science.unitn.it>
date Wed, 23 Nov 2011 19:22:13 +0100
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:f0afece42f48
1 % PFALLPSYMZ2 all pass filtering to stabilize TF poles and zeros.
2 %
3 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4 % DESCRIPTION:
5 %
6 % All pass filtering in order to stabilize transfer function poles and
7 % zeros. It inputs a partial fraction expanded discrete model and
8 % outputs a pole-zero minimum phase system
9 %
10 % CALL:
11 %
12 % resp= pfallpz2(ip,mresp,f,fs)
13 %
14 % INPUTS:
15 %
16 % ip: is a struct with fields named poles
17 % mresp: is a vector with functions response
18 % f: is the frequancies vector in (Hz)
19 % fs: is the sampling frequency in (Hz)
20 %
21 % OUTPUTS:
22 %
23 % resp: is the stable functions frequency response
24 %
25 % NOTE:
26 %
27 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
28 % VERSION: $Id: pfallpz.m,v 1.6 2009/06/10 15:47:00 luigi Exp $
29 %
30 % HISTORY: 12-09-2008 L Ferraioli
31 % Creation
32 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
33 function varargout = pfallpsymz2(ip,mresp,f,fs)
34
35 [a,b] = size(ip);
36 if a<b
37 ip = ip.'; % reshape as a column vector
38 end
39
40 [a,b] = size(f);
41 if a<b
42 f = f.'; % reshape as a column vector
43 end
44
45 if isempty(fs)
46 fs = 1;
47 end
48 [a,b] = size(fs);
49 if a ~= b
50 disp(' Fs has to be a number. Only first term will be considered! ')
51 fs = fs(1);
52 end
53
54
55 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
56 Nb = numel(ip);
57 for nn = 1:Nb
58
59 p = ip(nn).poles;
60
61 % stabilizing poles
62 sp = sym(p);
63 unst = abs(p) > 1;
64 sp(unst) = conj(sp(unst));
65
66 pp = sym(p(unst));
67 psp = sp(unst);
68 syms z
69 allpsym = 1;
70 for jj=1:numel(psp)
71 allpsym = allpsym.*((z-pp(jj))./(z*psp(jj)-1));
72 end
73 funcell{nn} = allpsym;
74
75 end
76
77 fullallprsp = 1;
78
79 for nn = 1:Nb
80 symexpr = funcell{nn};
81 nterm = subs(symexpr,z,cos((2*pi/fs).*f) + 1i.*sin((2*pi/fs).*f));
82 % willing to work with columns
83 if size(nterm,2)>1
84 nterm = nterm.';
85 end
86
87 fullallprsp = fullallprsp.*nterm;
88
89 end
90
91 % dballprsp = double(fullallprsp);
92 % rallp = real(dballprsp);
93 % iallp = imag(dballprsp);
94 % ang = angle(dballprsp);
95
96 for kk=1:Nb
97 sresp(:,kk) = mresp(:,kk).*fullallprsp;
98 % resp(:,kk) = mresp(:,kk).*(cos(ang)+1i.*sin(ang));
99 end
100
101 for kk=1:Nb
102 resp(:,kk) = double(sresp(:,kk));
103 end
104
105
106 % output
107 if nargout == 1
108 varargout{1} = resp;
109 else
110 error('Too many output arguments!')
111 end
112 end
113
114
115
116