0
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
1 % TEST_AO_LINLSQSVD tests the linlsqsvd method of the AO class.
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
3 % L Ferraioli 10-11-2010
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
4 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
5 % $Id: test_ao_linlsqsvd.m,v 1.2 2011/02/18 16:15:28 luigi Exp $
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
6 %
|
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 %% 1) Determine the coefficients of a linear combination of noises and
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
10 %% comapre with lscov:
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
11 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
12 % Make some data
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
13 fs = 10;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
14 nsecs = 10;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
15 B1 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'T'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
16 B1.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
17 B2 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'T'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
18 B2.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
19 B3 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'T'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
20 B3.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
21 n = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'm'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
22 c = [ao(1,plist('yunits','m/T')) ao(2,plist('yunits','m/T')) ao(3,plist('yunits','m T^-1'))];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
23 y = c(1)*B1 + c(2)*B2 + c(3)*B3 + n;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
24 y.simplifyYunits;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
25 %%% Get a fit with linlsqsvd
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
26 pobj1 = linlsqsvd(B1, B2, B3, y);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
27 % do linear combination: using lincom
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
28 yfit1 = lincom(B1, B2, B3, pobj1);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
29 yfit1.simplifyYunits;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
30 %%% do fit using lscov
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
31 pobj2 = lscov(B1, B2, B3, y);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
32 % do linear combination: using lincom
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
33 yfit2 = lincom(B1, B2, B3, pobj2);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
34 yfit2.simplifyYunits;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
35
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
36 %%% do linear combination: using eval
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
37 yfit3 = pobj1.eval(B1, B2, B3);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
38
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
39
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
40 % Plot (compare data with fit)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
41 iplot(y, yfit1, yfit2, yfit3)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
42
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
43 %% 2) Determine the coefficients of a linear combination of noises:
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
44 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
45 % Make some data
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
46 fs = 10;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
47 nsecs = 10;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
48 x1 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'T'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
49 x1.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
50 x2 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'm'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
51 x2.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
52 x3 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'C'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
53 x3.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
54 n = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'm'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
55 c = [ao(1,plist('yunits','m/T')) ao(2,plist('yunits','m/m')) ao(3,plist('yunits','m C^-1'))];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
56 y = c(1)*x1 + c(2)*x2 + c(3)*x3 + n;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
57 y.simplifyYunits;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
58
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
59 %%% Get a fit for c
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
60 pobj = linlsqsvd(x1, x2, x3, y);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
61 % do linear combination: using lincom
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
62 yfit1 = lincom(x1, x2, x3, pobj);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
63
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
64 %%% do linear combination: using eval
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
65 pl_split = plist('times', [1 5]);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
66 yfit2 = pobj.eval(plist('Xdata', {split(x1, pl_split), split(x2, pl_split), split(x3, pl_split)}));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
67
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
68 % Plot (compare data with fit)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
69 iplot(y, yfit1, yfit2, plist('Linestyles', {'-','--'}))
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
70
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
71
|