0
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
1 % Test file for ltpda_timedomain fit. This does a noise subtraction from
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
2 % engineering model data.
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
3 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
4 % M Hewitson 18-03-07
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
5 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
6 % $Id: test_timeDomain.m,v 1.3 2007/12/07 08:54:07 hewitson Exp $
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
7 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
8
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
9 clear all;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
10
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
11 %% Load data into analysis objects
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
12 disp('*** Loading data....')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
13 x12 = ao('x12.txt'); %x12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
14 phi1 = ao('phi1.txt'); %phi1
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
15 eta1 = ao('eta1.txt'); %eta1
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
16 phi12 = ao('phi12.txt'); %phi12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
17 eta12 = ao('eta12.txt'); %eta12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
18 disp('*** Done.')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
19
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
20 % get the sample frequency for designing bandpass filter
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
21 d1 = x12.data;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
22 fs = d1.fs;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
23
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
24 %% Load bandpass
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
25 pl = plist();
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
26 pl = append(pl, param('type', 'bandpass'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
27 % pl = append(pl, param('fs', fs));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
28 pl = append(pl, param('fc', [0.003 0.03]));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
29 pl = append(pl, param('order', 2));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
30 bp = miir(pl);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
31
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
32 %% Filter data and remove startup transients
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
33
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
34 % parameter list for split()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
35 % to remove startup transient of filters
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
36 spl = plist(param('times', [1000 2000]));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
37
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
38 % parameter list for filter()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
39 fpl = plist(param('filter', bp));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
40
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
41 % filter and split data
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
42 x12f = split(filter(x12, fpl), spl);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
43 phi1f = split(filter(phi1, fpl), spl);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
44 eta1f = split(filter(eta1, fpl), spl);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
45 phi12f = split(filter(phi12, fpl), spl);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
46 eta12f = split(filter(eta12, fpl), spl);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
47
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
48
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
49 %% Make Timedomainfit
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
50
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
51 coeffsAO = ltpda_timedomainfit(x12f, phi1f, eta1f, phi12f, eta12f);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
52
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
53 %% Make linear combination
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
54
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
55 x12ns = ltpda_lincom(x12, phi1, eta1, phi12, eta12, coeffsAO);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
56
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
57 %% Make LPSD of each
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
58
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
59 % Window function
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
60 w = specwin('Kaiser', 10, 150);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
61
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
62 % parameter list for lpsd
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
63 pl = plist();
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
64 pl = append(pl, param('Kdes', 100));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
65 pl = append(pl, param('Kmin', 2));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
66 pl = append(pl, param('Jdes', 500));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
67 pl = append(pl, param('Win', w));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
68
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
69 % LPSD of filtered x12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
70 x12fxx = ltpda_lpsd(x12f, pl);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
71 % LPSD of x12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
72 x12dt = ltpda_polydetrend(x12, plist(param('N', 1)));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
73 x12xx = ltpda_lpsd(x12dt, pl);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
74 % LPSD of noise-subtracted x12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
75 x12nsdt = ltpda_polydetrend(x12ns, plist(param('N', 3)));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
76 x12nsxx = ltpda_lpsd(x12ns, pl);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
77
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
78 %% Plot results
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
79
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
80 figure
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
81 plot([x12fxx x12xx x12nsxx])
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
82 ltpda_allxaxis(0.001, 1);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
83 subplot(3,1,1:2)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
84 legend('Bandpassed X12', 'X12', 'X12 with noise subtracted');
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
85 %%
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
86 figure
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
87 plot(x12ns.hist)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
88
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
89 %%
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
90
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
91 save(x12fxx, 'x12fxx.xml');
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
92 save(x12xx, 'x12xx.xml');
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
93 save(x12nsxx, 'x12nsxx.xml');
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
94
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
95 return
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
96
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
97 %% Reproduce from history
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
98
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
99 % Write an m-file from AO
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
100 ao2m(x12nsxx, 'test.m');
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
101
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
102 %% now run it
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
103 clear;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
104 test
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
105
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
106 figure
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
107 plot(a1)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
108 figure
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
109 plot(a1.hist)
|