source: sasview/park-1.2.1/park/fitresult.py @ 3570545

ESS_GUIESS_GUI_DocsESS_GUI_batch_fittingESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_iss879ESS_GUI_iss959ESS_GUI_openclESS_GUI_orderingESS_GUI_sync_sascalccostrafo411magnetic_scattrelease-4.1.1release-4.1.2release-4.2.2release_4.0.1ticket-1009ticket-1094-headlessticket-1242-2d-resolutionticket-1243ticket-1249ticket885unittest-saveload
Last change on this file since 3570545 was 3570545, checked in by Mathieu Doucet <doucetm@…>, 12 years ago

Adding park Part 2

  • Property mode set to 100644
File size: 7.5 KB
Line 
1import sys, math, time
2import numpy
3
4from formatnum import format_uncertainty
5
6class FitHandler(object):
7    """
8    Abstract interface for fit thread handler.
9
10    The methods in this class are called by the optimizer as the fit
11    progresses.
12
13    Note that it is up to the optimizer to call the fit handler correctly,
14    reporting all status changes and maintaining the 'done' flag.
15    """
16    done = False
17    """True when the fit job is complete"""
18    result = None
19    """The current best result of the fit"""
20
21    def improvement(self):
22        """
23        Called when a result is observed which is better than previous
24        results from the fit.
25
26        result is a FitResult object, with parameters, #calls and fitness.
27        """
28    def error(self, msg):
29        """
30        Model had an error; print traceback
31        """
32    def progress(self, current, expected):
33        """
34        Called each cycle of the fit, reporting the current and the
35        expected amount of work.   The meaning of these values is
36        optimizer dependent, but they can be converted into a percent
37        complete using (100*current)//expected.
38
39        Progress is updated each iteration of the fit, whatever that
40        means for the particular optimization algorithm.  It is called
41        after any calls to improvement for the iteration so that the
42        update handler can control I/O bandwidth by suppressing
43        intermediate improvements until the fit is complete.
44        """
45    def finalize(self):
46        """
47        Fit is complete; best results are reported
48        """
49    def abort(self):
50        """
51        Fit was aborted.
52        """
53
54class ConsoleUpdate(FitHandler):
55    """
56    Print progress to the console.
57    """
58    isbetter = False
59    """Record whether results improved since last update"""
60    progress_delta =  60
61    """Number of seconds between progress updates"""
62    improvement_delta = 5
63    """Number of seconds between improvement updates"""
64    def __init__(self,quiet=False,progress_delta=60,improvement_delta=5):
65        """
66        If quiet is true, only print out final summary, not progress and
67        improvements.
68        """
69        self.progress_time = time.time()
70        self.progress_percent = 0
71        self.improvement_time = self.progress_time
72        self.isbetter = False
73        self.quiet = quiet
74        self.progress_delta = progress_delta
75        self.improvement_delta = improvement_delta
76
77    def progress(self, k, n):
78        """
79        Report on progress.
80        """
81        if self.quiet: return
82        t = time.time()
83        p = int((100*k)//n)
84
85        # Show improvements if there are any
86        dt = t - self.improvement_time
87        if self.isbetter and dt > self.improvement_delta:
88            self.result.print_summary()
89            self.isbetter = False
90            self.improvement_time = t
91
92        # Update percent complete
93        dp = p-self.progress_percent
94        if dp < 1: return
95        dt = t - self.progress_time
96        if dt > self.progress_delta:
97            if 1 <= dp <= 2:
98                print "%d%% complete"%p
99                self.progress_percent = p
100                self.progress_time = t
101            elif 2 < dp <= 5:
102                if p//5 != self.progress_percent//5:
103                    print "%d%% complete"%(5*(p//5))
104                    self.progress_percent = p
105                    self.progress_time = t
106            else:
107                if p//10 != self.progress_percent//10:
108                    print "%d%% complete"%(10*(p//10))
109                    self.progress_percent = p
110                    self.progress_time = t
111
112    def improvement(self):
113        """
114        Called when a result is observed which is better than previous
115        results from the fit.
116        """
117        self.isbetter = True
118
119    def error(self, msg):
120        """
121        Model had an error; print traceback
122        """
123        if self.isbetter:
124            self.result.print_summary()
125        print msg
126
127    def finalize(self):
128        if self.isbetter:
129            self.result.print_summary()
130        print "Total function calls:",self.result.calls
131
132    def abort(self):
133        if self.isbetter:
134            self.result.print_summary()
135
136
137class FitParameter(object):
138    """
139    Fit result for an individual parameter.
140    """
141    def __init__(self, name, range, value):
142        self.name = name
143        self.range = range
144        self.value = value
145        self.stderr = None
146    def summarize(self):
147        """
148        Return parameter range string.
149
150        E.g.,  "       Gold .....|.... 5.2043 in [2,7]"
151        """
152        bar = ['.']*10
153        lo,hi = self.range
154        if numpy.isfinite(lo)and numpy.isfinite(hi):
155            portion = (self.value-lo)/(hi-lo)
156            if portion < 0: portion = 0.
157            elif portion >= 1: portion = 0.99999999
158            barpos = int(math.floor(portion*len(bar)))
159            bar[barpos] = '|'
160        bar = "".join(bar)
161        lostr = "[%g"%lo if numpy.isfinite(lo) else "(-inf"
162        histr = "%g]"%hi if numpy.isfinite(hi) else "inf)"
163        valstr = format_uncertainty(self.value, self.stderr)
164        return "%25s %s %s in %s,%s"  % (self.name,bar,valstr,lostr,histr)
165    def __repr__(self):
166        return "FitParameter('%s')"%self.name
167
168class FitResult(object):
169    """
170    Container for reporting fit results.
171    """
172    def __init__(self, parameters, calls, fitness):
173        self.parameters = parameters
174        """Fit parameter list, each with name, range and value attributes."""
175        self.calls = calls
176        """Number of function calls"""
177        self.fitness = fitness
178        """Value of the goodness of fit metric"""
179        self.pvec = numpy.array([p.value for p in self.parameters])
180        """Parameter vector"""
181        self.stderr = None
182        """Parameter uncertainties"""
183        self.cov = None
184        """Covariance matrix"""
185
186    def update(self, pvec, fitness, calls):
187        self.calls = calls
188        self.fitness = fitness
189        self.pvec = pvec.copy()
190        for i,p in enumerate(self.parameters):
191            p.value = pvec[i]
192
193    def calc_cov(self, fn):
194        """
195        Return the covariance matrix inv(J'J) at point p.
196        """
197        if hasattr(fn, 'jacobian'):
198            # Find cov of f at p
199            #     cov(f,p) = inv(J'J)
200            # Use SVD
201            #     J = U S V'
202            #     J'J = (U S V')' (U S V')
203            #         = V S' U' U S V'
204            #         = V S S V'
205            #     inv(J'J) = inv(V S S V')
206            #              = inv(V') inv(S S) inv(V)
207            #              = V inv (S S) V'
208            J = fn.jacobian(self.pvec)
209            u,s,vh = numpy.linalg.svd(J,0)
210            JTJinv = numpy.dot(vh.T.conj()/s**2,vh)
211            self.set_cov(JTJinv)
212
213    def set_cov(self, cov):
214        """
215        Return the covariance matrix inv(J'J) at point p.
216        """
217        self.cov = cov
218        if cov is not None:
219            self.stderr = numpy.sqrt(numpy.diag(self.cov))
220            # Set the uncertainties on the individual parameters
221            for k,p in enumerate(self.parameters):
222                p.stderr = self.stderr[k]
223        else:
224            self.stderr = None
225            # Reset the uncertainties on the individual parameters
226            for k,p in enumerate(self.parameters):
227                p.stderr = None
228
229    def __str__(self):
230        if self.parameters == None: return "No results"
231        L = ["P%-3d %s"%(n+1,p.summarize()) for n,p in enumerate(self.parameters)]
232        L.append("=== goodness of fit: %g"%(self.fitness))
233        return "\n".join(L)
234
235    def print_summary(self, fid=sys.stdout):
236        print >>fid, self
237
Note: See TracBrowser for help on using the repository browser.