0001 function [ao,bo] = pzm2ab(pzm, fs)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 disp(sprintf('$$$ converting %s', get(pzm, 'name')))
0013
0014 gain = get(pzm, 'gain');
0015 poles = get(pzm, 'poles');
0016 zeros = get(pzm, 'zeros');
0017 np = length(poles);
0018 nz = length(zeros);
0019
0020 ao = [];
0021 bo = [];
0022
0023
0024 cpoles = [];
0025 for j=1:np
0026 if get(poles(j), 'q') > 0.5
0027 cpoles = [cpoles poles(j)];
0028 end
0029 end
0030 czeros = [];
0031 for j=1:nz
0032 if get(zeros(j), 'q') > 0.5
0033 czeros = [czeros zeros(j)];
0034 end
0035 end
0036
0037 czi = 1;
0038 for j=1:length(cpoles)
0039 if czi <= length(czeros)
0040
0041 p = cpoles(j);
0042 z = czeros(czi);
0043
0044 [ai,bi] = cpolezero(get(p, 'f'), get(p, 'q'), get(z, 'f'), get(z, 'q'), fs);
0045 if ~isempty(ao)>0
0046 [ao,bo] = abcascade(ao,bo,ai,bi);
0047 else
0048 ao = ai;
0049 bo = bi;
0050 end
0051
0052
0053 czi = czi + 1;
0054 end
0055 end
0056
0057 if length(cpoles) > length(czeros)
0058
0059 for j=length(czeros)+1:length(cpoles)
0060 disp('$$$ computing complex pole');
0061 [ai,bi] = cp2iir(cpoles(j), fs);
0062 if ~isempty(ao)>0
0063 [ao,bo] = abcascade(ao,bo,ai,bi);
0064 else
0065 ao = ai;
0066 bo = bi;
0067 end
0068 end
0069 else
0070
0071 for j=length(cpoles)+1:length(czeros)
0072 disp('$$$ computing complex zero');
0073 [ai,bi] = cz2iir(czeros(j), fs);
0074 if ~isempty(ao)>0
0075 [ao,bo] = abcascade(ao,bo,ai,bi);
0076 else
0077 ao = ai;
0078 bo = bi;
0079 end
0080 end
0081 end
0082
0083
0084 for j=1:np
0085 pole = poles(j);
0086 if isnan(get(pole, 'q')) || get(pole, 'q') < 0.5
0087 disp('$$$ computing real pole');
0088 [ai,bi] = rp2iir(pole, fs);
0089 if ~isempty(ao)>0
0090 [ao,bo] = abcascade(ao,bo,ai,bi);
0091 else
0092 ao = ai;
0093 bo = bi;
0094 end
0095 end
0096 end
0097
0098 for j=1:nz
0099 zero = zeros(j);
0100 if isnan(get(zero, 'q')) || get(zero, 'q') < 0.5
0101 disp('$$$ computing real zero');
0102 [ai,bi] = rz2iir(zero, fs);
0103 if ~isempty(ao)>0
0104 [ao,bo] = abcascade(ao,bo,ai,bi);
0105 else
0106 ao = ai;
0107 bo = bi;
0108 end
0109 end
0110 end
0111
0112
0113
0114
0115
0116
0117
0118
0119
0120
0121
0122
0123
0124
0125
0126
0127
0128
0129
0130
0131
0132
0133
0134
0135
0136
0137
0138
0139
0140
0141
0142
0143 ao = ao.*gain;
0144
0145 function [a,b] = cpolezero(pf, pq, zf, zq, fs)
0146
0147
0148
0149
0150
0151
0152
0153
0154
0155
0156 disp('$$$ computing complex pole/zero pair');
0157
0158 wp = pf*2*pi;
0159 wp2 = wp^2;
0160 wz = zf*2*pi;
0161 wz2 = wz^2;
0162
0163 k = 4*fs*fs + 2*wp*fs/pq + wp2;
0164
0165 a(1) = (4*fs*fs + 2*wz*fs/zq + wz2)/k;
0166 a(2) = (2*wz2 - 8*fs*fs)/k;
0167 a(3) = (4*fs*fs - 2*wz*fs/zq + wz2)/k;
0168 b(1) = 1;
0169 b(2) = (2*wp2 - 8*fs*fs)/k;
0170 b(3) = (wp2 + 4*fs*fs - 2*wp*fs/pq)/k;
0171
0172
0173 g = iirdcgain(a,b);
0174 a = a / g;
0175
0176
0177 function g = iirdcgain(a,b)
0178
0179
0180
0181
0182
0183
0184
0185
0186
0187
0188
0189
0190
0191
0192
0193
0194
0195
0196
0197
0198 g = 0;
0199
0200 suma = sum(a);
0201 if(length(b)>1)
0202 sumb = sum(b);
0203 g = suma / sumb;
0204 else
0205 g = suma;
0206 end
0207
0208