Mercurial > hg > ltpda

comparison testing/utp_1.1/utp_fcns/get_test_csd_ao_noisegen2D.m @ 44:409a22968d5e default

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Add unit tests
author Daniele Nicolodi <nicolodi@science.unitn.it>
date Tue, 06 Dec 2011 18:42:11 +0100
parents
children
comparison
equal deleted inserted replaced
43:bc767aaa99a8 44:409a22968d5e
1 % test csd for ao/noisegen2D
2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3 %
4 % DESCRIPTION: Output a model for mdc1 cross spectral density
5 %
6 % CALL:
7 % [TF,CSD] = get_test_csd_ao_noisegen2D(pl)
8 %
9 % PARAMETERS: 'fs' - Sampling frequency
10 % 'f1' - the start frequency
11 % 'f2' - the stop frequency
12 % 'nf' - number of evaluation points
13 % 'scale' - spacing of frequencies: 'lin' or 'log'
14 % or
15 % 'f' - a vector of frequency values or an AO
16 % whereby the x-axis is taken for the frequency values
17 %
18 % OUTPUTS:
19 % CSD - a 2x2 matrix of analysis objects containing CSD
20 % models
21 % TF - a 2x2 matrix of analysis objects containing TF
22 % models
23 %
24 % M-FILE INFO: Get information about this methods by calling
25 % >> ao.getInfo('get_test_csd_ao_noisegen2D')
26 %
27 % Get information about a specified set-plist by calling:
28 % >> ao.getInfo('get_test_csd_ao_noisegen2D', 'none')
29 %
30 % VERSION: $Id: get_test_csd_ao_noisegen2D.m,v 1.4 2011/03/24 19:37:53 ingo Exp $
31 %
32 %
33 % HISTORY: 22-12-2008 L Ferraioli
34 % Creation
35 %
36 %
37 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
38 function varargout = get_test_csd_ao_noisegen2D(varargin)
39
40 % Check if this is a call for parameters
41 if utils.helper.isinfocall(varargin{:})
42 varargout{1} = getInfo(varargin{3});
43 return
44 end
45
46 import utils.const.*
47 utils.helper.msg(msg.MNAME, 'running %s/%s', mfilename('class'), mfilename);
48
49 % Collect input variable names
50 in_names = cell(size(varargin));
51 for ii = 1:nargin,in_names{ii} = inputname(ii);end
52
53 % Collect all AOs and plists
54 [as, ao_invars] = utils.helper.collect_objects(varargin(:), 'ao', in_names);
55 pl = utils.helper.collect_objects(varargin(:), 'plist', in_names);
56
57 % Get frequency vector
58 f = find(pl, 'f');
59 if isa(f, 'ao') && (isa(f.data, 'fsdata') || isa(f.data, 'xydata'))
60 f = f.data.getX;
61 end
62
63 % Get sampling frequency
64 fs = find(pl, 'fs');
65
66 if isempty(f)
67 % Compute from frequency range
68 f1 = find(pl, 'f1');
69 f2 = find(pl, 'f2');
70 ndata = find(pl, 'nf');
71 scale = find(pl, 'scale');
72
73 switch scale
74 case 'lin'
75 f = linspace(f1, f2, ndata);
76 case 'log'
77 f = logspace(log10(f1), log10(f2), ndata);
78 end
79 end
80
81
82 % Models response calculation
83 tf11 = calcTF11(f,fs);
84 tf12 = calcTF12(f,fs);
85 tf21 = calcTF21(f,fs);
86 tf22 = calcTF22(f,fs);
87
88 % CSD calculation
89 csd11 = tf11.*conj(tf11)+tf12.*conj(tf12);
90 csd12 = tf11.*conj(tf21)+tf12.*conj(tf22);
91 csd22 = tf22.*conj(tf22)+tf21.*conj(tf21);
92 csd21 = conj(csd12);
93 csd11.setName('CSD11');
94 csd12.setName('CSD12');
95 csd21.setName('CSD21');
96 csd22.setName('CSD22');
97
98 % Output data
99 if nargout == 1
100 varargout{1} = [csd11 csd12;csd21 csd22];
101 elseif nargout == 2
102 varargout{1} = [csd11 csd12;csd21 csd22];
103 varargout{2} = [tf11 tf12;tf21 tf22];
104 end
105
106 end
107
108 %--------------------------------------------------------------------------
109 % Get Info Object
110 %--------------------------------------------------------------------------
111 function ii = getInfo(varargin)
112 if nargin == 1 && strcmpi(varargin{1}, 'None')
113 sets = {};
114 pls = [];
115 elseif nargin == 1&& ~isempty(varargin{1}) && ischar(varargin{1})
116 sets{1} = varargin{1};
117 pls = getDefaultPlist(sets{1});
118 else
119 sets = {'List', 'Range'};
120 pls = [];
121 for kk=1:numel(sets)
122 pls = [pls getDefaultPlist(sets{kk})];
123 end
124 end
125 % Build info object
126 ii = minfo(mfilename, 'ao', '', utils.const.categories.mdc01, '$Id: get_test_csd_ao_noisegen2D.m,v 1.4 2011/03/24 19:37:53 ingo Exp $', sets, pls);
127 ii.setModifier(false);
128 end
129
130 %--------------------------------------------------------------------------
131 % Get Default Plist
132 %--------------------------------------------------------------------------
133 function plo = getDefaultPlist(set)
134 switch set
135 case 'List'
136 plo = plist('fs', 10, 'f', []);
137 case 'Range'
138 plo = plist('fs', 10, 'f1', -1, 'f2', -1, 'nf', 1000, 'scale', 'log');
139 otherwise
140 error('### Unknown set [%s] for a default list.', set);
141 end
142 end
143
144 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
145 % LOCAL FUNCTIONS
146 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
147 function tf11 = calcTF11(f,fs)
148
149
150
151 % Model Stefano TF11 coefficients
152 dRes11 = [2.44554138162509e-011 - 1.79482547894083e-011i;
153 2.44554138162509e-011 + 1.79482547894083e-011i;
154 2.66402334803101e-009 + 1.1025122049153e-009i;
155 2.66402334803101e-009 - 1.1025122049153e-009i;
156 -7.3560293387644e-009;
157 -1.82811618589835e-009 - 1.21803627800855e-009i;
158 -1.82811618589835e-009 + 1.21803627800855e-009i;
159 1.16258677367555e-009;
160 1.65216557639319e-016;
161 -1.78092396888606e-016;
162 -2.80420398962379e-017;
163 9.21305973049041e-013 - 8.24686706827269e-014i;
164 9.21305973049041e-013 + 8.24686706827269e-014i;
165 5.10730060739905e-010 - 3.76571756625722e-011i;
166 5.10730060739905e-010 + 3.76571756625722e-011i;
167 3.45893698149735e-009;
168 3.98139182134446e-014 - 8.25503935419059e-014i;
169 3.98139182134446e-014 + 8.25503935419059e-014i;
170 -1.40595719147164e-011];
171
172 dPoles11 = [0.843464045655194 - 0.0959986292915475i;
173 0.843464045655194 + 0.0959986292915475i;
174 0.953187595424927 - 0.0190043625473383i;
175 0.953187595424927 + 0.0190043625473383i;
176 0.967176277937188;
177 0.995012027005247 - 0.00268322602801729i;
178 0.995012027005247 + 0.00268322602801729i;
179 0.996564761885673;
180 0.999999366165445;
181 0.999981722418555;
182 0.999921882627659;
183 0.999624431675213 - 0.000813407848742761i;
184 0.999624431675213 + 0.000813407848742761i;
185 0.997312006278751 - 0.00265611346834941i;
186 0.997312006278751 + 0.00265611346834941i;
187 0.990516544257531;
188 0.477796923118318 - 0.311064085401834i;
189 0.477796923118318 + 0.311064085401834i;
190 0];
191
192 dDTerms11 = 0;
193
194 % response calculation
195 pfparams.type = 'disc';
196 pfparams.freq = f;
197 pfparams.fs = fs;
198 pfparams.res = dRes11;
199 pfparams.pol = dPoles11;
200 pfparams.dterm = dDTerms11;
201 pfr = utils.math.pfresp(pfparams);
202 mtf11 = pfr.resp;
203
204 % building AOs
205 tf11 = ao(plist('xvals', f, 'yvals', mtf11, 'fs', fs, 'type', 'fsdata'));
206 % name for this object
207 tf11.setName('TF11');
208
209 end
210
211 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
212 function tf12 = calcTF12(f,fs)
213
214 % Model Stefano TF12
215 dRes12 = [1.44258422208796e-017 + 7.07359428613009e-019i;
216 1.44258422208796e-017 - 7.07359428613009e-019i;
217 -3.4918408053655e-021 - 1.05662874569329e-021i;
218 -3.4918408053655e-021 + 1.05662874569329e-021i;
219 -7.61773292876976e-021;
220 4.84357724603939e-020 + 2.38824204294595e-019i;
221 4.84357724603939e-020 - 2.38824204294595e-019i;
222 -4.07088520945753e-020 - 2.31474543846105e-019i;
223 -4.07088520945753e-020 + 2.31474543846105e-019i;
224 8.73316588658882e-023;
225 -5.21840635377469e-020;
226 1.8461911504859e-023;
227 5.20105247464461e-020;
228 -4.68960092394415e-022;
229 -1.44261407664171e-017 + 6.8922564526833e-019i;
230 -1.44261407664171e-017 - 6.8922564526833e-019i;
231 3.13688133935426e-022];
232
233 dPoles12 = [0.477546340377332 - 0.310830571032376i;
234 0.477546340377332 + 0.310830571032376i;
235 0.99790715414307 - 0.0028490561287024i;
236 0.99790715414307 + 0.0028490561287024i;
237 0.998014205354671 ;
238 0.999585354543332 - 0.000780408757425194i;
239 0.999585354543332 + 0.000780408757425194i;
240 0.99966003029931 - 0.000830944038363768i;
241 0.99966003029931 + 0.000830944038363768i;
242 0.999962770401331 ;
243 0.999981881865521 ;
244 0.999999365763457 ;
245 0.999981706320212 ;
246 0.99992421574188 ;
247 0.477898460791003 + 0.311001926610074i;
248 0.477898460791003 - 0.311001926610074i;
249 0];
250 dDTerms12 = 0;
251
252 % response calculation
253 pfparams.type = 'disc';
254 pfparams.freq = f;
255 pfparams.fs = fs;
256 pfparams.res = dRes12;
257 pfparams.pol = dPoles12;
258 pfparams.dterm = dDTerms12;
259 pfr = utils.math.pfresp(pfparams);
260 mtf12 = pfr.resp;
261
262 % building AOs
263 tf12 = ao(plist('xvals', f, 'yvals', mtf12, 'fs', fs, 'type', 'fsdata'));
264 % name for this object
265 tf12.setName('TF12');
266
267 end
268
269 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
270 function tf21 = calcTF21(f,fs)
271
272 % Model Stefano Tf21
273 dRes21 = [-1.80035241582968e-016 + 1.99543917791863e-015i;
274 -1.80035241582968e-016 - 1.99543917791863e-015i;
275 -1.85590889333759e-013 - 1.23844418827409e-014i;
276 -1.85590889333759e-013 + 1.23844418827409e-014i;
277 5.03656596876842e-013 ;
278 -2.62470963499904e-013 + 2.30024232938878e-012i;
279 -2.62470963499904e-013 - 2.30024232938878e-012i;
280 -9.83780507870955e-018 ;
281 3.40426735130194e-021 ;
282 9.78322351492755e-018 ;
283 -1.65010934542937e-020 ;
284 2.60918565203438e-015 + 1.0546609464659e-015i;
285 2.60918565203438e-015 - 1.0546609464659e-015i;
286 3.18105585405455e-014 + 2.48839990780042e-013i;
287 3.18105585405455e-014 - 2.48839990780042e-013i;
288 3.23021641947666e-013 ;
289 4.81265000078114e-016 - 3.18269170053848e-017i;
290 4.81265000078114e-016 + 3.18269170053848e-017i;
291 5.16260024128201e-018];
292
293 dPoles21 = [0.872004077421604 - 0.110344282822693i;
294 0.872004077421604 + 0.110344282822693i;
295 0.956884129232757 - 0.0225532091775074i;
296 0.956884129232757 + 0.0225532091775074i;
297 0.966514825697177 ;
298 0.995140550419744 - 0.000755127639524413i;
299 0.995140550419744 + 0.000755127639524413i;
300 0.999981802194393 ;
301 0.99999936576546 ;
302 0.999981722418555 ;
303 0.999921882627659 ;
304 0.999624431675213 + 0.000813407848742761i;
305 0.999624431675213 - 0.000813407848742761i;
306 0.997312006278751 + 0.00265611346834941i;
307 0.997312006278751 - 0.00265611346834941i;
308 0.990516544257531 ;
309 0.477796923118318 + 0.311064085401834i;
310 0.477796923118318 - 0.311064085401834i;
311 0];
312
313 dDTerms21 = 0;
314
315 % response calculation
316 pfparams.type = 'disc';
317 pfparams.freq = f;
318 pfparams.fs = fs;
319 pfparams.res = dRes21;
320 pfparams.pol = dPoles21;
321 pfparams.dterm = dDTerms21;
322 pfr = utils.math.pfresp(pfparams);
323 mtf21 = pfr.resp;
324
325 % building AOs
326 tf21 = ao(plist('xvals', f, 'yvals', mtf21, 'fs', fs, 'type', 'fsdata'));
327 % name for this object
328 tf21.setName('TF21');
329
330 end
331
332 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
333 function tf22 = calcTF22(f,fs)
334
335 % Model Stefano Tf22
336 dRes22 = [1.1284521501259e-014;
337 -3.72133611555879e-014 - 2.08232683444075e-014i;
338 -3.72133611555879e-014 + 2.08232683444075e-014i;
339 9.84930639106637e-014 - 1.46640810672565e-013i;
340 9.84930639106637e-014 + 1.46640810672565e-013i;
341 2.51323684013671e-014 ;
342 -5.64078525288305e-014 ;
343 -1.62476406586366e-014 ;
344 -1.08424815979566e-011 + 8.32328079357669e-012i;
345 -1.08424815979566e-011 - 8.32328079357669e-012i;
346 2.41831559776112e-011];
347
348 dPoles22 = [0.988511243978897;
349 0.997305870640646 + 0.00211760900132725i;
350 0.997305870640646 - 0.00211760900132725i;
351 0.999626453270255 + 0.0008125673525946i;
352 0.999626453270255 - 0.0008125673525946i;
353 0.999999366366222 ;
354 0.999981706320212 ;
355 0.99992421574188 ;
356 0.477898460791003 - 0.311001926610074i;
357 0.477898460791003 + 0.311001926610074i;
358 0];
359
360 dDTerms22 = 0;
361
362 % response calculation
363 pfparams.type = 'disc';
364 pfparams.freq = f;
365 pfparams.fs = fs;
366 pfparams.res = dRes22;
367 pfparams.pol = dPoles22;
368 pfparams.dterm = dDTerms22;
369 pfr = utils.math.pfresp(pfparams);
370 mtf22 = pfr.resp;
371
372 % building AOs
373 tf22 = ao(plist('xvals', f, 'yvals', mtf22, 'fs', fs, 'type', 'fsdata'));
374 % name for this object
375 tf22.setName('TF22');
376
377 end