-                      r197ea24
+                      rfbc51ef
 #class Fitting
+import time
+import numpy
+import park
+from scipy import optimize
+from park import fit,fitresult
+from park import assembly
 from sans.guitools.plottables import Data1D
+#from sans.guitools import plottables
 from Loader import Load
+from scipy import optimize
+class SansParameter(park.Parameter):
+    """
+    SANS model parameters for use in the PARK fitting service.
+    The parameter attribute value is redirected to the underlying
+    parameter value in the SANS model.
+    """
+    def __init__(self, name, model):
+         self._model, self._name = model,name
+    def _getvalue(self): return self._model.getParam(self.name)
+    def _setvalue(self,value): self._model.setParam(self.name, value)
+    value = property(_getvalue,_setvalue)
+    def _getrange(self):
+        lo,hi = self._model.details[self.name][1:]
+        if lo is None: lo = -numpy.inf
+        if hi is None: hi = numpy.inf
+        return lo,hi
+    def _setrange(self,r):
+        self._model.details[self.name][1:] = r
+    range = property(_getrange,_setrange)
+class Model(object):
+    """
+        PARK wrapper for SANS models.
+    """
+    def __init__(self, sans_model):
+        self.model = sans_model
+        sansp = sans_model.getParamList()
+        parkp = [SansParameter(p,sans_model) for p in sansp]
+        self.parameterset = park.ParameterSet(sans_model.name,pars=parkp)
+    def eval(self,x):
+        return self.model.run(x)
+class Data(object):
+    """ Wrapper class  for SANS data """
+    def __init__(self, sans_data):
+        self.x= sans_data.x
+        self.y= sans_data.y
+        self.dx= sans_data.dx
+        self.dy= sans_data.dy
+        self.qmin=None
+        self.qmax=None
+    def setFitRange(self,mini=None,maxi=None):
+        """ to set the fit range"""
+        self.qmin=mini
+        self.qmax=maxi
+    def residuals(self, fn):
+        x,y,dy = [numpy.asarray(v) for v in (self.x,self.y,self.dy)]
+        if self.qmin==None and self.qmax==None:
+            return (y - fn(x))/dy
+        else:
+            idx = x>=self.qmin & x <= self.qmax
+            return (y[idx] - fn(x[idx]))/dy[idx]
+    def residuals_deriv(self, model, pars=[]):
+        """ Return residual derivatives .in this case just return empty array"""
+        return []
 class FitArrange:
     def __init__(self):
 …
             self.dList.remove(data)
 class Fitting:
+class ParkFit:
     """
         Performs the Fit.he user determine what kind of data
 …
         self.fitType =None
     def fit_engine(self,word):
+    def createProblem(self,pars={}):
         """
             Check the contraint value and specify what kind of fit to use
