<s> """<STR_LIT>""" <EOL> from numpy import ( abs , array , asarray , dot , eye , ndarray , ones_like , <EOL> sqrt , take , transpose , triu ) <EOL> from numpy . dual import inv <EOL> from numpy . linalg import LinAlgError <EOL> from scipy . optimize import _minpack <EOL> from scipy . optimize . minpack import _check_func <EOL> HAS_SCALAR_MIN = False <EOL> try : <EOL> from scipy . optimize import minimize as scipy_minimize <EOL> HAS_SCALAR_MIN = True <EOL> except ImportError : <EOL> pass <EOL> from . import uncertainties <EOL> from . parameter import isParameter <EOL> try : <EOL> from . . symboltable import Group <EOL> except : <EOL> Group = None <EOL> class MinimizerException ( Exception ) : <EOL> """<STR_LIT>""" <EOL> def __init__ ( self , msg ) : <EOL> Exception . __init__ ( self ) <EOL> self . msg = msg <EOL> def __str__ ( self ) : <EOL> return "<STR_LIT>" % ( self . msg ) <EOL> def larcheval_with_uncertainties ( * vals , ** kwargs ) : <EOL> """<STR_LIT>""" <EOL> _obj = kwargs . get ( '<STR_LIT>' , None ) <EOL> _pars = kwargs . get ( '<STR_LIT>' , None ) <EOL> _names = kwargs . get ( '<STR_LIT>' , None ) <EOL> _larch = kwargs . get ( '<STR_LIT>' , None ) <EOL> if ( _obj is None or _pars is None or <EOL> _names is None or _larch is None or <EOL> _obj . _ast is None ) : <EOL> return <NUM_LIT:0> <EOL> for val , name in zip ( vals , _names ) : <EOL> _pars [ name ] . _val = val <EOL> result = _larch . eval ( _obj . _ast ) <EOL> if isParameter ( result ) : <EOL> result = result . value <EOL> return result <EOL> wrap_ueval = uncertainties . wrap ( larcheval_with_uncertainties ) <EOL> def eval_stderr ( obj , uvars , _names , _pars , _larch ) : <EOL> """<STR_LIT>""" <EOL> if not isParameter ( obj ) : <EOL> return <EOL> if obj . _ast is None : <EOL> return <EOL> uval = wrap_ueval ( * uvars , _obj = obj , _names = _names , <EOL> _pars = _pars , _larch = _larch ) <EOL> try : <EOL> obj . stderr = uval . std_dev ( ) <EOL> obj . _uval = uval <EOL> except : <EOL> obj . stderr = <NUM_LIT:0> <EOL> obj . _uval = None <EOL> def leastsq ( func , x0 , args = ( ) , Dfun = None , ftol = <NUM_LIT> , xtol = <NUM_LIT> , <EOL> gtol = <NUM_LIT> , maxfev = <NUM_LIT:0> , epsfcn = None , factor = <NUM_LIT:100> , diag = None ) : <EOL> """<STR_LIT>""" <EOL> x0 = asarray ( x0 ) . flatten ( ) <EOL> n = len ( x0 ) <EOL> if not isinstance ( args , tuple ) : <EOL> args = ( args , ) <EOL> shape = _check_func ( '<STR_LIT>' , '<STR_LIT>' , func , x0 , args , n ) <EOL> if isinstance ( shape , tuple ) and len ( shape ) > <NUM_LIT:1> : <EOL> shape = shape [ <NUM_LIT:0> ] <EOL> m = shape [ <NUM_LIT:0> ] <EOL> if n > m : <EOL> raise TypeError ( '<STR_LIT>' % ( n , m ) ) <EOL> if maxfev == <NUM_LIT:0> : <EOL> maxfev = <NUM_LIT:200> * ( n + <NUM_LIT:1> ) <EOL> if epsfcn is None : <EOL> epsfcn = <NUM_LIT> <EOL> if Dfun is None : <EOL> retval = _minpack . _lmdif ( func , x0 , args , <NUM_LIT:1> , ftol , xtol , <EOL> gtol , maxfev , epsfcn , factor , diag ) <EOL> else : <EOL> _check_func ( '<STR_LIT>' , '<STR_LIT>' , Dfun , x0 , args , n , ( m , n ) ) <EOL> retval = _minpack . _lmder ( func , Dfun , x0 , args , <NUM_LIT:1> , <NUM_LIT:0> , ftol , xtol , <EOL> gtol , maxfev , factor , diag ) <EOL> errors = { <NUM_LIT:0> : [ "<STR_LIT>" , TypeError ] , <EOL> <NUM_LIT:1> : [ "<STR_LIT>" <EOL> "<STR_LIT>" % ftol , None ] , <EOL> <NUM_LIT:2> : [ "<STR_LIT>" <EOL> "<STR_LIT>" % xtol , None ] , <EOL> <NUM_LIT:3> : [ "<STR_LIT>" <EOL> "<STR_LIT>" <EOL> "<STR_LIT>" <EOL> "<STR_LIT>" % ( ftol , xtol ) , None ] , <EOL> <NUM_LIT:4> : [ "<STR_LIT>" <EOL> "<STR_LIT>" <EOL> "<STR_LIT>" % gtol , None ] , <EOL> <NUM_LIT:5> : [ "<STR_LIT>" <EOL> "<STR_LIT>" % maxfev , ValueError ] , <EOL> <NUM_LIT:6> : [ "<STR_LIT>" <EOL> "<STR_LIT>" "<STR_LIT>" % ftol , ValueError ] , <EOL> <NUM_LIT:7> : [ "<STR_LIT>" <EOL> "<STR_LIT>" % xtol , ValueError ] , <EOL> <NUM_LIT:8> : [ "<STR_LIT>" <EOL> "<STR_LIT>" <EOL> "<STR_LIT>" % gtol , ValueError ] , <EOL> '<STR_LIT>' : [ "<STR_LIT>" , TypeError ] } <EOL> info = retval [ - <NUM_LIT:1> ] <EOL> mesg = errors [ info ] [ <NUM_LIT:0> ] <EOL> cov_x = None <EOL> if info in [ <NUM_LIT:1> , <NUM_LIT:2> , <NUM_LIT:3> , <NUM_LIT:4> ] : <EOL> perm = take ( eye ( n ) , retval [ <NUM_LIT:1> ] [ '<STR_LIT>' ] - <NUM_LIT:1> , <NUM_LIT:0> ) <EOL> r = triu ( transpose ( retval [ <NUM_LIT:1> ] [ '<STR_LIT>' ] ) [ : n , : ] ) <EOL> R = dot ( r , perm ) <EOL> try : <EOL> cov_x = inv ( dot ( transpose ( R ) , R ) ) <EOL> except ( LinAlgError , ValueError ) : <EOL> pass <EOL> return ( retval [ <NUM_LIT:0> ] , cov_x ) + retval [ <NUM_LIT:1> : - <NUM_LIT:1> ] + ( mesg , info ) <EOL> class Minimizer ( object ) : <EOL> """<STR_LIT>""" <EOL> err_maxfev = """<STR_LIT>""" <EOL> def __init__ ( self , fcn , params , fcn_args = None , fcn_kws = None , <EOL> scale_covar = True , toler = <NUM_LIT> , <EOL> _larch = None , jacfcn = None , ** kws ) : <EOL> self . userfcn = fcn <EOL> self . paramgroup = params <EOL> self . userargs = fcn_args <EOL> if self . userargs is None : <EOL> self . userargs = [ ] <EOL> self . userkws = fcn_kws <EOL> if self . userkws is None : <EOL> self . userkws = { } <EOL> self . _larch = _larch <EOL> self . toler = toler <EOL> self . scale_covar = scale_covar <EOL> self . kws = kws <EOL> self . jacfcn = jacfcn <EOL> self . __prepared = False <EOL> self . prepare_fit ( ) <EOL> def __update_params ( self , fvars ) : <EOL> """<STR_LIT>""" <EOL> if not self . __prepared : <EOL> print ( '<STR_LIT>' ) <EOL> group = self . paramgroup <EOL> for name , val in zip ( self . var_names , fvars ) : <EOL> par = getattr ( group , name ) <EOL> par . _val = par . _from_internal ( val ) <EOL> def __residual ( self , fvars ) : <EOL> """<STR_LIT>""" <EOL> self . __update_params ( fvars ) <EOL> return self . userfcn ( self . paramgroup , * self . userargs , ** self . userkws ) <EOL> def __jacobian ( self , fvars ) : <EOL> """<STR_LIT>""" <EOL> self . __update_params ( fvars ) <EOL> return self . jacfcn ( self . paramgroup , * self . userargs , ** self . userkws ) <EOL> def prepare_fit ( self , force = False ) : <EOL> """<STR_LIT>""" <EOL> if self . __prepared and not force : <EOL> return <EOL> if self . _larch . symtable . isgroup ( self . paramgroup ) : <EOL> self . _larch . symtable . _sys . paramGroup = self . paramgroup <EOL> else : <EOL> self . _larch . write . write ( '<STR_LIT>' ) <EOL> return <EOL> self . var_names = [ ] <EOL> self . vars = [ ] <EOL> self . nvarys = <NUM_LIT:0> <EOL> for name in dir ( self . paramgroup ) : <EOL> par = getattr ( self . paramgroup , name ) <EOL> if isParameter ( par ) : <EOL> val0 = par . setup_bounds ( ) <EOL> if par . vary : <EOL> self . var_names . append ( name ) <EOL> self . vars . append ( val0 ) <EOL> elif par . expr is not None : <EOL> par . _getval ( ) <EOL> if not hasattr ( par , '<STR_LIT:name>' ) or par . name is None : <EOL> par . name = name <EOL> self . nvarys = len ( self . vars ) <EOL> self . __prepared = True <EOL> def leastsq ( self , ** kws ) : <EOL> """<STR_LIT>""" <EOL> self . prepare_fit ( force = True ) <EOL> toler = self . toler <EOL> lskws = dict ( xtol = toler , ftol = toler , <EOL> gtol = toler , maxfev = <NUM_LIT:1000> * ( self . nvarys + <NUM_LIT:1> ) , Dfun = None ) <EOL> lskws . update ( self . kws ) <EOL> lskws . update ( kws ) <EOL> if lskws [ '<STR_LIT>' ] is not None : <EOL> self . jacfcn = lskws [ '<STR_LIT>' ] <EOL> lskws [ '<STR_LIT>' ] = self . __jacobian <EOL> lsout = leastsq ( self . __residual , self . vars , ** lskws ) <EOL> del self . vars <EOL> _best , cov , infodict , errmsg , ier = lsout <EOL> resid = infodict [ '<STR_LIT>' ] <EOL> ndata = len ( resid ) <EOL> chisqr = ( resid ** <NUM_LIT:2> ) . sum ( ) <EOL> nfree = ndata - self . nvarys <EOL> redchi = chisqr / nfree <EOL> group = self . paramgroup <EOL> grad = ones_like ( _best ) <EOL> vbest = ones_like ( _best ) <EOL> named_params = { } <EOL> for ivar , name in enumerate ( self . var_names ) : <EOL> named_params [ name ] = par = getattr ( group , name ) <EOL> grad [ ivar ] = par . scale_gradient ( _best [ ivar ] ) <EOL> vbest [ ivar ] = par . value <EOL> par . stderr = <NUM_LIT:0> <EOL> par . correl = { } <EOL> par . _uval = None <EOL> infodict [ '<STR_LIT>' ] = transpose ( transpose ( infodict [ '<STR_LIT>' ] ) / <EOL> take ( grad , infodict [ '<STR_LIT>' ] - <NUM_LIT:1> ) ) <EOL> rvec = dot ( triu ( transpose ( infodict [ '<STR_LIT>' ] ) [ : self . nvarys , : ] ) , <EOL> take ( eye ( self . nvarys ) , infodict [ '<STR_LIT>' ] - <NUM_LIT:1> , <NUM_LIT:0> ) ) <EOL> try : <EOL> cov = inv ( dot ( transpose ( rvec ) , rvec ) ) <EOL> except ( LinAlgError , ValueError ) : <EOL> cov = None <EOL> if cov is not None : <EOL> if self . scale_covar : <EOL> cov = cov * chisqr / nfree <EOL> try : <EOL> uvars = uncertainties . correlated_values ( vbest , cov ) <EOL> except ( LinAlgError , ValueError ) : <EOL> cov , uvars = None , None <EOL> group . covar_vars = self . var_names <EOL> group . covar = cov <EOL> if uvars is not None : <EOL> for iv , name in enumerate ( self . var_names ) : <EOL> p = named_params [ name ] <EOL> p . stderr = uvars [ iv ] . std_dev ( ) <EOL> p . _uval = uvars [ iv ] <EOL> p . correl = { } <EOL> for jv , name2 in enumerate ( self . var_names ) : <EOL> if jv != iv : <EOL> p . correl [ name2 ] = ( cov [ iv , jv ] / <EOL> ( p . stderr * sqrt ( cov [ jv , jv ] ) ) ) <EOL> for nam in dir ( self . paramgroup ) : <EOL> obj = getattr ( self . paramgroup , nam ) <EOL> eval_stderr ( obj , uvars , self . var_names , <EOL> named_params , self . _larch ) <EOL> for uval , nam in zip ( uvars , self . var_names ) : <EOL> named_params [ nam ] . _val = uval . nominal_value <EOL> if self . _larch . error : <EOL> self . _larch . error = [ ] <EOL> message = '<STR_LIT>' <EOL> if ier == <NUM_LIT:0> : <EOL> message = '<STR_LIT>' <EOL> elif ier == <NUM_LIT:5> : <EOL> message = self . err_maxfev % lskws [ '<STR_LIT>' ] <EOL> elif ier > <NUM_LIT:5> : <EOL> message = '<STR_LIT>' <EOL> if cov is None : <EOL> message = '<STR_LIT>' % message <EOL> ofit = group <EOL> if Group is not None : <EOL> ofit = group . fit_details = Group ( ) <EOL> ofit . method = '<STR_LIT>' <EOL> ofit . fjac = infodict [ '<STR_LIT>' ] <EOL> ofit . fvec = infodict [ '<STR_LIT>' ] <EOL> ofit . qtf = infodict [ '<STR_LIT>' ] <EOL> ofit . ipvt = infodict [ '<STR_LIT>' ] <EOL> ofit . nfev = infodict [ '<STR_LIT>' ] <EOL> ofit . status = ier <EOL> ofit . message = errmsg <EOL> ofit . success = ier in [ <NUM_LIT:1> , <NUM_LIT:2> , <NUM_LIT:3> , <NUM_LIT:4> ] <EOL> ofit . toler = self . toler <EOL> group . residual = resid <EOL> group . message = message <EOL> group . chi_square = chisqr <EOL> group . chi_reduced = redchi <EOL> group . nvarys = self . nvarys <EOL> group . nfree = nfree <EOL> group . errorbars = cov is not None <EOL> return ier <EOL> def scalar_minimize ( self , method = '<STR_LIT>' , ** kws ) : <EOL> """<STR_LIT>""" <EOL> if not HAS_SCALAR_MIN : <EOL> raise NotImplementedError <EOL> self . prepare_fit ( ) <EOL> maxfev = <NUM_LIT:1000> * ( self . nvarys + <NUM_LIT:1> ) <EOL> opts = { '<STR_LIT>' : maxfev } <EOL> if method not in ( '<STR_LIT>' , '<STR_LIT>' , '<STR_LIT>' ) : <EOL> opts [ '<STR_LIT>' ] = maxfev <EOL> fmin_kws = dict ( method = method , tol = self . toler , options = opts ) <EOL> fmin_kws . update ( self . kws ) <EOL> fmin_kws . update ( kws ) <EOL> def penalty ( parvals ) : <EOL> "<STR_LIT>" <EOL> r = self . __residual ( parvals ) <EOL> if isinstance ( r , ndarray ) : <EOL> r = ( r * r ) . sum ( ) <EOL> return r <EOL> ret = scipy_minimize ( penalty , self . vars , ** fmin_kws ) <EOL> del self . vars <EOL> resid = self . __residual ( ret . x ) <EOL> ndata = len ( resid ) <EOL> chisqr = ( resid ** <NUM_LIT:2> ) . sum ( ) <EOL> nfree = ( ndata - self . nvarys ) <EOL> redchi = chisqr / nfree <EOL> ofit = group = self . paramgroup <EOL> if Group is not None : <EOL> ofit = group . fit_details = Group ( ) <EOL> ofit . method = method <EOL> ofit . nfev = ret . nfev <EOL> ofit . success = ret . success <EOL> ofit . status = ret . status <EOL> group . nvarys = self . nvarys <EOL> group . nfree = nfree <EOL> group . residual = resid <EOL> group . message = ret . message <EOL> group . chi_square = chisqr <EOL> group . chi_reduced = redchi <EOL> group . errorbars = False <EOL> def minimize ( fcn , group , args = None , kws = None , method = '<STR_LIT>' , <EOL> _larch = None , ** fit_kws ) : <EOL> """<STR_LIT>""" <EOL> if not _larch . symtable . isgroup ( group ) : <EOL> return '<STR_LIT>' <EOL> fit = Minimizer ( fcn , group , fcn_args = args , fcn_kws = kws , <EOL> scale_covar = True , _larch = _larch , ** fit_kws ) <EOL> _scalar_methods = { '<STR_LIT>' : '<STR_LIT>' , <EOL> '<STR_LIT>' : '<STR_LIT>' , <EOL> '<STR_LIT>' : '<STR_LIT>' , <EOL> '<STR_LIT>' : '<STR_LIT>' , <EOL> '<STR_LIT>' : '<STR_LIT>' , <EOL> '<STR_LIT>' : '<STR_LIT>' , <EOL> '<STR_LIT>' : '<STR_LIT>' , <EOL> '<STR_LIT>' : '<STR_LIT>' , <EOL> '<STR_LIT>' : '<STR_LIT>' , <EOL> '<STR_LIT>' : '<STR_LIT>' } <EOL> meth = method . lower ( ) <EOL> meth_found = False <EOL> if HAS_SCALAR_MIN : <EOL> for name , method in _scalar_methods . items ( ) : <EOL> if meth . startswith ( name ) : <EOL> meth_found = True <EOL> fit . scalar_minimize ( method = method ) <EOL> if not meth_found : <EOL> fit . leastsq ( ) <EOL> return fit <EOL> def fit_report ( group , show_correl = True , min_correl = <NUM_LIT:0.1> , precision = None , <EOL> _larch = None , ** kws ) : <EOL> """<STR_LIT>""" <EOL> if not _larch . symtable . isgroup ( group ) : <EOL> print ( '<STR_LIT>' ) <EOL> return <EOL> _larch . symtable . _sys . paramGroup = group <EOL> topline = '<STR_LIT>' <EOL> header = '<STR_LIT>' <EOL> exprformat = '<STR_LIT>' <EOL> out = [ topline ] <EOL> varformat = '<STR_LIT>' <EOL> fmt_sca = "<STR_LIT>" <EOL> fmt_err = "<STR_LIT>" <EOL> if precision is not None : <EOL> varformat = '<STR_LIT>' % precision <EOL> fmt_sca = "<STR_LIT>" % precision <EOL> fmt_err = "<STR_LIT>" % ( precision , precision ) <EOL> npts = len ( group . residual ) <EOL> ofit = getattr ( group , '<STR_LIT>' , None ) <EOL> if ofit is None : ofit = group <EOL> methodname = getattr ( ofit , '<STR_LIT>' , '<STR_LIT>' ) <EOL> success = getattr ( ofit , '<STR_LIT:success>' , False ) <EOL> if success : <EOL> subtitle = '<STR_LIT>' <EOL> else : <EOL> subtitle = '<STR_LIT>' <EOL> subtitle = "<STR_LIT>" % ( subtitle , methodname ) <EOL> out . append ( header % ( '<STR_LIT>' , subtitle ) ) <EOL> if hasattr ( group , '<STR_LIT:message>' ) : <EOL> out . append ( '<STR_LIT>' % ( group . message ) ) <EOL> out . append ( '<STR_LIT>' % ( npts , <EOL> group . nvarys , <EOL> group . nfree ) ) <EOL> if hasattr ( ofit , '<STR_LIT>' ) : <EOL> out . append ( '<STR_LIT>' % ( ofit . nfev ) ) <EOL> if hasattr ( group , '<STR_LIT>' ) : <EOL> out . append ( ( '<STR_LIT>' % fmt_sca ) % <EOL> ( group . chi_square ) ) <EOL> if hasattr ( group , '<STR_LIT>' ) : <EOL> out . append ( ( '<STR_LIT>' % fmt_sca ) % <EOL> ( group . chi_reduced ) ) <EOL> out . append ( '<STR_LIT:U+0020>' ) <EOL> out . append ( header % ( '<STR_LIT>' , '<STR_LIT>' ) ) <EOL> exprs = [ ] <EOL> for name in dir ( group ) : <EOL> var = getattr ( group , name ) <EOL> if len ( name ) < <NUM_LIT> : <EOL> name = ( name + '<STR_LIT:U+0020>' * <NUM_LIT> ) [ : <NUM_LIT> ] <EOL> if isParameter ( var ) : <EOL> sval = fmt_sca % var . value <EOL> if var . stderr is not None : <EOL> sval = fmt_err % ( var . value , var . stderr ) <EOL> if var . vary : <EOL> out . append ( varformat % ( name , sval , var . _initval ) ) <EOL> elif var . expr is not None : <EOL> exprs . append ( exprformat % ( name , sval , var . expr ) ) <EOL> if len ( exprs ) > <NUM_LIT:0> : <EOL> out . append ( header % ( '<STR_LIT>' , '<STR_LIT>' ) ) <EOL> out . extend ( exprs ) <EOL> covar_vars = getattr ( group , '<STR_LIT>' , [ ] ) <EOL> if show_correl and len ( covar_vars ) > <NUM_LIT:0> : <EOL> subtitle = '<STR_LIT>' % min_correl <EOL> out . append ( '<STR_LIT:U+0020>' ) <EOL> out . append ( header % ( '<STR_LIT>' , subtitle ) ) <EOL> correls = { } <EOL> for i , name in enumerate ( covar_vars ) : <EOL> par = getattr ( group , name ) <EOL> if not par . vary : <EOL> continue <EOL> if hasattr ( par , '<STR_LIT>' ) and par . correl is not None : <EOL> for name2 in covar_vars [ i + <NUM_LIT:1> : ] : <EOL> if name != name2 and name2 in par . correl : <EOL> correls [ "<STR_LIT>" % ( name , name2 ) ] = par . correl [ name2 ] <EOL> sort_correl = sorted ( correls . items ( ) , key = lambda it : abs ( it [ <NUM_LIT:1> ] ) ) <EOL> sort_correl . reverse ( ) <EOL> for name , val in sort_correl : <EOL> if abs ( val ) < min_correl : <EOL> break <EOL> if len ( name ) < <NUM_LIT:20> : <EOL> name = ( name + '<STR_LIT:U+0020>' * <NUM_LIT:20> ) [ : <NUM_LIT:20> ] <EOL> out . append ( '<STR_LIT>' % ( name , val ) ) <EOL> out . append ( '<STR_LIT:=>' * len ( topline ) ) <EOL> return '<STR_LIT:\n>' . join ( out ) </s>
