1 | from sans.pr.core.pr_inversion import Cinvertor |
---|
2 | import numpy |
---|
3 | |
---|
4 | class Invertor(Cinvertor): |
---|
5 | |
---|
6 | ## Chisqr of the last computation |
---|
7 | chisqr = 0 |
---|
8 | |
---|
9 | def __init__(self): |
---|
10 | Cinvertor.__init__(self) |
---|
11 | |
---|
12 | def __setattr__(self, name, value): |
---|
13 | """ |
---|
14 | Set the value of an attribute. |
---|
15 | Access the parent class methods for |
---|
16 | x, y, err and d_max. |
---|
17 | """ |
---|
18 | if name=='x': |
---|
19 | if 0.0 in value: |
---|
20 | raise ValueError, "Invertor: one of your q-values is zero. Delete that entry before proceeding" |
---|
21 | return self.set_x(value) |
---|
22 | elif name=='y': |
---|
23 | return self.set_y(value) |
---|
24 | elif name=='err': |
---|
25 | return self.set_err(value) |
---|
26 | elif name=='d_max': |
---|
27 | return self.set_dmax(value) |
---|
28 | elif name=='alpha': |
---|
29 | return self.set_alpha(value) |
---|
30 | |
---|
31 | return Cinvertor.__setattr__(self, name, value) |
---|
32 | |
---|
33 | def __getattr__(self, name): |
---|
34 | """ |
---|
35 | Return the value of an attribute |
---|
36 | For the moment x, y, err and d_max are write-only |
---|
37 | TODO: change that! |
---|
38 | """ |
---|
39 | import numpy |
---|
40 | if name=='x': |
---|
41 | out = numpy.ones(self.get_nx()) |
---|
42 | self.get_x(out) |
---|
43 | return out |
---|
44 | elif name=='y': |
---|
45 | out = numpy.ones(self.get_ny()) |
---|
46 | self.get_y(out) |
---|
47 | return out |
---|
48 | elif name=='err': |
---|
49 | out = numpy.ones(self.get_nerr()) |
---|
50 | self.get_err(out) |
---|
51 | return out |
---|
52 | elif name=='d_max': |
---|
53 | return self.get_dmax() |
---|
54 | elif name=='alpha': |
---|
55 | return self.get_alpha() |
---|
56 | elif name in self.__dict__: |
---|
57 | return self.__dict__[name] |
---|
58 | return None |
---|
59 | |
---|
60 | def invert(self, nfunc=5): |
---|
61 | """ |
---|
62 | Perform inversion to P(r) |
---|
63 | """ |
---|
64 | from scipy import optimize |
---|
65 | |
---|
66 | # First, check that the current data is valid |
---|
67 | if self.is_valid()<=0: |
---|
68 | raise RuntimeError, "Invertor.invert: Data array are of different length" |
---|
69 | |
---|
70 | p = numpy.ones(nfunc) |
---|
71 | out, cov_x, info, mesg, success = optimize.leastsq(self.residuals, p, full_output=1, warning=True) |
---|
72 | |
---|
73 | # Compute chi^2 |
---|
74 | res = self.residuals(out) |
---|
75 | chisqr = 0 |
---|
76 | for i in range(len(res)): |
---|
77 | chisqr += res[i] |
---|
78 | |
---|
79 | self.chi2 = chisqr |
---|
80 | |
---|
81 | return out, cov_x |
---|
82 | |
---|
83 | def pr_fit(self, nfunc=5): |
---|
84 | """ |
---|
85 | Perform inversion to P(r) |
---|
86 | """ |
---|
87 | from scipy import optimize |
---|
88 | |
---|
89 | # First, check that the current data is valid |
---|
90 | if self.is_valid()<=0: |
---|
91 | raise RuntimeError, "Invertor.invert: Data arrays are of different length" |
---|
92 | |
---|
93 | p = numpy.ones(nfunc) |
---|
94 | out, cov_x, info, mesg, success = optimize.leastsq(self.pr_residuals, p, full_output=1, warning=True) |
---|
95 | |
---|
96 | # Compute chi^2 |
---|
97 | res = self.pr_residuals(out) |
---|
98 | chisqr = 0 |
---|
99 | for i in range(len(res)): |
---|
100 | chisqr += res[i] |
---|
101 | |
---|
102 | self.chisqr = chisqr |
---|
103 | |
---|
104 | return out, cov_x |
---|
105 | |
---|
106 | def pr_err(self, c, c_cov, r): |
---|
107 | import math |
---|
108 | c_err = numpy.zeros(len(c)) |
---|
109 | for i in range(len(c)): |
---|
110 | try: |
---|
111 | c_err[i] = math.sqrt(math.fabs(c_cov[i][i])) |
---|
112 | except: |
---|
113 | import sys |
---|
114 | print sys.exc_value |
---|
115 | print "oups", c_cov[i][i] |
---|
116 | c_err[i] = c[i] |
---|
117 | |
---|
118 | return self.get_pr_err(c, c_err, r) |
---|
119 | |
---|
120 | |
---|
121 | if __name__ == "__main__": |
---|
122 | o = Invertor() |
---|
123 | |
---|
124 | |
---|
125 | |
---|
126 | |
---|
127 | |
---|