+
−%  test_ao_ltfe_variance_montecarlo
+
−%
+
−% Tests that the standard deviation returned by ao.dy in one
+
−% frequency bin is equivalent to the matlab's std taking
+
−% considering all realisations
+
−%
+
−% M Nofrarias 22-07-09
+
−%
+
−% $Id: test_ao_lcpsd_variance_montecarlo.m,v 1.1 2009/08/11 14:20:10 miquel Exp $
+
−
+
−% function test_ao_ltfe_variance_montecarlo()
+
−
+
−clear
+
−
+
−% data
+
−nsecs = 200;
+
−fs    = 10;
+
−pl = plist('nsecs', nsecs, 'fs', fs, 'tsfcn', 'sin(2*pi*2*t) + randn(size(t))');
+
−
+
−
+
−% Make a filter
+
−f1 = miir(plist('type', 'highpass', 'fc', 4, 'fs', fs));
+
−
+
−% Window
+
−Nfft = 1000;
+
−win  = specwin('Hanning', Nfft);
+
−pl2 = plist('Lmin',200,'win',win,'Olap',-1)
+
−
+
−% loop
+
−for  i = 1:100
+
−  a1 = ao(pl);
+
−  a2 = filter(a1,plist('filter', f1));
+
−  t1(i) = lcpsd(a1,a2,pl2);
+
−end
+
−
+
−%% mean
+
−index = 6;
+
−
+
−clear rel
+
−for i =1:len(t1(1))
+
−  mn(i) = [mean(t1(:).y(i))]; % both means are equal
+
−  rel(:,i) = [std(t1(:).y(i)) mean(t1(:).dy(i))]/abs(mn(i));
+
−end
+
−
+
−k = find(t1(1).procinfo,'k')
+
−
+
−
+
−figure
+
−loglog(t1(1).x,rel')
+
−figure
+
−loglog(t1(1).x,abs(rel(2,:)-rel(1,:)))
+
−