[SciPy-user] mpfit

Steve Schmerler elcorto at gmx.net
Mon Feb 21 09:54:46 CST 2005


Here the call-mpfit stuff. Maybe this helps.

Cheers, Steve

-- 
There are three types of people in this world: those who make things 
happen, those who watch things happen and those who wonder what 
happened. - Mary Kay Ash
-------------- next part --------------
#!/usr/bin/env python
# -*- coding: utf-8 -*-

# here are the essential parts of a small function which calls mpfit
# ==================================================================


# [...]


# assuming you have the data arrays 't' and 'y' in the namespace
def residuals(x, y = y, t = t, fjac = None, err = None):
    """
    y-func(t,x) - array with residuals 
    """
    status = 0
    return [status, array(y-func(t,x), 'd')]

residualFunction = residuals


def parametersInfos(initial_guess, limits, fixed, max_steps = None):
    """
    initial_guess -- array of len(x), the starting point
    limits -- n x 2 list of lists or array, n is the number of parameters to be fitted (== len(x)),
              limits = [[l, u], ..., [l, u]] where limits[i] contains the lower (l) and upper (u) bounds of the i-th parameter
    fixed -- a list of tuples of the form [(parameter_index,value),(),...,()] telling that parameter x[parameter_index] is  
             to be held constant at 'value'         
    max_steps -- list of len(x), if max_steps[i]  != 0: max iteration step length for x[i]
    """
    # no parameter bounds
    if (limits == None):
        parameter_infos = [{'mpside':0} for i in range(len(initial_guess))]
    else:
        parameter_infos = [{'limited':[1,1],'limits':limits[i],'mpside':0} for i in range(len(initial_guess))]        
    for j in range(len(parameter_infos)):
        # enter start values in dictionary 'parameter_infos'
        parameter_infos[j]['value'] = initial_guess[j]
        if max_steps:
            # set max iteration step for parameter j; max_step[j] = 0 means no max step
            parameter_infos[j]['mpmaxstep'] = max_steps[j]                
    # change 'parameter_infos' - entry if a parameter is fixed
    if (fixed != None):
        for f in fixed:
            # f is a tuple
            i = f[0]
            parameter_infos[i] = {'fixed':1}
            parameter_infos[i]['value'] = float(f[1])
    return parameter_infos


# for printing the termination status of mpfit
mpfit_messages = {0: "Improper input parameters in mpfit.",
                  1: "Both actual and predicted relative reductions in the sum of squares are at most ftol.",
                  2: "Relative error between two consecutive iterates is at most xtol.",
                  3: "Conditions for status = 1 and status = 2 both hold:\n\
                      Both actual and predicted relative reductions in the sum of squares are at most ftol.\n\
                      Relative error between two consecutive iterates is at most xtol.",
                  4: "The cosine of the angle between fvec and any column of the jacobian is at most gtol in absolute value.",
                  5: "The maximum number of iterations has been reached.",
                  6: "ftol is too small. No further reduction in the sum of squares is possible.",
                  7: "xtol is too small. No further improvement in the approximate solution x is possible.",
                  8: "gtol is too small. fvec is orthogonal to the columns of the jacobian to machine precision."}



parameter_infos = parametersInfos(initial_guess, limits, fixed, max_steps = max_steps)
fit = mpfit.mpfit(residualFunction, parinfo=parameter_infos, quiet=1, ftol=1e-15, xtol=1e-15, gtol=1e-15, maxiter=niter)    


More information about the SciPy-user mailing list