Mercurial > hg > ltpda

comparison testing/utp_1.1/utps/ao/utp_ao_firwhiten.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 % UTP_AO_PSD a set of UTPs for the ao/firwhiten method
2 %
3 % M Hewitson 06-08-08
4 %
5 % $Id: utp_ao_firwhiten.m,v 1.21 2011/11/16 08:16:08 mauro Exp $
6 %
7
8 % <MethodDescription>
9 %
10 % The firwhiten method of the ao class computes the spectral density of time-series AOs.
11 %
12 % </MethodDescription>
13
14 function results = utp_ao_firwhiten(varargin)
15
16 % Check the inputs
17 if nargin == 0
18
19 % Some keywords
20 class = 'ao';
21 mthd = 'firwhiten';
22
23 results = [];
24 disp('******************************************************');
25 disp(['**** Running UTPs for ' class '/' mthd]);
26 disp('******************************************************');
27
28 % Test AOs
29 [at1,at2,at3,at4,at5,at6,atvec,atmat] = eval(['get_test_objects_' class]);
30
31 % Exception list for the UTPs:
32 [ple1,ple2,ple3,ple4,ple5,ple6] = get_test_ples();
33
34 % Run the tests
35 results = [results utp_01]; % getInfo call
36 results = [results utp_02]; % Vector input
37 results = [results utp_03]; % Matrix input
38 results = [results utp_04]; % List input
39 results = [results utp_05]; % Test with mixed input
40 results = [results utp_06]; % Test history is working
41 results = [results utp_07]; % Test the modify call works
42 results = [results utp_08]; % Test the data shape
43 results = [results utp_09]; % Test with complex plist
44 results = [results utp_10]; % Test the spectral falttening
45 results = [results utp_11(mthd, at1, ple1)]; % Test plotinfo doesn't disappear
46 results = [results utp_12(mthd, at1, ple1)]; % Test errors are cleared
47
48 disp('Done.');
49 disp('******************************************************');
50
51 elseif nargin == 1 % Check for UTP functions
52 if strcmp(varargin{1}, 'isutp')
53 results = 1;
54 else
55 results = 0;
56 end
57 else
58 error('### Incorrect inputs')
59 end
60
61 %% UTP_01
62
63 % <TestDescription>
64 %
65 % Tests that the getInfo call works for this method.
66 %
67 % </TestDescription>
68 function result = utp_01
69
70
71 % <SyntaxDescription>
72 %
73 % Test that the getInfo call works for no sets, all sets, and each set
74 % individually.
75 %
76 % </SyntaxDescription>
77
78 try
79 % <SyntaxCode>
80 % Call for no sets
81 io(1) = eval([class '.getInfo(''' mthd ''', ''None'')']);
82 % Call for all sets
83 io(2) = eval([class '.getInfo(''' mthd ''')']);
84 % Call for each set
85 for kk=1:numel(io(2).sets)
86 io(kk+2) = eval([class '.getInfo(''' mthd ''', ''' io(2).sets{kk} ''')']);
87 end
88 % </SyntaxCode>
89 stest = true;
90 catch err
91 disp(err.message)
92 stest = false;
93 end
94
95 % <AlgoDescription>
96 %
97 % 1) Check that getInfo call returned an minfo object in all cases.
98 % 2) Check that all plists have the correct parameters.
99 %
100 % </AlgoDescription>
101
102 atest = true;
103 if stest
104 % <AlgoCode>
105 % check we have minfo objects
106 if isa(io, 'minfo')
107 prefs = getappdata(0, 'LTPDApreferences');
108 defaultWinType = char(prefs.getMiscPrefs.getDefaultWindow);
109
110 %%% SET 'None'
111 pn = 1;
112 if ~isempty(io(pn).sets), atest = false; end
113 if ~isempty(io(pn).plists), atest = false; end
114 %%% Check all Sets
115 pn = 2;
116 if ~any(strcmpi(io(pn).sets, 'Default')), atest = false; end
117 if numel(io(pn).plists) ~= numel(io(pn).sets), atest = false; end
118 %%%%%%%%%% SET 'Default'
119 pn = 3;
120 if io(pn).plists.nparams ~= 7, atest = false; end
121 % Check key
122 if ~io(pn).plists.isparam('nfft'), atest = false; end
123 if ~io(pn).plists.isparam('bw'), atest = false; end
124 if ~io(pn).plists.isparam('hc'), atest = false; end
125 if ~io(pn).plists.isparam('win'), atest = false; end
126 if ~io(pn).plists.isparam('order'), atest = false; end
127 if ~io(pn).plists.isparam('firwin'), atest = false; end
128 if ~io(pn).plists.isparam('ntaps'), atest = false; end
129 % Check default value
130 if ~isequal(io(pn).plists.find('nfft'), -1), atest = false; end
131 if ~isequal(io(pn).plists.find('bw'), 20), atest = false; end
132 if ~isequal(io(pn).plists.find('hc'), .8), atest = false; end
133 if ~strcmpi(io(pn).plists.find('win'), defaultWinType), atest = false; end
134 if ~isequal(io(pn).plists.find('order'), -1), atest = false; end
135 if ~strcmpi(io(pn).plists.find('win'), defaultWinType), atest = false; end
136 if ~isequal(io(pn).plists.find('ntaps'), 256), atest = false; end
137 % Check options
138 if ~isequal(io(pn).plists.getOptionsForParam('nfft'), {-1}), atest = false; end
139 if ~isequal(io(pn).plists.getOptionsForParam('bw'), {20}), atest = false; end
140 if ~isequal(io(pn).plists.getOptionsForParam('hc'), {.8}), atest = false; end
141 if ~isequal(io(pn).plists.getOptionsForParam('win'), specwin.getTypes), atest = false; end
142 if ~isequal(io(pn).plists.getOptionsForParam('order'), {-1 0 1 2 3 4 5 6 7 8 9}), atest = false; end
143 if ~isequal(io(pn).plists.getOptionsForParam('firwin'), specwin.getTypes), atest = false; end
144 if ~isequal(io(pn).plists.getOptionsForParam('ntaps'), {256}), atest = false; end
145 end
146 % </AlgoCode>
147 else
148 atest = false;
149 end
150
151 % Return a result structure
152 result = utp_prepare_result(atest, stest, dbstack, mfilename);
153 end % END UTP_01
154
155 %% UTP_02
156
157 % <TestDescription>
158 %
159 % Tests that the firwhiten method works with a vector of AOs as input.
160 %
161 % </TestDescription>
162 function result = utp_02
163
164 % <SyntaxDescription>
165 %
166 % Test that the firwhiten method works for a vector of AOs as input.
167 %
168 % </SyntaxDescription>
169
170 try
171 % <SyntaxCode>
172 avec = [at5 at5 at6];
173 N = 512;
174 fwin = specwin('Hanning', N+1);
175 [out, outf, outxx] = firwhiten(avec, plist('Ntaps', N, 'FIRwin', fwin));
176 % </SyntaxCode>
177 stest = true;
178 catch err
179 disp(err.message)
180 stest = false;
181 end
182
183 % <AlgoDescription>
184 %
185 % 1) Check the number of elements in 'out'
186 % 2) Check the number of filters (outf) and noise-floor estimates (outxx)
187 % 3) Check that each output AO contains the correct data.
188 %
189 % </AlgoDescription>
190
191 atest = true;
192 if stest
193 % <AlgoCode>
194 % Check we have the correct number of outputs
195 if numel(out) ~= numel(avec), atest = false; end
196 if numel(outf) ~= numel(avec), atest = false; end
197 if numel(outxx) ~= numel(avec), atest = false; end
198 % Check the output data
199 prefs = getappdata(0, 'LTPDApreferences');
200 swin = char(prefs.getMiscPrefs.getDefaultWindow);
201 for kk = 1:numel(avec)
202 pl_psd = plist('Nfft', length(avec(kk).y), 'Win', specwin(swin), 'Order', -1, 'Scale', 'ASD');
203 axx = psd(avec(kk), pl_psd);
204 nxx = smoother(axx, plist('width', 20, 'hc', 0.8));
205 w = 1./nxx;
206 ff = mfir(w, plist('Win', fwin, 'N', N));
207 res = filter(avec(kk), ff);
208 if ~isequal(res.x, out(kk).x), atest = false; end
209 if ~isequal(res.y, out(kk).y), atest = false; end
210 end
211 % </AlgoCode>
212 else
213 atest = false;
214 end
215
216 % Return a result structure
217 result = utp_prepare_result(atest, stest, dbstack, mfilename);
218 end % END UTP_02
219
220 %% UTP_03
221
222 % <TestDescription>
223 %
224 % Tests that the firwhiten method works with a matrix of AOs as input.
225 %
226 % </TestDescription>
227 function result = utp_03
228
229 % <SyntaxDescription>
230 %
231 % Test that the firwhiten method works for a matrix of AOs as input.
232 %
233 % </SyntaxDescription>
234
235 try
236 % <SyntaxCode>
237 amat = [at1 at5 at6; at5 at6 at1];
238 N = 512;
239 fwin = specwin('Hanning', N+1);
240 [out, outf, outxx] = firwhiten(amat, plist('Ntaps', N, 'FIRwin', fwin));
241 % </SyntaxCode>
242 stest = true;
243 catch err
244 disp(err.message)
245 stest = false;
246 end
247
248 % <AlgoDescription>
249 %
250 % 1) Check the number of elements in 'out'
251 % 2) Check the number of filters (outf) and noise-floor estimates (outxx)
252 % 3) Check that each output AO contains the correct data.
253 %
254 % </AlgoDescription>
255
256 atest = true;
257 if stest
258 % <AlgoCode>
259 % Check we have the correct number of outputs
260 if numel(out) ~= numel(amat), atest = false; end
261 if numel(outf) ~= numel(amat), atest = false; end
262 if numel(outxx) ~= numel(amat), atest = false; end
263 % Check the output data
264 prefs = getappdata(0, 'LTPDApreferences');
265 swin = char(prefs.getMiscPrefs.getDefaultWindow);
266 for kk = 1:numel(amat)
267 pl_psd = plist('Nfft', length(amat(kk).y), 'Win', specwin(swin), 'Order', -1, 'Scale', 'ASD');
268 axx = psd(amat(kk), pl_psd);
269 nxx = smoother(axx, plist('width', 20, 'hc', 0.8));
270 w = 1./nxx;
271 ff = mfir(w, plist('Win', fwin, 'N', N));
272 res = filter(amat(kk), ff);
273 if ~isequal(res.x, out(kk).x), atest = false; end
274 if ~isequal(res.y, out(kk).y), atest = false; end
275 end
276 % </AlgoCode>
277 else
278 atest = false;
279 end
280
281 % Return a result structure
282 result = utp_prepare_result(atest, stest, dbstack, mfilename);
283 end % END UTP_03
284
285 %% UTP_04
286
287 % <TestDescription>
288 %
289 % Tests that the firwhiten method works with a list of AOs as input.
290 %
291 % </TestDescription>
292 function result = utp_04
293
294 % <SyntaxDescription>
295 %
296 % Test that the firwhiten method works for a list of AOs as input.
297 %
298 % </SyntaxDescription>
299
300 try
301 % <SyntaxCode>
302 N = 512;
303 fwin = specwin('Hanning', N+1);
304 [out, outf, outxx] = firwhiten(at1,at5,at6, plist('Ntaps', N, 'FIRwin', fwin));
305 % </SyntaxCode>
306 stest = true;
307 catch err
308 disp(err.message)
309 stest = false;
310 end
311
312 % <AlgoDescription>
313 %
314 % 1) Check the number of elements in 'out'
315 % 2) Check the number of filters (outf) and noise-floor estimates (outxx)
316 % 3) Check that each output AO contains the correct data.
317 %
318 % </AlgoDescription>
319
320 atest = true;
321 aoin = [at1, at5, at6];
322 if stest
323 % <AlgoCode>
324 % Check we have the correct number of outputs
325 if numel(out) ~= 3, atest = false; end
326 if numel(outf) ~= 3, atest = false; end
327 if numel(outxx) ~= 3, atest = false; end
328 % Check the output data
329 prefs = getappdata(0, 'LTPDApreferences');
330 swin = char(prefs.getMiscPrefs.getDefaultWindow);
331 for kk = 1:numel(aoin)
332 pl_psd = plist('Nfft', length(aoin(kk).y), 'Win', specwin(swin), 'Order', -1, 'Scale', 'ASD');
333 axx = psd(aoin(kk), pl_psd);
334 nxx = smoother(axx, plist('width', 20, 'hc', 0.8));
335 w = 1./nxx;
336 ff = mfir(w, plist('Win', fwin, 'N', N));
337 res = filter(aoin(kk), ff);
338 if ~isequal(res.x, out(kk).x), atest = false; end
339 if ~isequal(res.y, out(kk).y), atest = false; end
340 end
341 % </AlgoCode>
342 else
343 atest = false;
344 end
345
346 % Return a result structure
347 result = utp_prepare_result(atest, stest, dbstack, mfilename);
348 end % END UTP_04
349
350 %% UTP_05
351
352 % <TestDescription>
353 %
354 % Tests that the firwhiten method works with a mix of different shaped AOs as
355 % input.
356 %
357 % </TestDescription>
358 function result = utp_05
359
360 % <SyntaxDescription>
361 %
362 % Test that the firwhiten method works with an input of matrices and vectors
363 % and single AOs.
364 %
365 % </SyntaxDescription>
366
367 try
368 % <SyntaxCode>
369 N = 512;
370 fwin = specwin('Hanning', N+1);
371 pl = plist('Ntaps', N, 'FIRwin', fwin);
372 [out, outf, outxx] = firwhiten(at1,[at5 at6],at5,[at5 at1; at6 at1],at6, pl);
373 % </SyntaxCode>
374 stest = true;
375 catch err
376 disp(err.message)
377 stest = false;
378 end
379
380 % <AlgoDescription>
381 %
382 % 1) Check the number of elements in 'out'
383 % 2) Check the number of filters (outf) and noise-floor estimates (outxx)
384 % 3) Check that each output AO contains the correct data.
385 %
386 % </AlgoDescription>
387
388 atest = true;
389 aoin = [at1, reshape([at5 at6], 1, []), at5, reshape([at5 at1; at6 at1], 1, []), at6];
390 if stest
391 % <AlgoCode>
392 % Check we have the correct number of outputs
393 if numel(out) ~= 9, atest = false; end
394 if numel(outf) ~= 9, atest = false; end
395 if numel(outxx) ~= 9, atest = false; end
396 % Check the output data
397 prefs = getappdata(0, 'LTPDApreferences');
398 swin = char(prefs.getMiscPrefs.getDefaultWindow);
399 for kk = 1:numel(aoin)
400 pl_psd = plist('Nfft', length(aoin(kk).y), 'Win', specwin(swin), 'Order', -1, 'Scale', 'ASD');
401 axx = psd(aoin(kk), pl_psd);
402 nxx = smoother(axx, plist('width', 20, 'hc', 0.8));
403 w = 1./nxx;
404 ff = mfir(w, plist('Win', fwin, 'N', N));
405 res = filter(aoin(kk), ff);
406 if ~isequal(res.x, out(kk).x), atest = false; end
407 if ~isequal(res.y, out(kk).y), atest = false; end
408 end
409 % </AlgoCode>
410 else
411 atest = false;
412 end
413
414 % Return a result structure
415 result = utp_prepare_result(atest, stest, dbstack, mfilename);
416 end % END UTP_05
417
418 %% UTP_06
419
420 % <TestDescription>
421 %
422 % Tests that the firwhiten method properly applies history.
423 %
424 % </TestDescription>
425 function result = utp_06
426
427 % <SyntaxDescription>
428 %
429 % Test that the result of applying the firwhiten method can be processed back
430 % to an m-file.
431 %
432 % </SyntaxDescription>
433
434 try
435 % <SyntaxCode>
436 N = 512;
437 fwin = specwin('Hanning', N+1);
438 pl = plist('Ntaps', N, 'FIRwin', fwin);
439 out = firwhiten(at5, pl);
440 mout = rebuild(out);
441 % </SyntaxCode>
442 stest = true;
443 catch err
444 disp(err.message)
445 stest = false;
446 end
447
448 % <AlgoDescription>
449 %
450 % 1) Check that the last entry in the history of 'out' corresponds to
451 % 'firwhiten'.
452 % 2) Check that the re-built object is the same object as 'out'.
453 %
454 % </AlgoDescription>
455
456 atest = true;
457 if stest
458 % <AlgoCode>
459 % Check the last step in the history of 'out'
460 if ~strcmp(out.hist.methodInfo.mname, 'firwhiten'), atest = false; end
461 % Check the re-built object
462 if ~eq(mout, out, ple4), atest = false; end
463 % </AlgoCode>
464 else
465 atest = false;
466 end
467
468 % Return a result structure
469 result = utp_prepare_result(atest, stest, dbstack, mfilename);
470 end % END UTP_06
471
472 %% UTP_07
473
474 % <TestDescription>
475 %
476 % Tests that the firwhiten method can modify the input AO.
477 %
478 % </TestDescription>
479 function result = utp_07
480
481 % <SyntaxDescription>
482 %
483 % Test that the firwhiten method can modify the input AO by calling
484 % with no output.
485 %
486 % </SyntaxDescription>
487
488 try
489 % <SyntaxCode>
490 N = 512;
491 fwin = specwin('Hanning', N+1);
492 pl = plist('Ntaps', N, 'FIRwin', fwin);
493 % copy at1 to work with
494 ain = ao(at1);
495 % modify ain
496 aout = ain.firwhiten(pl);
497 ain.firwhiten(pl);
498 % </SyntaxCode>
499 stest = true;
500 catch err
501 disp(err.message)
502 stest = false;
503 end
504
505 % <AlgoDescription>
506 %
507 % 1) Check that 'at1' and 'ain' are now different.
508 % 2) Check that 'ain' is firwhiten(at1).
509 %
510 % </AlgoDescription>
511
512 atest = true;
513 if stest
514 % <AlgoCode>
515 % Check that firwhiten modified the input by comparing to the copy
516 if eq(ao(at1), ain, ple1), atest = false; end
517 % Check that firwhiten doesn't modified the input for the function notation
518 if ~eq(aout, ain, ple1), atest = false; end
519 % Check that the modified input is the firwhiten of the copy
520 % Check the output data of ain
521 prefs = getappdata(0, 'LTPDApreferences');
522 swin = char(prefs.getMiscPrefs.getDefaultWindow);
523 pl_psd = plist('Nfft', length(at1.y), 'Win', specwin(swin), 'Order', -1, 'Scale', 'ASD');
524 axx = psd(at1, pl_psd);
525 nxx = smoother(axx, plist('width', 20, 'hc', 0.8));
526 w = 1./nxx;
527 ff = mfir(w, plist('Win', fwin, 'N', N));
528 res = filter(at1, ff);
529 if ~isequal(res.x, ain.x), atest = false; end
530 if ~isequal(res.y, ain.y), atest = false; end
531 % </AlgoCode>
532 else
533 atest = false;
534 end
535
536 % Return a result structure
537 result = utp_prepare_result(atest, stest, dbstack, mfilename);
538 end % END UTP_07
539
540 %% UTP_08
541
542 % <TestDescription>
543 %
544 % Tests that the firwhiten method keeps the data shape of the input object.
545 %
546 % </TestDescription>
547 function result = utp_08
548
549 % <SyntaxDescription>
550 %
551 % Test that the firwhiten method keeps the data shape of the input object. The
552 % input AO must be an AO with row data and an AO with column data.
553 %
554 % </SyntaxDescription>
555
556 try
557 % <SyntaxCode>
558 N = 512;
559 fwin = specwin('Hanning', N+1);
560 pl = plist('Ntaps', N, 'FIRwin', fwin);
561 out1 = firwhiten(at5, pl);
562 out2 = firwhiten(at6, pl);
563 % </SyntaxCode>
564 stest = true;
565 catch err
566 disp(err.message)
567 stest = false;
568 end
569
570 % <AlgoDescription>
571 %
572 % 1) Check that the shpe of the data doesn't change.
573 %
574 % </AlgoDescription>
575
576 atest = true;
577 if stest
578 % <AlgoCode>
579 % Check the shape of the output data
580 if size(out1.data.y,2) ~= 1, atest = false; end
581 if size(out2.data.y,1) ~= 1, atest = false; end
582 % </AlgoCode>
583 else
584 atest = false;
585 end
586
587 % Return a result structure
588 result = utp_prepare_result(atest, stest, dbstack, mfilename);
589 end % END UTP_08
590
591 %% UTP_09
592
593 % <TestDescription>
594 %
595 % Tests that the firwhiten method with a complex plist.
596 %
597 % </TestDescription>
598 function result = utp_09
599
600 % <SyntaxDescription>
601 %
602 % Test that the result of applying the firwhiten method with a complex plist
603 % can be processed back to a m-file.
604 %
605 % </SyntaxDescription>
606
607 try
608 % <SyntaxCode>
609 N = 512;
610 Nfft = 100;
611 order = 1;
612 fwin = specwin('Hanning', N+1);
613 swin = specwin('Kaiser', 10, 100);
614 pl = plist('Ntaps', N, 'FIRwin', fwin, 'win', swin, 'Nfft', Nfft, 'order', order);
615 [out, outf, outxx] = firwhiten(at5, pl);
616 mout = rebuild(out);
617 % </SyntaxCode>
618 stest = true;
619 catch err
620 disp(err.message)
621 stest = false;
622 end
623
624 % <AlgoDescription>
625 %
626 % 1) Check the output data
627 % 2) Check the output filter
628 % 3) Check the noise-floor estimates
629 % 4) Check that the re-built object is the same object as 'out'.
630 %
631 % </AlgoDescription>
632
633 atest = true;
634 if stest
635 % <AlgoCode>
636 % Check the last step in the history of 'out'
637 if ~strcmp(out.hist.methodInfo.mname, 'firwhiten'), atest = false; end
638 % Check the re-built object
639 if ~eq(mout, out, ple2), atest = false; end
640
641 % </AlgoCode>
642 else
643 atest = false;
644 end
645
646 % Return a result structure
647 result = utp_prepare_result(atest, stest, dbstack, mfilename);
648 end % END UTP_09
649
650 %% UTP_10
651
652 % <TestDescription>
653 %
654 % Test the spectral falttening capability of firwhiten method.
655 %
656 % </TestDescription>
657 function result = utp_10
658
659 % <SyntaxDescription>
660 %
661 % Test that the application of the firwhiten method enhances the
662 % spectral flatness of input data.
663 %
664 % </SyntaxDescription>
665
666 try
667 % <SyntaxCode>
668 % Making test data
669 fs = at5.fs;
670 pzm = pzmodel(2, [pz(0.1, 2) pz(0.5, 1)] , [pz(1) pz(2.55, 2)]);
671 ft = miir(pzm,plist('fs',fs));
672 af = filter(at5, ft); % Colored noise
673
674 % Whitening
675 N = 512;
676 Nfft = 100;
677 order = 1;
678 fwin = specwin('Hanning', N+1);
679 swin = specwin('Kaiser', 10, 100);
680 pl = plist('Ntaps', N, 'FIRwin', fwin, 'win', swin, 'Nfft', Nfft, 'order', order);
681 out = firwhiten(af, pl);
682
683 % </SyntaxCode>
684 stest = true;
685 catch err
686 disp(err.message)
687 stest = false;
688 end
689
690 % <AlgoDescription>
691 %
692 % 1) Calculate PSD of input and whitened data
693 % 2) Compare relative spectral flatness coefficients
694 %
695 % </AlgoDescription>
696
697 atest = true;
698 if stest
699 % <AlgoCode>
700 % Check the output data
701 % pl_psd = plist('Nfft', Nfft, 'Win', swin, 'Order', order, 'Scale', 'PSD');
702 axx = af.psd;
703 awxx = out.psd;
704
705 % Claculating flatness
706 sf1 = utils.math.spflat(axx.data.y);
707 sf2 = utils.math.spflat(awxx.data.y);
708
709 % Checking flatness
710 if sf1>sf2, atest = false; end
711
712 % </AlgoCode>
713 else
714 atest = false;
715 end
716
717 % Return a result structure
718 result = utp_prepare_result(atest, stest, dbstack, mfilename);
719 end % END UTP_10
720
721
722
723 end