+        """
+        self.fitType = word
+        return True
+            return (M1,D1)
+        """
+        mylist=[]
+        for k,value in self.fitArrangeList.iteritems():
+            couple=()
+            model=value.get_model()
+            parameters= self.set_param(model, pars)
+            model = Model(model)
+            #print "model created",model.parameterset[0].value,model.parameterset[1].value
+            # Make all parameters fitting parameters
+            for p in model.parameterset:
+                p.set([-numpy.inf,numpy.inf])
+                #p.set([-10,10])
+            Ldata=value.get_data()
+            data=self._concatenateData(Ldata)
+            #print "this data",data
+            #print "data.residuals in createProblem",Ldata[0].residuals
+            #print "data.residuals in createProblem",data.residuals
+            #couple1=(model,Ldata[0])
+            #mylist.append(couple1)
+            couple=(model,data)
+            mylist.append(couple)
+        #print mylist
+        return mylist
+        #return model,data
     def fit(self,pars, qmin=None, qmax=None):
 …
              Do the fit
         """
+        #for item in self.fitArrangeList.:
+        fitproblem=self.fitArrangeList.values()[0]
+        listdata=[]
+        model = fitproblem.get_model()
+        listdata = fitproblem.get_data()
+        parameters = self.set_param(model,pars)
+        modelList=self.createProblem(pars)
+        #model,data=self.createProblem()
+        #fitness=assembly.Fitness(model,data)
+        problem =  park.Assembly(modelList)
+        #print "problem :",problem[0].parameterset,problem[0].parameterset.fitted
+        #problem[0].parameterset['A'].set([0,1000])
+        #print "problem :",problem[0].parameterset,problem[0].parameterset.fitted
+        fit.fit(problem, handler= fitresult.ConsoleUpdate(improvement_delta=0.1))
+        #return fit.fit(problem)
+        #fit.fit(problem, handler= fitresult.ConsoleUpdate(improvement_delta=0.1))
+    def set_model(self,model,Uid):
+        """ Set model """
+        if self.fitArrangeList.has_key(Uid):
+            self.fitArrangeList[Uid].set_model(model)
+        else:
+            fitproblem= FitArrange()
+            fitproblem.set_model(model)
+            self.fitArrangeList[Uid]=fitproblem
+    def set_data(self,data,Uid):
+        """ Receive plottable and create a list of data to fit"""
+        data=Data(data)
+        if self.fitArrangeList.has_key(Uid):
+            self.fitArrangeList[Uid].add_data(data)
+        else:
+            fitproblem= FitArrange()
+            fitproblem.add_data(data)
+            self.fitArrangeList[Uid]=fitproblem
+    def get_model(self,Uid):
+        """ return list of data"""
+        return self.fitArrangeList[Uid]
+    def set_param(self,model, pars):
+        """ Recieve a dictionary of parameter and save it """
+        parameters=[]
+        if model==None:
+            raise ValueError, "Cannot set parameters for empty model"
+        else:
+            #for key ,value in pars:
+            for key, value in pars.iteritems():
+                param = Parameter(model, key, value)
+                parameters.append(param)
+        return parameters
+    def add_constraint(self, constraint):
+        """ User specify contraint to fit """
+        self.constraint = str(constraint)
+    def get_constraint(self):
+        """ return the contraint value """
+        return self.constraint
+    def set_constraint(self,constraint):
+        """
+            receive a string as a constraint
+            @param constraint: a string used to constraint some parameters to get a
+                specific value
+        """
+        self.constraint= constraint
+    def _concatenateData(self, listdata=[]):
+        """ concatenate each fields of all Data contains ins listdata
+         return data
+        """
         if listdata==[]:
             raise ValueError, " data list missing"
         else:
-            # Do the fit with more than one data set and one model
             xtemp=[]
             ytemp=[]
             dytemp=[]
+            resid=[]
+            resid_deriv=[]
             for data in listdata:
                 for i in range(len(data.x)):
 …
                     if not data.dy[i] in dytemp:
                         dytemp.append(data.dy[i])
+            if qmin==None:
+                qmin= min(xtemp)
+            if qmax==None:
+                qmax= max(xtemp)
+            chisqr, out, cov = fitHelper(model,parameters, xtemp,ytemp, dytemp ,qmin,qmax)
+            return chisqr, out, cov
+    def set_model(self,model,Uid):
+        """ Set model """
+        if self.fitArrangeList.has_key(Uid):
+            self.fitArrangeList[Uid].set_model(model)
+        else:
+            fitproblem= FitArrange()
+            fitproblem.set_model(model)
+            self.fitArrangeList[Uid]=fitproblem
+    def set_data(self,data,Uid):
+        """ Receive plottable and create a list of data to fit"""
+        if self.fitArrangeList.has_key(Uid):
+            self.fitArrangeList[Uid].add_data(data)
+        else:
+            fitproblem= FitArrange()
+            fitproblem.add_data(data)
+            self.fitArrangeList[Uid]=fitproblem
+    def get_model(self,Uid):
+        """ return list of data"""
+        return self.fitArrangeList[Uid]
+    def set_param(self,model, pars):
+        """ Recieve a dictionary of parameter and save it """
+        parameters=[]
+        if model==None:
+            raise ValueError, "Cannot set parameters for empty model"
+        else:
+            #for key ,value in pars:
+            for key, value in pars.iteritems():
+                param = Parameter(model, key, value)
+                parameters.append(param)
+        return parameters
+    def add_constraint(self, constraint):
+        """ User specify contraint to fit """
+        self.constraint = str(constraint)
+    def get_constraint(self):
+        """ return the contraint value """
+        return self.constraint
+    def set_constraint(self,constraint):
+        """
+            receive a string as a constraint
+            @param constraint: a string used to constraint some parameters to get a
+                specific value
+        """
+        self.constraint= constraint
+            newplottable= Data1D(xtemp,ytemp,None,dytemp)
+            newdata=Data(newplottable)
+            #print "this is new data",newdata.dy
+            return newdata
 class Parameter:
     """
 …
         """
         return self.model.getParam(self.name)
+class FitHelper:
+    def __init__(self,model, pars, x, y, err_y ,qmin=None, qmax=None):
+        self.x = x
+        self.y = y
+        self.model = model
+        self.err_y = err_y
+        self.qmin = qmin
+        self.qmax= qmax
+        self.pars = pars
+    def __call__(self, params):
+        i = 0
+        for p in self.pars:
+            p.set(params[i])
+            i += 1
+        residuals = []
+        for j in range(len(self.x)):
+            if self.x[j]>self.qmin and self.x[j]<self.qmax:
+                residuals.append( ( self.y[j] - self.model.runXY(self.x[j]) ) / self.err_y[j] )
+        return residuals
+def fitHelper(model, pars, x, y, err_y ,qmin=None, qmax=None):
+    """
+        Fit function
+        @param model: sans model object
+        @param pars: list of parameters
+        @param x: vector of x data
+        @param y: vector of y data
+        @param err_y: vector of y errors
+    """
+    f = FitHelper(model, pars, x, y, err_y ,qmin, qmax)
+    def ff(params):
+        """
+            Calculates the vector of residuals for each point
+            in y for a given set of input parameters.
+            @param params: list of parameter values
+            @return: vector of residuals
+        """
+        i = 0
+        for p in pars:
+            p.set(params[i])
+            i += 1
+        residuals = []
+        for j in range(len(x)):
+            if x[j]>qmin and x[j]<qmax:
+                residuals.append( ( y[j] - model.runXY(x[j]) ) / err_y[j] )
+        return residuals
+    def chi2(params):
+        """
+            Calculates chi^2
+            @param params: list of parameter values
+            @return: chi^2
+        """
+        sum = 0
+        res = f(params)
+        for item in res:
+            sum += item*item
+        return sum
+    p = [param() for param in pars]
+    out, cov_x, info, mesg, success = optimize.leastsq(f, p, full_output=1, warning=True)
+    print info, mesg, success
+    # Calculate chi squared
+    if len(pars)>1:
+        chisqr = chi2(out)
+    elif len(pars)==1:
+        chisqr = chi2([out])
+    return chisqr, out, cov_x
 …
     data1.name = "data1"
     load.load_data(data1)
     Fit =Fitting()
     from LineModel import LineModel
+    fitter =ParkFit()
+    from sans.guitools.LineModel import LineModel
     model  = LineModel()
+    Fit.set_model(model,1)
+    Fit.set_data(data1,1)
+    chisqr, out, cov=Fit.fit({'A':2,'B':1},None,None)
+    print"fit only one data",chisqr, out, cov
+    # test fit with 2 data and one model
+    Fit =Fitting()
+    Fit.set_model(model,2 )
+    load.set_filename("testdata1.txt")
+    load.set_values()
+    data2 = Data1D(x=[], y=[], dx=None,dy=None)
+    data2.name = "data2"
+    load.load_data(data2)
+    Fit.set_data(data2,2)
+    load.set_filename("testdata2.txt")
+    load.set_values()
+    data3 = Data1D(x=[], y=[], dx=None,dy=None)
+    data3.name = "data2"
+    load.load_data(data3)
+    Fit.set_data(data3,2)
+    chisqr, out, cov=Fit.fit({'A':2,'B':1},None,None)
+    print"fit two data",chisqr, out, cov
+    fitter.set_model(model,1)
+    fitter.set_data(data1,1)
+    print"PARK fit result \n",fitter.fit({'A':2,'B':1},None,None)

Note: See TracChangeset for help on using the changeset viewer.

SasView

Changeset fbc51ef in sasview

Legend:

park_integration/test/ParkFitting.py

Download in other formats: