diff m-toolbox/test/test_ao_firwhiten.m @ 0:f0afece42f48

Import.
author Daniele Nicolodi <nicolodi@science.unitn.it>
date Wed, 23 Nov 2011 19:22:13 +0100
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/m-toolbox/test/test_ao_firwhiten.m	Wed Nov 23 19:22:13 2011 +0100
@@ -0,0 +1,93 @@
+% function test_ao_firwhiten
+% A test script for ao/firwhiten
+% 
+% M Hewitson 21-04-08
+% 
+% $Id: test_ao_firwhiten.m,v 1.5 2009/08/26 18:27:06 luigi Exp $
+% 
+
+mc
+
+%% Make an AO
+
+fs    = 10;
+nsecs = 1e3;
+
+% pole zero model
+pzm = pzmodel(2, [pz(0.1, 2) pz(0.5, 1)] , [pz(0.8) pz(1, 2)]);
+
+% parameter list for ltpda_noisegen
+pl = plist();
+pl = append(pl, param('nsecs', nsecs));
+pl = append(pl, param('fs', fs));
+
+% calling the noisegenerator
+a = fngen(pzm, pl);
+% a = set(a, 'yunits', 'V');
+
+% Make series of sine waves
+
+% s = ao(plist('waveform', 'sine wave', ...
+%              'A', [0.1 1 0.1], ...
+%              'f', [3.3  0.42 1.23], ...
+%              'phi', [0 0 0], ...
+%              'fs', fs, 'nsecs', nsecs));
+
+s = ao(plist('tsfcn', '0.05*sin(2*pi*1.23*t)', ...
+             'fs', fs, 'nsecs', nsecs));
+% s = ao(plist('waveform', 'sine wave', ...
+%              'A', [0.05], ...
+%              'f', [1.23], ...
+%              'phi', [0], ...
+%              'fs', fs, 'nsecs', nsecs));
+
+% Add signal and noise
+a = a + s;
+
+%% Whiten
+Nfft = fs*100;
+[aw, ffs, nfs] = firwhiten(a, plist('Ntaps', 1024, 'Nfft', Nfft));
+
+% plot filter response
+iplot(resp(ffs))
+
+
+%% Remove end effects
+
+spl = plist('split_type', 'times', 'times', [10 -20]);
+a   = timeshift(split(a, spl));
+aw  = timeshift(split(aw, spl));
+
+%% Compare time-series
+iplot(a,aw, plist('Arrangement', 'subplots', 'XRanges', {'All', [0 10]}))
+
+%% Compare spectra
+
+pxx = pwelch(a, aw, plist('Nfft', fs*100));
+iplot(pxx)
+
+%% Help test
+
+fs = 1;
+a = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', 1e5, 'yunits','m'));
+
+% filter
+pl = plist('type','bandpass', ...
+  'order',3, ...
+  'gain',1, ...
+  'fc',[0.03 0.1], ...
+  'fs',fs);
+ft = miir(pl);
+
+% coloring white noise
+af = filter(a, ft);
+
+% firwhiten
+[aw, ffs, nfs] = firwhiten(af, plist('Ntaps', 5000, 'Nfft', 1e5, 'BW', 5));
+
+awxx = aw.psd;
+afxx = af.psd;
+
+iplot(afxx,awxx)
+
+% END
\ No newline at end of file