Mercurial > hg > ltpda

comparison testing/utp_1.1/utps/ao/utp_ao_lpsd.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_LPSD a set of UTPs for the ao/lpsd method
2 %
3 % M Hewitson 06-08-08
4 %
5 % $Id: utp_ao_lpsd.m,v 1.30 2011/05/22 21:21:30 mauro Exp $
6 %
7
8 % <MethodDescription>
9 %
10 % The lpsd method of the ao class computes the spectral density of
11 % time-series AOs on a log-frequency axis.
12 %
13 % </MethodDescription>
14
15 function results = utp_ao_lpsd(varargin)
16
17 % Check the inputs
18 if nargin == 0
19
20 % Some keywords
21 class = 'ao';
22 mthd = 'lpsd';
23
24 results = [];
25 disp('******************************************************');
26 disp(['**** Running UTPs for ' class '/' mthd]);
27 disp('******************************************************');
28
29 % Test AOs
30 [at1,at2,at3,at4,at5,at6,atvec,atmat] = eval(['get_test_objects_' class]);
31
32 % Exception list for the UTPs:
33 [ple1,ple2,ple3,ple4,ple5,ple6] = get_test_ples();
34
35 % Get default window from the preferences
36 prefs = getappdata(0, 'LTPDApreferences');
37 defaultWinType = char(prefs.getMiscPrefs.getDefaultWindow);
38
39 % Run the tests
40 results = [results utp_01]; % getInfo call
41 results = [results utp_02]; % Vector input
42 results = [results utp_03]; % Matrix input
43 results = [results utp_04]; % List input
44 results = [results utp_05]; % Test with mixed input
45 results = [results utp_06]; % Test history is working
46 results = [results utp_07]; % Test output of the data
47 results = [results utp_08]; % Test against MATLAB's fft
48 results = [results utp_11(mthd, at1, ple1)]; % Test plotinfo doesn't disappear
49 results = [results utp_17]; % Test units handling: PSD
50 results = [results utp_18]; % Test units handling: ASD
51 results = [results utp_19]; % Test units handling: PS
52 results = [results utp_20]; % Test units handling: AS
53 results = [results utp_41]; % Test different windows: Rectangular
54 results = [results utp_42]; % Test different windows: BH92
55 results = [results utp_43]; % Test different windows: Hamming
56 results = [results utp_44]; % Test different windows: Hanning
57 results = [results utp_45]; % Test different windows: Bartlett
58 results = [results utp_46]; % Test different windows: Nuttall3
59 results = [results utp_47]; % Test different windows: Kaiser psll = [random]
60 results = [results utp_48]; % Test different windows: Kaiser psll = [default]
61 results = [results utp_49]; % Test different windows: Nuttall4
62 results = [results utp_50]; % Test different windows: SFT3F
63
64 disp('Done.');
65 disp('******************************************************');
66
67 elseif nargin == 1 % Check for UTP functions
68 if strcmp(varargin{1}, 'isutp')
69 results = 1;
70 else
71 results = 0;
72 end
73 else
74 error('### Incorrect inputs')
75 end
76
77 %% UTP_01
78
79 % <TestDescription>
80 %
81 % Tests that the getInfo call works for this method.
82 %
83 % </TestDescription>
84 function result = utp_01
85
86
87 % <SyntaxDescription>
88 %
89 % Test that the getInfo call works for no sets, all sets, and each set
90 % individually.
91 %
92 % </SyntaxDescription>
93
94 % <SyntaxCode>
95 try
96 % Call for no sets
97 io(1) = eval([class '.getInfo(''' mthd ''', ''None'')']);
98 % Call for all sets
99 io(2) = eval([class '.getInfo(''' mthd ''')']);
100 % Call for each set
101 for kk=1:numel(io(2).sets)
102 io(kk+2) = eval([class '.getInfo(''' mthd ''', ''' io(2).sets{kk} ''')']);
103 end
104 stest = true;
105 catch err
106 disp(err.message)
107 stest = false;
108 end
109 % </SyntaxCode>
110
111 % <AlgoDescription>
112 %
113 % 1) Check that getInfo call returned an minfo object in all cases.
114 % 2) Check that all plists have the correct parameters.
115 %
116 % </AlgoDescription>
117
118 % <AlgoCode>
119 atest = true;
120 if stest
121 % check we have minfo objects
122 if isa(io, 'minfo')
123
124 %%% SET 'None'
125 if ~isempty(io(1).sets), atest = false; end
126 if ~isempty(io(1).plists), atest = false; end
127 %%% Check all Sets
128 if ~any(strcmpi(io(2).sets, 'Default')), atest = false; end
129 if numel(io(2).plists) ~= numel(io(2).sets), atest = false; end
130 %%%%%%%%%% SET 'Default'
131 if io(3).plists.nparams ~= 10, atest = false; end
132 % Check key
133 if ~io(3).plists.isparam('kdes'), atest = false; end
134 if ~io(3).plists.isparam('jdes'), atest = false; end
135 if ~io(3).plists.isparam('lmin'), atest = false; end
136 if ~io(3).plists.isparam('win'), atest = false; end
137 if ~io(3).plists.isparam('olap'), atest = false; end
138 if ~io(3).plists.isparam('scale'), atest = false; end
139 if ~io(3).plists.isparam('order'), atest = false; end
140 if ~io(3).plists.isparam('psll'), atest = false; end
141 if ~io(3).plists.isparam('times'), atest = false; end
142 if ~io(3).plists.isparam('split'), atest = false; end
143 % Check default value
144 if ~isequal(io(3).plists.find('kdes'), 100), atest = false; end
145 if ~isequal(io(3).plists.find('jdes'), 1000), atest = false; end
146 if ~isequal(io(3).plists.find('lmin'), 0), atest = false; end
147 if ~strcmpi(io(3).plists.find('win'), defaultWinType), atest = false; end
148 if ~isequal(io(3).plists.find('olap'), -1), atest = false; end
149 if ~isequal(io(3).plists.find('scale'), 'PSD'), atest = false; end
150 if ~isequal(io(3).plists.find('order'), 0), atest = false; end
151 if ~isequal(io(3).plists.find('psll'), 200), atest = false; end
152 if ~isEmptyDouble(io(3).plists.find('times')), atest = false; end
153 if ~isEmptyDouble(io(3).plists.find('split')), atest = false; end
154 % Check options
155 if ~isequal(io(3).plists.getOptionsForParam('kdes'), {100}), atest = false; end
156 if ~isequal(io(3).plists.getOptionsForParam('jdes'), {1000}), atest = false; end
157 if ~isequal(io(3).plists.getOptionsForParam('lmin'), {0}), atest = false; end
158 if ~isequal(io(3).plists.getOptionsForParam('win'), specwin.getTypes), atest = false; end
159 if ~isequal(io(3).plists.getOptionsForParam('olap'), {-1}), atest = false; end
160 if ~isequal(io(3).plists.getOptionsForParam('scale'), {'PSD', 'ASD', 'PS', 'AS'}), atest = false; end
161 if ~isequal(io(3).plists.getOptionsForParam('order'), {-1 0 1 2 3 4 5 6 7 8 9}), atest = false; end
162 if ~isequal(io(3).plists.getOptionsForParam('psll'), {200}), atest = false; end
163 if ~isequal(io(3).plists.getOptionsForParam('times'), {[]}), atest = false; end
164 if ~isequal(io(3).plists.getOptionsForParam('split'), {[]}), atest = false; end
165 end
166 else
167 atest = false;
168 end
169 % </AlgoCode>
170
171 % Return a result structure
172 result = utp_prepare_result(atest, stest, dbstack, mfilename);
173 end % END UTP_01
174
175 %% UTP_02
176
177 % <TestDescription>
178 %
179 % Tests that the lpsd method works with a vector of AOs as input.
180 %
181 % </TestDescription>
182 function result = utp_02
183
184 % <SyntaxDescription>
185 %
186 % Test that the lpsd method works for a vector of AOs as input.
187 %
188 % </SyntaxDescription>
189
190 % <SyntaxCode>
191 try
192 avec = [at1 at5 at6];
193 % Vector output
194 out = lpsd(avec);
195 % List output
196 [out1, out2, out3] = lpsd(avec);
197 stest = true;
198 catch err
199 disp(err.message)
200 stest = false;
201 end
202 % </SyntaxCode>
203
204 % <AlgoDescription>
205 %
206 % 1) Check that the number of elements in 'out' is the same as in the input.
207 % 2) Check that each output object contains the correct values.
208 %
209 % </AlgoDescription>
210
211 % <AlgoCode>
212 atest = true;
213 if stest
214 % Check we have the correct number of outputs
215 if numel(out) ~= numel(avec), atest = false; end
216 % Check we have the correct values in the outputs
217 % Vector output
218 for kk = 1:numel(out)
219 if ~eq(out(kk), lpsd(avec(kk)), ple1), atest = false; end
220 end
221 % List output
222 if ~eq(out1, lpsd(avec(1)), ple1), atest = false; end
223 if ~eq(out2, lpsd(avec(2)), ple1), atest = false; end
224 if ~eq(out3, lpsd(avec(3)), ple1), atest = false; end
225 else
226 atest = false;
227 end
228 % </AlgoCode>
229
230 % Return a result structure
231 result = utp_prepare_result(atest, stest, dbstack, mfilename);
232 end % END UTP_02
233
234 %% UTP_03
235
236 % <TestDescription>
237 %
238 % Tests that the lpsd method works with a matrix of AOs as input.
239 %
240 % </TestDescription>
241 function result = utp_03
242
243 % <SyntaxDescription>
244 %
245 % Test that the lpsd method works for a matrix of AOs as input.
246 %
247 % </SyntaxDescription>
248
249 % <SyntaxCode>
250 try
251 amat = [at1 at5 at6; at5 at6 at1];
252 % Vector output
253 out = lpsd(amat);
254 % List output
255 [out1, out2, out3, out4, out5, out6] = lpsd(amat);
256 stest = true;
257 catch err
258 disp(err.message)
259 stest = false;
260 end
261 % </SyntaxCode>
262
263 % <AlgoDescription>
264 %
265 % 1) Check that the number of elements in 'out' is the same as in the input.
266 % 2) Check that each output object contains the correct values.
267 %
268 % </AlgoDescription>
269
270 % <AlgoCode>
271 atest = true;
272 if stest
273 % Check we have the correct number of outputs
274 if numel(out) ~= numel(amat), atest = false; end
275 % Check we have the correct values in the outputs
276 % Vector output
277 for kk = 1:numel(out)
278 if ~eq(out(kk), lpsd(amat(kk)), ple1), atest = false; end
279 end
280 % List output
281 if ~eq(out1, lpsd(amat(1)), ple1), atest = false; end
282 if ~eq(out2, lpsd(amat(2)), ple1), atest = false; end
283 if ~eq(out3, lpsd(amat(3)), ple1), atest = false; end
284 if ~eq(out4, lpsd(amat(4)), ple1), atest = false; end
285 if ~eq(out5, lpsd(amat(5)), ple1), atest = false; end
286 if ~eq(out6, lpsd(amat(6)), ple1), atest = false; end
287 else
288 atest = false;
289 end
290 % </AlgoCode>
291
292 % Return a result structure
293 result = utp_prepare_result(atest, stest, dbstack, mfilename);
294 end % END UTP_03
295
296 %% UTP_04
297
298 % <TestDescription>
299 %
300 % Tests that the lpsd method works with a list of AOs as input.
301 %
302 % </TestDescription>
303 function result = utp_04
304
305 % <SyntaxDescription>
306 %
307 % Test that the lpsd method works for a list of AOs as input.
308 %
309 % </SyntaxDescription>
310
311 % <SyntaxCode>
312 try
313 % Vector output
314 out = lpsd(at1,at5,at6);
315 % List output
316 [out1, out2, out3] = lpsd(at1,at5,at6);
317 stest = true;
318 catch err
319 disp(err.message)
320 stest = false;
321 end
322 % </SyntaxCode>
323
324 % <AlgoDescription>
325 %
326 % 1) Check that the number of elements in 'out' is the same as in the input.
327 % 2) Check that each output AO contains the correct data.
328 %
329 % </AlgoDescription>
330
331 % <AlgoCode>
332 atest = true;
333 if stest
334 % Check we have the correct number of outputs
335 if numel(out) ~= 3, atest = false; end
336 % Check we have the correct values in the outputs
337 % Vector output
338 if ~eq(out(1), lpsd(at1), ple1), atest = false; end
339 if ~eq(out(2), lpsd(at5), ple1), atest = false; end
340 if ~eq(out(3), lpsd(at6), ple1), atest = false; end
341 % List output
342 if ~eq(out1, lpsd(at1), ple1), atest = false; end
343 if ~eq(out2, lpsd(at5), ple1), atest = false; end
344 if ~eq(out3, lpsd(at6), ple1), atest = false; end
345
346 else
347 atest = false;
348 end
349 % </AlgoCode>
350
351 % Return a result structure
352 result = utp_prepare_result(atest, stest, dbstack, mfilename);
353 end % END UTP_04
354
355 %% UTP_05
356
357 % <TestDescription>
358 %
359 % Tests that the lpsd method works with a mix of different shaped AOs as
360 % input.
361 %
362 % </TestDescription>
363 function result = utp_05
364
365 % <SyntaxDescription>
366 %
367 % Test that the lpsd method works with an input of matrices and vectors
368 % and single AOs.
369 %
370 % </SyntaxDescription>
371
372 % <SyntaxCode>
373 try
374 % Vector output
375 out = lpsd(at1,[at5 at6],at5,[at5 at1; at6 at1],at6);
376 % List output
377 [out1, out2, out3, out4, out5, out6, out7, out8, out9] = ...
378 lpsd(at1,[at5 at6],at5,[at5 at1; at6 at1],at6);
379 stest = true;
380 catch err
381 disp(err.message)
382 stest = false;
383 end
384 % </SyntaxCode>
385
386 % <AlgoDescription>
387 %
388 % 1) Check that the number of elements in 'out' is the same as in
389 % input.
390 % 2) Check that each output AO contains the correct data.
391 %
392 % </AlgoDescription>
393
394 % <AlgoCode>
395 atest = true;
396 if stest
397 % Check we have the correct number of outputs
398 if numel(out) ~= 9, atest = false; end
399 % Check we have the correct values in the outputs
400 % Vector output
401 if ~eq(out(1), lpsd(at1), ple1), atest = false; end
402 if ~eq(out(2), lpsd(at5), ple1), atest = false; end
403 if ~eq(out(3), lpsd(at6), ple1), atest = false; end
404 if ~eq(out(4), lpsd(at5), ple1), atest = false; end
405 if ~eq(out(5), lpsd(at5), ple1), atest = false; end
406 if ~eq(out(6), lpsd(at6), ple1), atest = false; end
407 if ~eq(out(7), lpsd(at1), ple1), atest = false; end
408 if ~eq(out(8), lpsd(at1), ple1), atest = false; end
409 if ~eq(out(9), lpsd(at6), ple1), atest = false; end
410 % List output
411 if ~eq(out1, lpsd(at1), ple1), atest = false; end
412 if ~eq(out2, lpsd(at5), ple1), atest = false; end
413 if ~eq(out3, lpsd(at6), ple1), atest = false; end
414 if ~eq(out4, lpsd(at5), ple1), atest = false; end
415 if ~eq(out5, lpsd(at5), ple1), atest = false; end
416 if ~eq(out6, lpsd(at6), ple1), atest = false; end
417 if ~eq(out7, lpsd(at1), ple1), atest = false; end
418 if ~eq(out8, lpsd(at1), ple1), atest = false; end
419 if ~eq(out9, lpsd(at6), ple1), atest = false; end
420
421 else
422 atest = false;
423 end
424 % </AlgoCode>
425
426 % Return a result structure
427 result = utp_prepare_result(atest, stest, dbstack, mfilename);
428 end % END UTP_05
429
430 %% UTP_06
431
432 % <TestDescription>
433 %
434 % Tests that the lpsd method properly applies history.
435 %
436 % </TestDescription>
437 function result = utp_06
438
439 % <SyntaxDescription>
440 %
441 % Test that the result of applying the lpsd method can be processed back
442 % to an m-file.
443 %
444 % </SyntaxDescription>
445
446 % <SyntaxCode>
447 try
448 out = lpsd(at5);
449 mout = rebuild(out);
450 stest = true;
451 catch err
452 disp(err.message)
453 stest = false;
454 end
455 % </SyntaxCode>
456
457 % <AlgoDescription>
458 %
459 % 1) Check that the last entry in the history of 'out' corresponds to
460 % 'lpsd'.
461 % 2) Check that the re-built object is the same object as 'out'.
462 %
463 % </AlgoDescription>
464
465 % <AlgoCode>
466 atest = true;
467 if stest
468 % Check the last step in the history of 'out'
469 if ~strcmp(out.hist.methodInfo.mname, 'lpsd'), atest = false; end
470 % Check the re-built object
471 if ~eq(mout, out, ple2), atest = false; end
472 else
473 atest = false;
474 end
475 % </AlgoCode>
476
477 % Return a result structure
478 result = utp_prepare_result(atest, stest, dbstack, mfilename);
479 end % END UTP_06
480
481 %% UTP_07
482
483 % <TestDescription>
484 %
485 % Check that the lpsd method pass back the output objects to a list of
486 % output variables or to a single variable.
487 %
488 % </TestDescription>
489 function result = utp_07
490
491 % <SyntaxDescription>
492 %
493 % Call the method with a list of output variables and with a single output
494 % variable. Additionaly check that the rebuild method works on the output.
495 %
496 % </SyntaxDescription>
497
498 try
499 % <SyntaxCode>
500 [o1, o2] = lpsd(at5, at6);
501 o3 = lpsd(at5, at6);
502 mout1 = rebuild(o1);
503 mout2 = rebuild(o2);
504 mout3 = rebuild(o3);
505 % </SyntaxCode>
506 stest = true;
507 catch err
508 disp(err.message)
509 stest = false;
510 end
511
512 % <AlgoDescription>
513 %
514 % 1) Check that the output contains the right number of objects
515 % 2) Check that the 'rebuild' method produces the same object as 'out'.
516 %
517 % </AlgoDescription>
518
519 atest = true;
520 if stest
521 % <AlgoCode>
522 % Check the number of outputs
523 if numel(o1) ~=1, atest = false; end
524 if numel(o2) ~=1, atest = false; end
525 if numel(o3) ~=2, atest = false; end
526 % Check the rebuilding of the object
527 if ~eq(o1, mout1, ple2), atest = false; end
528 if ~eq(o2, mout2, ple2), atest = false; end
529 if ~eq(o3, mout3, ple2), atest = false; end
530 % </AlgoCode>
531 else
532 atest = false;
533 end
534
535 % Return a result structure
536 result = utp_prepare_result(atest, stest, dbstack, mfilename);
537 end % END UTP_07
538 %% UTP_08
539
540 % <TestDescription>
541 %
542 % Tests that the lpsd method agrees with MATLAB's pwelch if the frequency
543 % resolution of the latter is changed to match those used in lpsd
544 %
545 %
546 % % </TestDescription>
547 function result = utp_08
548
549 % <SyntaxDescription>
550 %
551 % Test that the applying psd works on a single AO.
552 %
553 % </SyntaxDescription>
554
555 % <SyntaxCode>
556 try
557 % Construct test AO
558 nsecs = 10;
559 fs = 1000;
560 Nfft = 2*fs;
561 pl = plist('nsecs', nsecs, 'fs', fs, 'tsfcn', 'randn(size(t))');
562 a1 = ao(pl);
563 % parameter list for lpsd
564 win = specwin('Hanning', Nfft);
565 pl = plist('Kdes', 10, 'Lmin', 2, 'Jdes', 100, 'Win', win.type, 'Olap', 50, 'Order', -1);
566 % compute lpsd
567 a2 = lpsd(a1, pl);
568 % </SyntaxCode>
569 stest = true;
570 catch err
571 disp(err.message)
572 stest = false;
573 end
574
575 % <AlgoDescription>
576 %
577 % 1) Check that the DFT equations (eq. (3) and (4) in [1]) are fulfilled.
578 % 2) Check that we get the same outputs for each frequency bin when
579 % computing the psd using matlab's fft with the frequency resolution
580 % values retrieved by lpsd (and already tested in (1))
581 %
582 % [1] G. Heinzel, lpsd revisited: ltf, S2-AEI-TN-3052
583 %
584 % </AlgoDescription>
585 atest = true;
586 if stest
587 % <AlgoCode>
588 % get freq. resolution, bin number and seg. length from procinfo
589 b = find(a2.procinfo,'m');
590 r = find(a2.procinfo,'r');
591 L = find(a2.procinfo,'L');
592 % test eq.(3): f = r*b up to numerical accuracy
593 if all(abs(a2.data.x - r.*b) > max(eps(a2.data.x),eps(r.*b)))
594 atest = false;
595 end
596 % test eq.(4): fs = r*L up to numerical accuracy
597 if all(abs(a2.data.fs - r.*L) > max(eps(a2.data.fs),eps(r.*L)))
598 atest = false;
599 end
600
601 % </AlgoCode>
602 else
603 atest = false;
604 end
605 % Return a result structure
606 result = utp_prepare_result(atest, stest, dbstack, mfilename);
607 end % END UTP_08
608
609 %% UTP_17
610
611 % <TestDescription>
612 %
613 % Tests handling of units:
614 % 1) white noise produced from normal pdf, with a given mean value and
615 % sigma (distribution's 1st and 2nd orders)
616 % 2) LPSD of the white noise
617 % 3) compares the units of the input and output
618 %
619 % </TestDescription>
620 function result = utp_17
621
622 % <SyntaxDescription>
623 %
624 % 1) Prepare the test tsdata:
625 % white noise from normal distribution + offset
626 % 2) Assign a random unit
627 % 3) LPSD of the white noise
628 %
629 % </SyntaxDescription>
630
631 % <SyntaxCode>
632 try
633
634 noise_type = 'Normal';
635 win_type = 'BH92';
636 scale = 'PSD';
637
638 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, plist);
639
640 stest = true;
641
642 catch err
643 disp(err.message)
644 stest = false;
645 end
646 % </SyntaxCode>
647
648 % <AlgoDescription>
649 %
650 % 1) Check that (calculated LPSD yunits) equals (input units)^2 / Hz
651 % 2) Check that (calculated LPSD xunits) equals Hz
652 %
653 % </AlgoDescription>
654
655 % <AlgoCode>
656 atest = true;
657
658 if stest
659 if ne(S.yunits, (a.yunits).^2 * unit('Hz^-1')) || ne(S.xunits, unit('Hz'))
660 atest = false;
661 end
662 else
663 atest = false;
664 end
665 % </AlgoCode>
666
667 % Return a result structure
668 result = utp_prepare_result(atest, stest, dbstack, mfilename);
669 end % END UTP_17
670
671 %% UTP_18
672
673 % <TestDescription>
674 %
675 % Tests handling of units:
676 % 1) white noise produced from uniform pdf, with a given mean value and
677 % sigma (distribution's 1st and 2nd orders)
678 % 2) LASD of the white noise
679 % 3) compares the units of the input and output
680 %
681
682 % </TestDescription>
683 function result = utp_18
684
685 % <SyntaxDescription>
686 %
687 % 1) Prepare the test tsdata:
688 % white noise from uniform distribution + offset
689 % 2) Assign a random unit
690 % 3) LASD of the white noise
691 %
692 % </SyntaxDescription>
693
694 % <SyntaxCode>
695 try
696
697 noise_type = 'Uniform';
698 win_type = 'Hamming';
699 scale = 'ASD';
700
701 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, plist);
702
703 stest = true;
704
705 catch err
706 disp(err.message)
707 stest = false;
708 end
709 % </SyntaxCode>
710
711 % <AlgoDescription>
712 %
713 % 1) Check that (calculated LASD yunits) equals (input units) / Hz^(1/2)
714 % 2) Check that (calculated LASD xunits) equals Hz
715 %
716 % </AlgoDescription>
717
718 % <AlgoCode>
719 atest = true;
720
721 if stest
722 if ne(S.yunits, (a.yunits) * unit('Hz^-1/2')) || ne(S.xunits, unit('Hz'))
723 atest = false;
724 end
725 else
726 atest = false;
727 end
728 % </AlgoCode>
729
730 % Return a result structure
731 result = utp_prepare_result(atest, stest, dbstack, mfilename);
732 end % END UTP_18
733
734 %% UTP_19
735
736 % <TestDescription>
737 %
738 % Tests handling of units:
739 % 1) white noise produced from normal pdf, with a given mean value and
740 % sigma (distribution's 1st and 2nd orders)
741 % 2) LPS of the white noise
742 % 3) compares the units of the input and output
743 %
744
745 % </TestDescription>
746 function result = utp_19
747
748 % <SyntaxDescription>
749 %
750 % 1) Prepare the test tsdata:
751 % white noise from normal distribution + offset
752 % 2) Assign a random unit
753 % 3) LPS of the white noise
754 %
755 % </SyntaxDescription>
756
757 % <SyntaxCode>
758 try
759
760 noise_type = 'Normal';
761 win_type = 'Hanning';
762 scale = 'PS';
763
764 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, plist);
765
766 stest = true;
767
768 catch err
769 disp(err.message)
770 stest = false;
771 end
772 % </SyntaxCode>
773
774 % <AlgoDescription>
775 %
776 % 1) Check that (calculated LPS yunits) equals (input units)^2
777 % 2) Check that (calculated LPS xunits) equals Hz
778 %
779 % </AlgoDescription>
780
781 % <AlgoCode>
782 atest = true;
783
784 if stest
785 if ne(S.yunits, (a.yunits).^2) || ne(S.xunits, unit('Hz'))
786 atest = false;
787 end
788 else
789 atest = false;
790 end
791 % </AlgoCode>
792
793 % Return a result structure
794 result = utp_prepare_result(atest, stest, dbstack, mfilename);
795 end % END UTP_19
796
797 %% UTP_20
798
799 % <TestDescription>
800 %
801 % Tests handling of units:
802 % 1) white noise produced from uniform distribution, with a given mean value and
803 % sigma (distribution's 1st and 2nd orders)
804 % 2) LAS of the white noise
805 % 3) compares the units of the input and output
806 %
807
808 % </TestDescription>
809 function result = utp_20
810
811 % <SyntaxDescription>
812 %
813 % 1) Prepare the test tsdata:
814 % white noise from uniform distribution + offset
815 % 2) Assign a random unit
816 % 3) LAS of the white noise
817 %
818 % </SyntaxDescription>
819
820 % <SyntaxCode>
821 try
822
823 noise_type = 'Uniform';
824 win_type = 'Rectangular';
825 scale = 'AS';
826
827 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, plist);
828
829 stest = true;
830
831 catch err
832 disp(err.message)
833 stest = false;
834 end
835 % </SyntaxCode>
836
837 % <AlgoDescription>
838 %
839 % 1) Check that (calculated LAS yunits) equals (input units)
840 % 2) Check that (calculated LAS xunits) equals Hz
841 %
842 % </AlgoDescription>
843
844 % <AlgoCode>
845 atest = true;
846
847 if stest
848 if ne(S.yunits, a.yunits) || ne(S.xunits, unit('Hz'))
849 atest = false;
850 end
851 else
852 atest = false;
853 end
854 % </AlgoCode>
855
856 % Return a result structure
857 result = utp_prepare_result(atest, stest, dbstack, mfilename);
858 end % END UTP_20
859
860 %% UTP_41
861
862 % <TestDescription>
863 %
864 % Tests the effect of windowing:
865 % 1) white noise produced from normal pdf, with:
866 % a given mean value and sigma (distribution's 1st and 2nd order)
867 % 2) lpsd passing the window name (Rectangular)
868 % 3) lpsd passing the window object (Rectangular type)
869 % 4) compares the two psds
870 %
871
872 % </TestDescription>
873 function result = utp_41
874
875 % <SyntaxDescription>
876 %
877 % 1) Prepare the test tsdata:
878 % white noise from normal distribution + offset
879 % 2) Calculate the statistical parameters
880 % 3) Estimate the psd without detrending, Rectangular window (name)
881 % 4) Estimate the psd without detrending, Rectangular window (object)
882 % 5) Compare results
883 %
884 % </SyntaxDescription>
885
886 % <SyntaxCode>
887 try
888
889 noise_type = 'Normal';
890 scale = 'PSD';
891 win_type = 'Rectangular';
892
893 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, ...
894 plist('win_test', true));
895
896 stest = true;
897
898 catch err
899 disp(err.message)
900 stest = false;
901 end
902 % </SyntaxCode>
903
904 % <AlgoDescription>
905 %
906 % 1) Check that calculated psds are identical
907 %
908 % </AlgoDescription>
909
910 % <AlgoCode>
911 atest = true;
912
913 if stest
914 % Compare the psd evaluated with the window name / object
915 if ~eq(S, S1, ple1)
916 atest = false;
917 end
918 else
919 atest = false;
920 end
921 % </AlgoCode>
922
923 % Return a result structure
924 result = utp_prepare_result(atest, stest, dbstack, mfilename);
925 end % END UTP_41
926
927 %% UTP_42
928
929 % <TestDescription>
930 %
931 % Tests the effect of windowing:
932 % 1) white noise produced from normal pdf, with:
933 % a given mean value and sigma (distribution's 1st and 2nd order)
934 % 2) lpsd passing the window name (BH92)
935 % 3) lpsd passing the window object (BH92 type)
936 % 4) compares the two psds
937 %
938
939 % </TestDescription>
940 function result = utp_42
941
942 % <SyntaxDescription>
943 %
944 % 1) Prepare the test tsdata:
945 % white noise from normal distribution + offset
946 % 2) Calculate the statistical parameters
947 % 3) Estimate the psd without detrending, BH92 window (name)
948 % 4) Estimate the psd without detrending, BH92 window (object)
949 % 5) Compare results
950 %
951 % </SyntaxDescription>
952
953 % <SyntaxCode>
954 try
955
956 noise_type = 'Normal';
957 scale = 'PSD';
958 win_type = 'BH92';
959
960 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, ...
961 plist('win_test', true));
962
963 stest = true;
964
965 catch err
966 disp(err.message)
967 stest = false;
968 end
969 % </SyntaxCode>
970
971 % <AlgoDescription>
972 %
973 % 1) Check that calculated psds are identical
974 %
975 % </AlgoDescription>
976
977 % <AlgoCode>
978 atest = true;
979
980 if stest
981 % Compare the psd evaluated with the window name / object
982 if ~eq(S, S1, ple1)
983 atest = false;
984 end
985 else
986 atest = false;
987 end
988 % </AlgoCode>
989
990 % Return a result structure
991 result = utp_prepare_result(atest, stest, dbstack, mfilename);
992 end % END UTP_42
993
994
995 %% UTP_43
996
997 % <TestDescription>
998 %
999 % Tests the effect of windowing:
1000 % 1) white noise produced from normal pdf, with:
1001 % a given mean value and sigma (distribution's 1st and 2nd order)
1002 % 2) lpsd passing the window name (Hamming)
1003 % 3) lpsd passing the window object (Hamming type)
1004 % 4) compares the two psds
1005 %
1006
1007 % </TestDescription>
1008 function result = utp_43
1009
1010 % <SyntaxDescription>
1011 %
1012 % 1) Prepare the test tsdata:
1013 % white noise from normal distribution + offset
1014 % 2) Calculate the statistical parameters
1015 % 3) Estimate the psd without detrending, Hamming window (name)
1016 % 4) Estimate the psd without detrending, Hamming window (object)
1017 % 5) Compare results
1018 %
1019 % </SyntaxDescription>
1020
1021 % <SyntaxCode>
1022 try
1023
1024 noise_type = 'Normal';
1025 scale = 'PSD';
1026 win_type = 'Hamming';
1027
1028 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, ...
1029 plist('win_test', true));
1030
1031 stest = true;
1032
1033 catch err
1034 disp(err.message)
1035 stest = false;
1036 end
1037 % </SyntaxCode>
1038
1039 % <AlgoDescription>
1040 %
1041 % 1) Check that calculated psds are identical
1042 %
1043 % </AlgoDescription>
1044
1045 % <AlgoCode>
1046 atest = true;
1047
1048 if stest
1049 % Compare the psd evaluated with the window name / object
1050 if ~eq(S, S1, ple1)
1051 atest = false;
1052 end
1053 else
1054 atest = false;
1055 end
1056 % </AlgoCode>
1057
1058 % Return a result structure
1059 result = utp_prepare_result(atest, stest, dbstack, mfilename);
1060 end % END UTP_43
1061
1062 %% UTP_44
1063
1064 % <TestDescription>
1065 %
1066 % Tests the effect of windowing:
1067 % 1) white noise produced from normal pdf, with:
1068 % a given mean value and sigma (distribution's 1st and 2nd order)
1069 % 2) lpsd passing the window name (Hanning)
1070 % 3) lpsd passing the window object (Hanning type)
1071 % 4) compares the two psds
1072 %
1073
1074 % </TestDescription>
1075 function result = utp_44
1076
1077 % <SyntaxDescription>
1078 %
1079 % 1) Prepare the test tsdata:
1080 % white noise from normal distribution + offset
1081 % 2) Calculate the statistical parameters
1082 % 3) Estimate the psd without detrending, Hanning window (name)
1083 % 4) Estimate the psd without detrending, Hanning window (object)
1084 % 5) Compare results
1085 %
1086 % </SyntaxDescription>
1087
1088 % <SyntaxCode>
1089 try
1090
1091 noise_type = 'Normal';
1092 scale = 'PSD';
1093 win_type = 'Hanning';
1094
1095 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, ...
1096 plist('win_test', true));
1097
1098 stest = true;
1099
1100 catch err
1101 disp(err.message)
1102 stest = false;
1103 end
1104 % </SyntaxCode>
1105
1106 % <AlgoDescription>
1107 %
1108 % 1) Check that calculated psds are identical
1109 %
1110 % </AlgoDescription>
1111
1112 % <AlgoCode>
1113 atest = true;
1114
1115 if stest
1116 % Compare the psd evaluated with the window name / object
1117 if ~eq(S, S1, ple1)
1118 atest = false;
1119 end
1120 else
1121 atest = false;
1122 end
1123 % </AlgoCode>
1124
1125 % Return a result structure
1126 result = utp_prepare_result(atest, stest, dbstack, mfilename);
1127 end % END UTP_44
1128
1129 %% UTP_43
1130
1131 % <TestDescription>
1132 %
1133 % Tests the effect of windowing:
1134 % 1) white noise produced from normal pdf, with:
1135 % a given mean value and sigma (distribution's 1st and 2nd order)
1136 % 2) lpsd passing the window name (Bartlett)
1137 % 3) lpsd passing the window object (Bartlett type)
1138 % 4) compares the two psds
1139 %
1140
1141 % </TestDescription>
1142 function result = utp_45
1143
1144 % <SyntaxDescription>
1145 %
1146 % 1) Prepare the test tsdata:
1147 % white noise from normal distribution + offset
1148 % 2) Calculate the statistical parameters
1149 % 3) Estimate the psd without detrending, Bartlett window (name)
1150 % 4) Estimate the psd without detrending, Bartlett window (object)
1151 % 5) Compare results
1152 %
1153 % </SyntaxDescription>
1154
1155 % <SyntaxCode>
1156 try
1157
1158 noise_type = 'Normal';
1159 scale = 'PSD';
1160 win_type = 'Bartlett';
1161
1162 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, ...
1163 plist('win_test', true));
1164
1165 stest = true;
1166
1167 catch err
1168 disp(err.message)
1169 stest = false;
1170 end
1171 % </SyntaxCode>
1172
1173 % <AlgoDescription>
1174 %
1175 % 1) Check that calculated psds are identical
1176 %
1177 % </AlgoDescription>
1178
1179 % <AlgoCode>
1180 atest = true;
1181
1182 if stest
1183 % Compare the psd evaluated with the window name / object
1184 if ~eq(S, S1, ple1)
1185 atest = false;
1186 end
1187 else
1188 atest = false;
1189 end
1190 % </AlgoCode>
1191
1192 % Return a result structure
1193 result = utp_prepare_result(atest, stest, dbstack, mfilename);
1194 end % END UTP_45
1195
1196 %% UTP_46
1197
1198 % <TestDescription>
1199 %
1200 % Tests the effect of windowing:
1201 % 1) white noise produced from normal pdf, with:
1202 % a given mean value and sigma (distribution's 1st and 2nd order)
1203 % 2) lpsd passing the window name (Nuttall3)
1204 % 3) lpsd passing the window object (Nuttall3 type)
1205 % 4) compares the two psds
1206 %
1207
1208 % </TestDescription>
1209 function result = utp_46
1210
1211 % <SyntaxDescription>
1212 %
1213 % 1) Prepare the test tsdata:
1214 % white noise from normal distribution + offset
1215 % 2) Calculate the statistical parameters
1216 % 3) Estimate the psd without detrending, Nuttall3 window (name)
1217 % 4) Estimate the psd without detrending, Nuttall3 window (object)
1218 % 5) Compare results
1219 %
1220 % </SyntaxDescription>
1221
1222 % <SyntaxCode>
1223 try
1224
1225 noise_type = 'Normal';
1226 scale = 'PSD';
1227 win_type = 'Nuttall3';
1228
1229 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, ...
1230 plist('win_test', true));
1231
1232 stest = true;
1233
1234 catch err
1235 disp(err.message)
1236 stest = false;
1237 end
1238 % </SyntaxCode>
1239
1240 % <AlgoDescription>
1241 %
1242 % 1) Check that calculated psds are identical
1243 %
1244 % </AlgoDescription>
1245
1246 % <AlgoCode>
1247 atest = true;
1248
1249 if stest
1250 % Compare the psd evaluated with the window name / object
1251 if ~eq(S, S1, ple1)
1252 atest = false;
1253 end
1254 else
1255 atest = false;
1256 end
1257 % </AlgoCode>
1258
1259 % Return a result structure
1260 result = utp_prepare_result(atest, stest, dbstack, mfilename);
1261 end % END UTP_46
1262
1263 %% UTP_47
1264
1265 % <TestDescription>
1266 %
1267 % Tests the effect of windowing:
1268 % 1) white noise produced from normal pdf, with:
1269 % a given mean value and sigma (distribution's 1st and 2nd order)
1270 % 2) lpsd passing the window name (Kaiser, psll = 150)
1271 % 3) lpsd passing the window object (Kaiser type, psll = 150)
1272 % 4) compares the two psds
1273 %
1274
1275 % </TestDescription>
1276 function result = utp_47
1277
1278 % <SyntaxDescription>
1279 %
1280 % 1) Prepare the test tsdata:
1281 % white noise from normal distribution + offset
1282 % 2) Calculate the statistical parameters
1283 % 3) Estimate the psd without detrending, Kaiser window (name)
1284 % 4) Estimate the psd without detrending, Kaiser window (object)
1285 % 5) Compare results
1286 %
1287 % </SyntaxDescription>
1288
1289 % <SyntaxCode>
1290 try
1291
1292 noise_type = 'Normal';
1293 scale = 'PSD';
1294 win_type = 'Kaiser';
1295 psll = utils.math.randelement([10:10:200],1);
1296
1297 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, ...
1298 plist('win_test', true));
1299
1300 stest = true;
1301
1302 catch err
1303 disp(err.message)
1304 stest = false;
1305 end
1306 % </SyntaxCode>
1307
1308 % <AlgoDescription>
1309 %
1310 % 1) Check that calculated psds are identical
1311 %
1312 % </AlgoDescription>
1313
1314 % <AlgoCode>
1315 atest = true;
1316
1317 if stest
1318 % Compare the psd evaluated with the window name / object
1319 if ~eq(S, S1, ple1)
1320 atest = false;
1321 end
1322 else
1323 atest = false;
1324 end
1325 % </AlgoCode>
1326
1327 % Return a result structure
1328 result = utp_prepare_result(atest, stest, dbstack, mfilename);
1329 end % END UTP_47
1330
1331 %% UTP_48
1332
1333 % <TestDescription>
1334 %
1335 % Tests the effect of windowing:
1336 % 1) white noise produced from normal pdf, with:
1337 % a given mean value and sigma (distribution's 1st and 2nd order)
1338 % 2) lpsd passing the window name (Kaiser, psll = default)
1339 % 3) lpsd passing the window object (Kaiser type, psll = default)
1340 % 4) compares the two psds
1341 %
1342
1343 % </TestDescription>
1344 function result = utp_48
1345
1346 % <SyntaxDescription>
1347 %
1348 % 1) Prepare the test tsdata:
1349 % white noise from normal distribution + offset
1350 % 2) Calculate the statistical parameters
1351 % 3) Estimate the psd without detrending, Kaiser window (name)
1352 % 4) Estimate the psd without detrending, Kaiser window (object)
1353 % 5) Compare results
1354 %
1355 % </SyntaxDescription>
1356
1357 % <SyntaxCode>
1358 try
1359
1360 noise_type = 'Normal';
1361 scale = 'PSD';
1362 win_type = 'Kaiser';
1363
1364 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, ...
1365 plist('win_test', true));
1366
1367 stest = true;
1368
1369 catch err
1370 disp(err.message)
1371 stest = false;
1372 end
1373 % </SyntaxCode>
1374
1375 % <AlgoDescription>
1376 %
1377 % 1) Check that calculated psds are identical
1378 %
1379 % </AlgoDescription>
1380
1381 % <AlgoCode>
1382 atest = true;
1383
1384 if stest
1385 % Compare the psd evaluated with the window name / object
1386 if ~eq(S, S1, ple1)
1387 atest = false;
1388 end
1389 else
1390 atest = false;
1391 end
1392 % </AlgoCode>
1393
1394 % Return a result structure
1395 result = utp_prepare_result(atest, stest, dbstack, mfilename);
1396 end % END UTP_48
1397
1398
1399 %% UTP_49
1400
1401 % <TestDescription>
1402 %
1403 % Tests the effect of windowing:
1404 % 1) white noise produced from normal pdf, with:
1405 % a given mean value and sigma (distribution's 1st and 2nd order)
1406 % 2) lpsd passing the window name (Nuttall4)
1407 % 3) lpsd passing the window object (Nuttall4 type)
1408 % 4) compares the two psds
1409 %
1410
1411 % </TestDescription>
1412 function result = utp_49
1413
1414 % <SyntaxDescription>
1415 %
1416 % 1) Prepare the test tsdata:
1417 % white noise from normal distribution + offset
1418 % 2) Calculate the statistical parameters
1419 % 3) Estimate the psd without detrending, Nuttall4 window (name)
1420 % 4) Estimate the psd without detrending, Nuttall4 window (object)
1421 % 5) Compare results
1422 %
1423 % </SyntaxDescription>
1424
1425 % <SyntaxCode>
1426 try
1427
1428 noise_type = 'Normal';
1429 scale = 'PSD';
1430 win_type = 'Nuttall4';
1431
1432 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, ...
1433 plist('win_test', true));
1434
1435 stest = true;
1436
1437 catch err
1438 disp(err.message)
1439 stest = false;
1440 end
1441 % </SyntaxCode>
1442
1443 % <AlgoDescription>
1444 %
1445 % 1) Check that calculated psds are identical
1446 %
1447 % </AlgoDescription>
1448
1449 % <AlgoCode>
1450 atest = true;
1451
1452 if stest
1453 % Compare the psd evaluated with the window name / object
1454 if ~eq(S, S1, ple1)
1455 atest = false;
1456 end
1457 else
1458 atest = false;
1459 end
1460 % </AlgoCode>
1461
1462 % Return a result structure
1463 result = utp_prepare_result(atest, stest, dbstack, mfilename);
1464 end % END UTP_49
1465
1466 %% UTP_50
1467
1468 % <TestDescription>
1469 %
1470 % Tests the effect of windowing:
1471 % 1) white noise produced from normal pdf, with:
1472 % a given mean value and sigma (distribution's 1st and 2nd order)
1473 % 2) lpsd passing the window name (SFT3F)
1474 % 3) lpsd passing the window object (SFT3F type)
1475 % 4) compares the two psds
1476 %
1477
1478 % </TestDescription>
1479 function result = utp_50
1480
1481 % <SyntaxDescription>
1482 %
1483 % 1) Prepare the test tsdata:
1484 % white noise from normal distribution + offset
1485 % 2) Calculate the statistical parameters
1486 % 3) Estimate the psd without detrending, SFT3F window (name)
1487 % 4) Estimate the psd without detrending, SFT3F window (object)
1488 % 5) Compare results
1489 %
1490 % </SyntaxDescription>
1491
1492 % <SyntaxCode>
1493 try
1494
1495 noise_type = 'Normal';
1496 scale = 'PSD';
1497 win_type = 'SFT3F';
1498
1499 [a, S, S1] = prepare_analyze_noise(win_type, noise_type, scale, ...
1500 plist('win_test', true));
1501
1502 stest = true;
1503
1504 catch err
1505 disp(err.message)
1506 stest = false;
1507 end
1508 % </SyntaxCode>
1509
1510 % <AlgoDescription>
1511 %
1512 % 1) Check that calculated psds are identical
1513 %
1514 % </AlgoDescription>
1515
1516 % <AlgoCode>
1517 atest = true;
1518
1519 if stest
1520 % Compare the psd evaluated with the window name / object
1521 if ~eq(S, S1, ple1)
1522 atest = false;
1523 end
1524 else
1525 atest = false;
1526 end
1527 % </AlgoCode>
1528
1529 % Return a result structure
1530 result = utp_prepare_result(atest, stest, dbstack, mfilename);
1531 end % END UTP_50
1532
1533
1534 %% Helper function for window call construction
1535
1536 function [a, S_name, S_obj] = prepare_analyze_noise(win_type, noise_type, scale, pli)
1537 % Array of parameters to pick from
1538 fs_list = [0.1;1;2;5;10];
1539 nsecs_list = [20 100 1000:1000:10000]';
1540 sigma_distr_list = [1e-6 2e-3 0.25 1:0.1:10]';
1541 trend_0_list = [1e-6 2e-3 0.25 1:0.1:10]';
1542
1543
1544 % Build time-series test data
1545
1546 % Picks the values at random from the list
1547 fs = utils.math.randelement(fs_list, 1);
1548 nsecs = utils.math.randelement(nsecs_list, 1);
1549 sigma_distr = utils.math.randelement(sigma_distr_list, 1);
1550 trend_0 = utils.math.randelement(trend_0_list, 1);
1551
1552 % Pick units and prefix from those supported
1553 unit_list = unit.supportedUnits;
1554 % remove the first empty unit '' from the list, because then is it
1555 % possible that we add a prefix to an empty unit
1556 unit_list = unit_list(2:end);
1557 prefix_list = unit.supportedPrefixes;
1558
1559 % White noise
1560 a_n = ao(plist('waveform', 'noise', ...
1561 'type', noise_type, 'fs', fs, 'nsecs', nsecs, 'sigma', sigma_distr));
1562
1563 % Constant signal
1564 a_c = ao(trend_0);
1565
1566 % Total signal
1567 a = a_n + a_c;
1568
1569 % Set units
1570 a.setYunits(unit([cell2mat(utils.math.randelement(prefix_list,1)) cell2mat(utils.math.randelement(unit_list,1))]));
1571
1572 % Evaluate the lpsd of the white noise time-series data
1573 olap = find(pli, 'olap');
1574 if isempty(olap)
1575 olap = 0; %% Should we take the rov instead?
1576 end
1577 detrend_order = 0;
1578
1579 switch lower(win_type)
1580 case 'kaiser'
1581 psll = find(pli, 'psll');
1582 if isempty(psll)
1583 psll = find(ao.getInfo('lpsd').plists, 'psll');
1584 end
1585 pl_spec = plist('Win', win_type, 'psll', psll, 'olap', olap, ...
1586 'order', detrend_order, 'scale', scale);
1587 otherwise
1588 pl_spec = plist('Win', win_type, 'olap', olap, ...
1589 'order', detrend_order, 'scale', scale);
1590 end
1591
1592 if find(pli, 'win_test')
1593 % Calls the lpsd applying the detrend and window internally
1594 % (passig window object)
1595 S_obj = a.lpsd(pl_spec);
1596 else
1597 S_obj = ao;
1598 end
1599
1600 % Calls the lpsd applying the detrend and window internally
1601 % (passig window name)
1602 pl_spec.pset('Win', win_type);
1603 S_name = a.lpsd(pl_spec);
1604
1605 end
1606
1607
1608 end