0
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
1 mc
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
2
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
3 a = 2;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
4 b = 3;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
5 x = -1:0.01:1;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
6 xa = ao(plist('dtype', 'cdata', 'yvals', x));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
7
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
8 n = 1*randn(size(x));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
9 sn = 2;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
10 na = sn.*ao(plist('dtype', 'cdata', 'yvals', n));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
11
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
12 y0 = a.*(xa + b).^3;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
13 ya = y0 + na;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
14
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
15 %% linearisation
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
16
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
17
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
18 % ylin = dyda.*da + dydb.*db + y0;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
19
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
20 da = 0.2;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
21 db = -0.3;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
22 ag = a+da;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
23 bg = b+db;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
24 y0 = ag.*(xa + bg).^3;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
25
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
26 dyda = (xa+bg).^3;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
27 dydb = ag.*(xa+bg).^2;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
28
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
29 tgtData = ya-y0;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
30 fitData = lscov(dyda, dydb, tgtData);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
31
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
32 dfit = fitData - tgtData;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
33 iplot(dfit)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
34 iplot(hist(dfit))
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
35
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
36 params = find(fitData.procinfo, 'P')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
37 af = ag + params(1)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
38 bf = bg + params(2)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
39
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
40 errs = find(fitData.procinfo, 'STD');
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
41
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
42 afe = 100*errs(1)./af
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
43 bfe = 100*errs(2)./bf
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
44
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
45 %% Non-linear fit
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
46
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
47 xy = ao(plist('xvals', x, 'yvals', ya.y));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
48 pl = plist('Function', 'P(1).*(Xdata + P(2)).^3', 'P0', [ag bg]);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
49 nlfit = curvefit(xy, pl);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
50
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
51
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
52
|