source: sasview/park_integration/AbstractFitEngine.py @ c6036f5

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 c6036f5 was 83195f7, checked in by Jae Cho <jhjcho@…>, 15 years ago

updated according to the changes in 2d inputs

  • Property mode set to 100644
File size: 22.6 KB
Line 
1import logging, sys
2import park,numpy,math, copy
3from DataLoader.data_info import Data1D
4from DataLoader.data_info import Data2D
5class SansParameter(park.Parameter):
6    """
7        SANS model parameters for use in the PARK fitting service.
8        The parameter attribute value is redirected to the underlying
9        parameter value in the SANS model.
10    """
11    def __init__(self, name, model):
12        """
13            @param name: the name of the model parameter
14            @param model: the sans model to wrap as a park model
15        """
16        self._model, self._name = model,name
17        #set the value for the parameter of the given name
18        self.set(model.getParam(name))
19         
20    def _getvalue(self):
21        """
22            override the _getvalue of park parameter
23            @return value the parameter associates with self.name
24        """
25        return self._model.getParam(self.name)
26   
27    def _setvalue(self,value):
28        """
29            override the _setvalue pf park parameter
30            @param value: the value to set on a given parameter
31        """
32        self._model.setParam(self.name, value)
33       
34    value = property(_getvalue,_setvalue)
35   
36    def _getrange(self):
37        """
38            Override _getrange of park parameter
39            return the range of parameter
40        """
41        #if not  self.name in self._model.getDispParamList():
42        lo,hi = self._model.details[self.name][1:3]
43        if lo is None: lo = -numpy.inf
44        if hi is None: hi = numpy.inf
45        #else:
46            #lo,hi = self._model.details[self.name][1:]
47            #if lo is None: lo = -numpy.inf
48            #if hi is None: hi = numpy.inf
49        if lo >= hi:
50            raise ValueError,"wrong fit range for parameters"
51       
52        return lo,hi
53   
54    def _setrange(self,r):
55        """
56            override _setrange of park parameter
57            @param r: the value of the range to set
58        """
59        self._model.details[self.name][1:3] = r
60    range = property(_getrange,_setrange)
61   
62class Model(park.Model):
63    """
64        PARK wrapper for SANS models.
65    """
66    def __init__(self, sans_model, **kw):
67        """
68            @param sans_model: the sans model to wrap using park interface
69        """
70        park.Model.__init__(self, **kw)
71        self.model = sans_model
72        self.name = sans_model.name
73        #list of parameters names
74        self.sansp = sans_model.getParamList()
75        #list of park parameter
76        self.parkp = [SansParameter(p,sans_model) for p in self.sansp]
77        #list of parameterset
78        self.parameterset = park.ParameterSet(sans_model.name,pars=self.parkp)
79        self.pars=[]
80 
81 
82    def getParams(self,fitparams):
83        """
84            return a list of value of paramter to fit
85            @param fitparams: list of paramaters name to fit
86        """
87        list=[]
88        self.pars=[]
89        self.pars=fitparams
90        for item in fitparams:
91            for element in self.parkp:
92                 if element.name ==str(item):
93                     list.append(element.value)
94        return list
95   
96   
97    def setParams(self,paramlist, params):
98        """
99            Set value for parameters to fit
100            @param params: list of value for parameters to fit
101        """
102        try:
103            for i in range(len(self.parkp)):
104                for j in range(len(paramlist)):
105                    if self.parkp[i].name==paramlist[j]:
106                        self.parkp[i].value = params[j]
107                        self.model.setParam(self.parkp[i].name,params[j])
108        except:
109            raise
110 
111    def eval(self,x):
112        """
113            override eval method of park model.
114            @param x: the x value used to compute a function
115        """
116        try:
117            return self.model.evalDistribution(x)
118        except:
119            raise
120
121   
122class FitData1D(Data1D):
123    """
124        Wrapper class  for SANS data
125        FitData1D inherits from DataLoader.data_info.Data1D. Implements
126        a way to get residuals from data.
127    """
128    def __init__(self,x, y,dx= None, dy=None, smearer=None):
129        Data1D.__init__(self, x=numpy.array(x), y=numpy.array(y), dx=dx, dy=dy)
130        """
131            @param smearer: is an object of class QSmearer or SlitSmearer
132            that will smear the theory data (slit smearing or resolution
133            smearing) when set.
134           
135            The proper way to set the smearing object would be to
136            do the following:
137           
138            from DataLoader.qsmearing import smear_selection
139            smearer = smear_selection(some_data)
140            fitdata1d = FitData1D( x= [1,3,..,],
141                                    y= [3,4,..,8],
142                                    dx=None,
143                                    dy=[1,2...], smearer= smearer)
144           
145            Note that some_data _HAS_ to be of class DataLoader.data_info.Data1D
146           
147            Setting it back to None will turn smearing off.
148           
149        """
150       
151        self.smearer = smearer
152        if dy ==None or dy==[]:
153            self.dy= numpy.zeros(len(self.y)) 
154        else:
155            self.dy= numpy.asarray(dy)
156     
157        # For fitting purposes, replace zero errors by 1
158        #TODO: check validity for the rare case where only
159        # a few points have zero errors
160        self.dy[self.dy==0]=1
161       
162        ## Min Q-value
163        #Skip the Q=0 point, especially when y(q=0)=None at x[0].
164        if min (self.x) ==0.0 and self.x[0]==0 and not numpy.isfinite(self.y[0]):
165            self.qmin = min(self.x[self.x!=0])
166        else:                             
167            self.qmin= min (self.x)
168        ## Max Q-value
169        self.qmax = max (self.x)
170       
171        # Range used for input to smearing
172        self._qmin_unsmeared = self.qmin
173        self._qmax_unsmeared = self.qmax
174        # Identify the bin range for the unsmeared and smeared spaces
175        self.idx = (self.x>=self.qmin) & (self.x <= self.qmax)
176        self.idx_unsmeared = (self.x>=self._qmin_unsmeared) & (self.x <= self._qmax_unsmeared)
177 
178       
179       
180    def setFitRange(self,qmin=None,qmax=None):
181        """ to set the fit range"""
182        # Skip Q=0 point, (especially for y(q=0)=None at x[0]).
183        #ToDo: Fix this.
184        if qmin==0.0 and not numpy.isfinite(self.y[qmin]):
185            self.qmin = min(self.x[self.x!=0])
186        elif qmin!=None:                       
187            self.qmin = qmin           
188
189        if qmax !=None:
190            self.qmax = qmax
191           
192        # Determine the range needed in unsmeared-Q to cover
193        # the smeared Q range
194        self._qmin_unsmeared = self.qmin
195        self._qmax_unsmeared = self.qmax   
196       
197        self._first_unsmeared_bin = 0
198        self._last_unsmeared_bin  = len(self.x)-1
199       
200        if self.smearer!=None:
201            self._first_unsmeared_bin, self._last_unsmeared_bin = self.smearer.get_bin_range(self.qmin, self.qmax)
202            self._qmin_unsmeared = self.x[self._first_unsmeared_bin]
203            self._qmax_unsmeared = self.x[self._last_unsmeared_bin]
204           
205        # Identify the bin range for the unsmeared and smeared spaces
206        self.idx = (self.x>=self.qmin) & (self.x <= self.qmax)
207        self.idx_unsmeared = (self.x>=self._qmin_unsmeared) & (self.x <= self._qmax_unsmeared)
208 
209       
210    def getFitRange(self):
211        """
212            @return the range of data.x to fit
213        """
214        return self.qmin, self.qmax
215       
216    def residuals(self, fn):
217        """
218            Compute residuals.
219           
220            If self.smearer has been set, use if to smear
221            the data before computing chi squared.
222           
223            @param fn: function that return model value
224            @return residuals
225        """
226        # Compute theory data f(x)
227        fx= numpy.zeros(len(self.x))
228        fx[self.idx_unsmeared] = fn(self.x[self.idx_unsmeared])
229       
230        ## Smear theory data
231        if self.smearer is not None:
232            fx = self.smearer(fx, self._first_unsmeared_bin, self._last_unsmeared_bin)
233       
234        ## Sanity check
235        if numpy.size(self.dy)!= numpy.size(fx):
236            raise RuntimeError, "FitData1D: invalid error array %d <> %d" % (numpy.shape(self.dy),
237                                                                              numpy.size(fx))
238                                                                             
239        return (self.y[self.idx]-fx[self.idx])/self.dy[self.idx]
240     
241 
242       
243    def residuals_deriv(self, model, pars=[]):
244        """
245            @return residuals derivatives .
246            @note: in this case just return empty array
247        """
248        return []
249   
250   
251class FitData2D(Data2D):
252    """ Wrapper class  for SANS data """
253    def __init__(self,sans_data2d ,data=None, err_data=None):
254        Data2D.__init__(self, data= data, err_data= err_data)
255        """
256            Data can be initital with a data (sans plottable)
257            or with vectors.
258        """
259        self.res_err_image=[]
260        self.index_model=[]
261        self.qmin= None
262        self.qmax= None
263        self.set_data(sans_data2d )
264       
265    def set_data(self, sans_data2d, qmin=None, qmax=None ):
266        """
267            Determine the correct qx_data and qy_data within range to fit
268        """
269        self.data     = sans_data2d.data
270        self.err_data = sans_data2d.err_data
271        self.qx_data = sans_data2d.qx_data
272        self.qy_data = sans_data2d.qy_data
273        self.mask       = sans_data2d.mask
274
275        x_max = max(math.fabs(sans_data2d.xmin), math.fabs(sans_data2d.xmax))
276        y_max = max(math.fabs(sans_data2d.ymin), math.fabs(sans_data2d.ymax))
277       
278        ## fitting range
279        if qmin == None:
280            self.qmin = 1e-16
281        if qmax == None:
282            self.qmax = math.sqrt(x_max*x_max +y_max*y_max)
283        ## new error image for fitting purpose
284        if self.err_data== None or self.err_data ==[]:
285            self.res_err_data= numpy.zeros(len(self.data))
286        else:
287            self.res_err_data = copy.deepcopy(self.err_data)
288        self.res_err_data[self.res_err_data==0]=numpy.sqrt(self.data[self.res_err_data==0])
289       
290        self.radius= numpy.sqrt(self.qx_data**2 + self.qy_data**2)
291       
292        # Note: mask = True: for MASK while mask = False for NOT to mask
293        self.index_model = ((self.qmin <= self.radius)&(self.radius<= self.qmax))
294        self.mask = (self.index_model) & (self.mask)
295           
296    def setFitRange(self,qmin=None,qmax=None):
297        """ to set the fit range"""
298        if qmin==0.0:
299            self.qmin = 1e-16
300        elif qmin!=None:                       
301            self.qmin = qmin           
302        if qmax!=None:
303            self.qmax= qmax       
304        self.radius= numpy.sqrt(self.qx_data**2 + self.qy_data**2)
305        self.index_model = ((self.qmin <= self.radius)&(self.radius<= self.qmax))
306        self.mask = (self.index_model) &(self.mask)
307       
308       
309    def getFitRange(self):
310        """
311            @return the range of data.x to fit
312        """
313        return self.qmin, self.qmax
314     
315    def residuals(self, fn): 
316        """
317            @return the residuals
318        """       
319        # use only the data point within ROI range
320        res=(self.data[self.mask] - fn([self.qx_data[self.mask],
321                             self.qy_data[self.mask]]))/self.res_err_data[self.mask]
322        return res
323       
324 
325    def residuals_deriv(self, model, pars=[]):
326        """
327            @return residuals derivatives .
328            @note: in this case just return empty array
329        """
330        return []
331   
332class FitAbort(Exception):
333    """
334        Exception raise to stop the fit
335    """
336    print"Creating fit abort Exception"
337
338
339class SansAssembly:
340    """
341         Sans Assembly class a class wrapper to be call in optimizer.leastsq method
342    """
343    def __init__(self,paramlist,Model=None , Data=None, curr_thread= None):
344        """
345            @param Model: the model wrapper fro sans -model
346            @param Data: the data wrapper for sans data
347        """
348        self.model = Model
349        self.data  = Data
350        self.paramlist=paramlist
351        self.curr_thread= curr_thread
352        self.res=[]
353        self.func_name="Functor"
354    def chisq(self, params):
355        """
356            Calculates chi^2
357            @param params: list of parameter values
358            @return: chi^2
359        """
360        sum = 0
361        for item in self.res:
362            sum += item*item
363        if len(self.res)==0:
364            return None
365        return sum/ len(self.res)
366   
367    def __call__(self,params):
368        """
369            Compute residuals
370            @param params: value of parameters to fit
371        """
372        #import thread
373        self.model.setParams(self.paramlist,params)
374        self.res= self.data.residuals(self.model.eval)
375        #if self.curr_thread != None :
376        #    try:
377        #        self.curr_thread.isquit()
378        #    except:
379        #        raise FitAbort,"stop leastsqr optimizer"   
380        return self.res
381   
382class FitEngine:
383    def __init__(self):
384        """
385            Base class for scipy and park fit engine
386        """
387        #List of parameter names to fit
388        self.paramList=[]
389        #Dictionnary of fitArrange element (fit problems)
390        self.fitArrangeDict={}
391       
392    def _concatenateData(self, listdata=[]):
393        """ 
394            _concatenateData method concatenates each fields of all data contains ins listdata.
395            @param listdata: list of data
396            @return Data: Data is wrapper class for sans plottable. it is created with all parameters
397             of data concatenanted
398            @raise: if listdata is empty  will return None
399            @raise: if data in listdata don't contain dy field ,will create an error
400            during fitting
401        """
402        #TODO: we have to refactor the way we handle data.
403        # We should move away from plottables and move towards the Data1D objects
404        # defined in DataLoader. Data1D allows data manipulations, which should be
405        # used to concatenate.
406        # In the meantime we should switch off the concatenation.
407        #if len(listdata)>1:
408        #    raise RuntimeError, "FitEngine._concatenateData: Multiple data files is not currently supported"
409        #return listdata[0]
410       
411        if listdata==[]:
412            raise ValueError, " data list missing"
413        else:
414            xtemp=[]
415            ytemp=[]
416            dytemp=[]
417            self.mini=None
418            self.maxi=None
419               
420            for item in listdata:
421                data=item.data
422                mini,maxi=data.getFitRange()
423                if self.mini==None and self.maxi==None:
424                    self.mini=mini
425                    self.maxi=maxi
426                else:
427                    if mini < self.mini:
428                        self.mini=mini
429                    if self.maxi < maxi:
430                        self.maxi=maxi
431                       
432                   
433                for i in range(len(data.x)):
434                    xtemp.append(data.x[i])
435                    ytemp.append(data.y[i])
436                    if data.dy is not None and len(data.dy)==len(data.y):   
437                        dytemp.append(data.dy[i])
438                    else:
439                        raise RuntimeError, "Fit._concatenateData: y-errors missing"
440            data= Data(x=xtemp,y=ytemp,dy=dytemp)
441            data.setFitRange(self.mini, self.maxi)
442            return data
443       
444       
445    def set_model(self,model,Uid,pars=[], constraints=[]):
446        """
447            set a model on a given uid in the fit engine.
448            @param model: sans.models type
449            @param Uid :is the key of the fitArrange dictionnary where model is saved as a value
450            @param pars: the list of parameters to fit
451            @param constraints: list of
452                tuple (name of parameter, value of parameters)
453                the value of parameter must be a string to constraint 2 different
454                parameters.
455                Example:
456                we want to fit 2 model M1 and M2 both have parameters A and B.
457                constraints can be:
458                 constraints = [(M1.A, M2.B+2), (M1.B= M2.A *5),...,]
459            @note : pars must contains only name of existing model's paramaters
460        """
461        if model == None:
462            raise ValueError, "AbstractFitEngine: Need to set model to fit"
463       
464        new_model= model
465        if not issubclass(model.__class__, Model):
466            new_model= Model(model)
467       
468        if len(constraints)>0:
469            for constraint in constraints:
470                name, value = constraint
471                try:
472                    new_model.parameterset[ str(name)].set( str(value) )
473                except:
474                    msg= "Fit Engine: Error occurs when setting the constraint"
475                    msg += " %s for parameter %s "%(value, name)
476                    raise ValueError, msg
477               
478        if len(pars) >0:
479            temp=[]
480            for item in pars:
481                if item in new_model.model.getParamList():
482                    temp.append(item)
483                    self.paramList.append(item)
484                else:
485                   
486                    msg = "wrong parameter %s used"%str(item)
487                    msg += "to set model %s. Choose"%str(new_model.model.name)
488                    msg += "parameter name within %s"%str(new_model.model.getParamList())
489                    raise ValueError,msg
490             
491            #A fitArrange is already created but contains dList only at Uid
492            if self.fitArrangeDict.has_key(Uid):
493                self.fitArrangeDict[Uid].set_model(new_model)
494                self.fitArrangeDict[Uid].pars= pars
495            else:
496            #no fitArrange object has been create with this Uid
497                fitproblem = FitArrange()
498                fitproblem.set_model(new_model)
499                fitproblem.pars= pars
500                self.fitArrangeDict[Uid] = fitproblem
501               
502        else:
503            raise ValueError, "park_integration:missing parameters"
504   
505    def set_data(self,data,Uid,smearer=None,qmin=None,qmax=None):
506        """ Receives plottable, creates a list of data to fit,set data
507            in a FitArrange object and adds that object in a dictionary
508            with key Uid.
509            @param data: data added
510            @param Uid: unique key corresponding to a fitArrange object with data
511        """
512        if data.__class__.__name__=='Data2D':
513            fitdata=FitData2D(sans_data2d=data, data=data.data, err_data= data.err_data)
514        else:
515            fitdata=FitData1D(x=data.x, y=data.y , dx= data.dx,dy=data.dy,smearer=smearer)
516       
517        fitdata.setFitRange(qmin=qmin,qmax=qmax)
518        #A fitArrange is already created but contains model only at Uid
519        if self.fitArrangeDict.has_key(Uid):
520            self.fitArrangeDict[Uid].add_data(fitdata)
521        else:
522        #no fitArrange object has been create with this Uid
523            fitproblem= FitArrange()
524            fitproblem.add_data(fitdata)
525            self.fitArrangeDict[Uid]=fitproblem   
526   
527    def get_model(self,Uid):
528        """
529            @param Uid: Uid is key in the dictionary containing the model to return
530            @return  a model at this uid or None if no FitArrange element was created
531            with this Uid
532        """
533        if self.fitArrangeDict.has_key(Uid):
534            return self.fitArrangeDict[Uid].get_model()
535        else:
536            return None
537   
538    def remove_Fit_Problem(self,Uid):
539        """remove   fitarrange in Uid"""
540        if self.fitArrangeDict.has_key(Uid):
541            del self.fitArrangeDict[Uid]
542           
543    def select_problem_for_fit(self,Uid,value):
544        """
545            select a couple of model and data at the Uid position in dictionary
546            and set in self.selected value to value
547            @param value: the value to allow fitting. can only have the value one or zero
548        """
549        if self.fitArrangeDict.has_key(Uid):
550             self.fitArrangeDict[Uid].set_to_fit( value)
551             
552             
553    def get_problem_to_fit(self,Uid):
554        """
555            return the self.selected value of the fit problem of Uid
556           @param Uid: the Uid of the problem
557        """
558        if self.fitArrangeDict.has_key(Uid):
559             self.fitArrangeDict[Uid].get_to_fit()
560   
561class FitArrange:
562    def __init__(self):
563        """
564            Class FitArrange contains a set of data for a given model
565            to perform the Fit.FitArrange must contain exactly one model
566            and at least one data for the fit to be performed.
567            model: the model selected by the user
568            Ldata: a list of data what the user wants to fit
569           
570        """
571        self.model = None
572        self.dList =[]
573        self.pars=[]
574        #self.selected  is zero when this fit problem is not schedule to fit
575        #self.selected is 1 when schedule to fit
576        self.selected = 0
577       
578    def set_model(self,model):
579        """
580            set_model save a copy of the model
581            @param model: the model being set
582        """
583        self.model = model
584       
585    def add_data(self,data):
586        """
587            add_data fill a self.dList with data to fit
588            @param data: Data to add in the list 
589        """
590        if not data in self.dList:
591            self.dList.append(data)
592           
593    def get_model(self):
594        """ @return: saved model """
595        return self.model   
596     
597    def get_data(self):
598        """ @return:  list of data dList"""
599        #return self.dList
600        return self.dList[0] 
601     
602    def remove_data(self,data):
603        """
604            Remove one element from the list
605            @param data: Data to remove from dList
606        """
607        if data in self.dList:
608            self.dList.remove(data)
609    def set_to_fit (self, value=0):
610        """
611           set self.selected to 0 or 1  for other values raise an exception
612           @param value: integer between 0 or 1
613        """
614        self.selected= value
615       
616    def get_to_fit(self):
617        """
618            @return self.selected value
619        """
620        return self.selected
Note: See TracBrowser for help on using the repository browser.