Python GaussianModel.make_params方法代码示例

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def GaussConst(signal, guess):
    
    amp, centre, stdev, offset = guess
    
    data = np.array([range(len(signal)), signal]).T
    X = data[:,0]
    Y = data[:,1]

    gauss_mod = GaussianModel(prefix='gauss_')
    const_mod = ConstantModel(prefix='const_')
    
    pars = gauss_mod.make_params(center=centre, sigma=stdev, amplitude=amp)
    pars += const_mod.guess(Y, x=X)
    pars['gauss_center'].min = centre - 5.
    pars['gauss_center'].max = centre + 5.
    pars['gauss_sigma'].max = stdev + 5.
    
    mod = gauss_mod + const_mod
    result = mod.fit(Y, pars, x=X)
    
    fwhm = result.best_values['gauss_sigma'] #* 2.3548
    centr = result.best_values['gauss_center']
    
    # Values within two stdevs i.e. 95%
    pl.plot(np.repeat(centr - fwhm * 2, len(Y)),
            np.arange(len(Y)), 'b-')
    pl.plot(np.repeat(centr + fwhm * 2, len(Y)),
            np.arange(len(Y)), 'b-')
    
    return X, result.best_fit, result.best_values['gauss_sigma'] * 4

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def lmfit_ngauss(x,y, *params):
  params = params[0]
  mods = []
  prefixes = []
  for i in range(0, len(params), 3):
    pref = "g%02i_" % (i/3)
    gauss_i = GaussianModel(prefix=pref)

    if i == 0:
      pars = gauss_i.guess(y, x=x)
    else:
      pars.update(gauss_i.make_params())

    A = params[i]
    l_cen = params[i+1]
    sigma = params[i+2]

    pars[pref+'amplitude'].set(A)
    pars[pref+'center'].set(l_cen)
    pars[pref+'sigma'].set(sigma)

    mods.append(gauss_i)
    prefixes.append(pref)

  mod = mods[0]

  if len(mods) > 1:
    for m in mods[1:]:
      mod += m

  print mod

  init = mod.eval(pars, x=x)
  out = mod.fit(y, pars, x=x)
  return mod, out, init

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def test_2gaussians():

  x = np.linspace(0.0, 10.0, num=1000)
  y = gaussian(x, -1, 3, 0.75) + gaussian(x, -0.5, 5, 0.8) + np.random.normal(0, 0.01, x.shape[0])

  gauss1  = GaussianModel(prefix='g1_')
  
  pars = gauss1.guess(y, x=x)
  pars['g1_amplitude'].set(-0.9)
  pars['g1_center'].set(2.5)
  pars['g1_sigma'].set(0.5)

  gauss2  = GaussianModel(prefix='g2_')
  pars.update(gauss2.make_params())
  pars['g2_amplitude'].set(-0.4)
  pars['g2_center'].set(5)
  pars['g2_sigma'].set(0.5)

  mod = gauss1 + gauss2

  init = mod.eval(pars, x=x)

  plt.plot(x, y)
  plt.plot(x, init, 'k--')

  out = mod.fit(y, pars, x=x)

  print(out.fit_report(min_correl=0.5))

  plt.plot(x, out.best_fit, 'r-')
  plt.show()

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def call_gauss(x, y, cen, count, pars):
	label='g'+str(count)+'_'
	gauss = GaussianModel(prefix=label)
	pars.update(gauss.make_params())
	pars[label+'center'].set(cen, min=cen-0.01, max=cen+0.01)
	pars[label+'amplitude'].set(-0.5, min=-10., max=0.0001)
	pars[label+'sigma'].set(0.1, min=0.005, max=0.25)
	return gauss

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def xrf_calib_init_roi(mca, roiname):
    """initial calibration step for MCA:
    find energy locations for one ROI
    """
    if not isLarchMCAGroup(mca):
        print( 'Not a valid MCA')
        return
    energy = 1.0*mca.energy
    chans = 1.0*np.arange(len(energy))
    counts = mca.counts
    bgr = getattr(mca, 'bgr', None)
    if bgr is not None:
        counts = counts - bgr
    if not hasattr(mca, 'init_calib'):
        mca.init_calib = OrderedDict()

    roi = None
    for xroi in mca.rois:
        if xroi.name == roiname:
            roi = xroi
            break
    if roi is None:
        return
    words = roiname.split()
    elem = words[0].title()
    family = 'Ka'
    if len(words) > 1:
        family = words[1].title()
    if family == 'Lb':
        family = 'Lb1'
    eknown = xray_line(elem, family)[0]/1000.0
    llim = max(0, roi.left - roi.bgr_width)
    hlim = min(len(chans)-1, roi.right + roi.bgr_width)
    segcounts = counts[llim:hlim]
    maxcounts = max(segcounts)
    ccen = llim + np.where(segcounts==maxcounts)[0]
    ecen = ccen * mca.slope + mca.offset
    bkgcounts = counts[llim] + counts[hlim]
    if maxcounts < 2*bkgcounts:
        mca.init_calib[roiname] = (eknown, ecen, 0.0, ccen, None)
    else:
        model = GaussianModel() + ConstantModel()
        params = model.make_params(amplitude=maxcounts,
                                   sigma=(chans[hlim]-chans[llim])/2.0,
                                   center=ccen-llim, c=0.00)
        params['center'].min = -10
        params['center'].max = hlim - llim + 10
        params['c'].min = -10
        out = model.fit(counts[llim:hlim], params, x=chans[llim:hlim])
        ccen = llim + out.params['center'].value
        ecen = ccen * mca.slope + mca.offset
        fwhm = out.params['fwhm'].value * mca.slope
        mca.init_calib[roiname] = (eknown, ecen, fwhm, ccen, out)

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def create_model_params(x, y):
    """Create the model and parameters."""
    exp_mod = ExponentialModel(prefix='exp_')
    params = exp_mod.guess(y, x=x)

    gauss1 = GaussianModel(prefix='g1_')
    params.update(gauss1.make_params())

    gauss2 = GaussianModel(prefix='g2_')
    params.update(gauss2.make_params())

    params['g1_center'].set(value=105, min=75, max=125)
    params['g1_sigma'].set(value=15, min=3)
    params['g1_amplitude'].set(value=2000, min=10)

    params['g2_center'].set(value=155, min=125, max=175)
    params['g2_sigma'].set(value=15, min=3)
    params['g2_amplitude'].set(value=2000, min=10)

    model = gauss1 + gauss2 + exp_mod
    return model, params

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def get_gaussianmodel(amplitude=1.0, center=5.0, sigma=1.0, noise=0.1):
    # create data to be fitted
    np.random.seed(7392)
    x = np.linspace(-20, 20, 201)
    y = gaussian(x, amplitude, center=center, sigma=sigma)
    y = y + np.random.normal(size=len(x), scale=noise)

    model = GaussianModel()
    params = model.make_params(amplitude=amplitude/5.0,
                               center=center-1.0,
                               sigma=sigma*2.0)
    return x, y, model, params

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def fitTwoGaussians(x,y):
	background  = PolynomialModel(2)
	pars = background.make_params()
	peak1 = GaussianModel(prefix='p1_')
	pars.update( peak1.make_params())
	peak2 = GaussianModel(prefix='p2_')
	pars.update( peak2.make_params())
	# Guess some parameters from data to help the fitting
	span = max(x)-min(x)
	c1Guess = (y[-1]-y[0])/(x[-1]-x[0])
	c0Guess = y[0]-c1Guess*x[0]
	bgGuess = background.func(x=x,c0=c0Guess,c1=c1Guess,c2=0.)
	signalGuess=min(y-bgGuess)
	sigmaGuess = span/30.
	amplitudeGuess = signalGuess*(sigmaGuess*np.sqrt(2.0*np.pi))
	# Fit variables initialization
	
	# pars.add('splitting',0.0001,max=span)
	
	pars['c2'].set(0.,min=-0.000001,max=0.001)
	pars['c1'].set(c1Guess)
	pars['c0'].set(c0Guess)
	pars['p1_center'].set(min(x)+span*0.35,min=min(x),max=max(x))
	pars['p2_center'].set(min(x)+span*0.55,min=min(x),max=max(x))
	# pars['p2_center'].set(min(x)+span*0.65,expr='p1_center+splitting')
	pars['p1_amplitude'].set(amplitudeGuess,max=amplitudeGuess/10000.)
	pars['p2_amplitude'].set(amplitudeGuess,max=amplitudeGuess/10000.)
	pars['p1_sigma'].set(sigmaGuess, min=sigmaGuess/100.,max=sigmaGuess*10000.)
	pars['p2_sigma'].set(sigmaGuess, min=sigmaGuess/100.,max=sigmaGuess*10000.)
	#Add some useful parameters to evaluate
	pars.add('p1_signal', expr='p1_amplitude/(p1_sigma*sqrt(2.0*pi))')
	pars.add('p2_signal', expr='p2_amplitude/(p2_sigma*sqrt(2.0*pi))')
	pars.add('p1_contrast', expr='-p1_amplitude/(p1_sigma*sqrt(2.0*pi)*(c0+c1*p1_center+c2*p1_center**2))')
	pars.add('p2_contrast', expr='-p2_amplitude/(p2_sigma*sqrt(2.0*pi)*(c0+c1*p2_center+c2*p2_center**2))')
	pars.add('splitting',pars['p2_center']-pars['p1_center'],expr='p2_center-p1_center',min=0.00001)
	model = peak1 + peak2 + background
	init = model.eval(pars, x=x)
	out = model.fit(y, pars, x=x)
	# print out.fit_report()
	return init,out

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def test_example_2_Gaussians_1_exp():
  dat = np.loadtxt('NIST_Gauss2.dat')
  x = dat[:, 1]
  y = dat[:, 0]

  exp_mod = ExponentialModel(prefix='exp_')
  pars = exp_mod.guess(y, x=x)

  gauss1  = GaussianModel(prefix='g1_')
  pars.update(gauss1.make_params())

  pars['g1_center'].set(105, min=75, max=125)
  pars['g1_sigma'].set(15, min=3)
  pars['g1_amplitude'].set(2000, min=10)

  gauss2  = GaussianModel(prefix='g2_')

  pars.update(gauss2.make_params())

  pars['g2_center'].set(155, min=125, max=175)
  pars['g2_sigma'].set(15, min=3)
  pars['g2_amplitude'].set(2000, min=10)

  mod = gauss1 + gauss2 + exp_mod


  init = mod.eval(pars, x=x)
  plt.plot(x, y)
  plt.plot(x, init, 'k--')

  out = mod.fit(y, pars, x=x)

  print(out.fit_report(min_correl=0.5))

  plt.plot(x, out.best_fit, 'r-')
  plt.show()

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def lmfit_ngauss_constrains(x,y, params, constrains):
  """
  INPUT:
  x - is the wavelength array
  y - is the normalized flux
  params - is a list/array of initial guess values for the parameters
  		   (this controls the number of gaussians to be fitted
  		   	number of gaussians: len(params)/3 - 3 parameters per Gaussian)
  contrains - the limits of the constrains for the fit of the parameters
  OUTPUT:
  mod - the lmfit model object used for the fit
  out - the lmfit fit object that contains all the results of the fit
  init- array with the initial guess model (usefull to see the initial guess when plotting)
  """

  mods = []
  prefixes = []
  for i in range(0, len(params), 3):
    pref = "g%02i_" % (i/3)
    gauss_i = GaussianModel(prefix=pref)

    if i == 0:
      pars = gauss_i.guess(y, x=x)
    else:
      pars.update(gauss_i.make_params())
    A = params[i]
    limA = constrains[i]
    l_cen = params[i+1]
    limL = constrains[i+1]
    sigma = params[i+2]
    limS = constrains[i+2]

    pars[pref+'amplitude'].set(A, min=limA[0], max=limA[1])
    pars[pref+'center'].set(l_cen, min=limL[0], max=limL[1])
    pars[pref+'sigma'].set(sigma, min=limS[0], max=limS[1])

    mods.append(gauss_i)
    prefixes.append(pref)

  mod = mods[0]

  if len(mods) > 1:
    for m in mods[1:]:
      mod += m
  init = mod.eval(pars, x=x)
  out = mod.fit(y, pars, x=x)

  return mod, out, init

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def test_itercb():
    x = np.linspace(0, 20, 401)
    y = gaussian(x, amplitude=24.56, center=7.6543, sigma=1.23)
    y = y  - .20*x + 3.333 + np.random.normal(scale=0.23,  size=len(x))
    mod = GaussianModel(prefix='peak_') + LinearModel(prefix='bkg_')

    pars = mod.make_params(peak_amplitude=21.0,
                           peak_center=7.0,
                           peak_sigma=2.0,
                           bkg_intercept=2,
                           bkg_slope=0.0)

    out = mod.fit(y, pars, x=x, iter_cb=per_iteration)

    assert(out.nfev == 23)
    assert(out.aborted)
    assert(not out.errorbars)
    assert(not out.success)

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def measure_line_index_recover_spectrum(wave, params, norm=False):
    """ recover the fitted line profile from params

    Parameters
    ----------
    wave: array-like
        the wavelength to which the recovered flux correspond

    params: 5-element tuple
        the 1 to 5 elements are:
        mod_linear_slope
        mod_linear_intercept
        mod_gauss_amplitude
        mod_gauss_center
        mod_gauss_sigma

    norm: bool
        if True, linear model (continuum) is deprecated
        else linear + Gaussian model is used

    """
    from lmfit.models import LinearModel, GaussianModel
    mod_linear = LinearModel(prefix='mod_linear_')
    mod_gauss = GaussianModel(prefix='mod_gauss_')
    par_linear = mod_linear.make_params()
    par_gauss = mod_gauss.make_params()
    par_linear['mod_linear_slope'].value = params[0]
    par_linear['mod_linear_intercept'].value = params[1]
    par_gauss['mod_gauss_amplitude'].value = params[2]
    par_gauss['mod_gauss_center'].value = params[3]
    par_gauss['mod_gauss_sigma'].value = params[4]
    if not norm:
        flux = 1 - mod_gauss.eval(params=par_gauss, x=wave)
    else:
        flux = \
            (1 - mod_gauss.eval(params=par_gauss, x=wave)) * \
            mod_linear.eval(params=par_linear, x=wave)
    return flux

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def lmfit_ngauss_constrains(x,y, params, constrains):
  #params = params[0]
  #constrains = constrains[0]
  mods = []
  prefixes = []
  for i in range(0, len(params), 3):
    pref = "g%02i_" % (i/3)
    gauss_i = GaussianModel(prefix=pref)

    if i == 0:
      pars = gauss_i.guess(y, x=x)
    else:
      pars.update(gauss_i.make_params())

    A = params[i]
    limA = constrains[i]
    l_cen = params[i+1]
    limL = constrains[i+1]
    sigma = params[i+2]
    limS = constrains[i+2]

    pars[pref+'amplitude'].set(A, min=limA[0], max=limA[1])
    pars[pref+'center'].set(l_cen, min=limL[0], max=limL[1])
    pars[pref+'sigma'].set(sigma, min=limS[0], max=limS[1])

    mods.append(gauss_i)
    prefixes.append(pref)

  mod = mods[0]

  if len(mods) > 1:
    for m in mods[1:]:
      mod += m

  init = mod.eval(pars, x=x)
  out = mod.fit(y, pars, x=x)
  return mod, out, init

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def test_reports_created():
    """do a simple Model fit but with all the bells-and-whistles
    and verify that the reports are created
    """
    x = np.linspace(0, 12, 601)
    data = gaussian(x, amplitude=36.4, center=6.70, sigma=0.88)
    data = data + np.random.normal(size=len(x), scale=3.2)
    model = GaussianModel()
    params = model.make_params(amplitude=50, center=5, sigma=2)

    params['amplitude'].min = 0
    params['sigma'].min = 0
    params['sigma'].brute_step = 0.001

    result = model.fit(data, params, x=x)

    report = result.fit_report()
    assert(len(report) > 500)

    html_params = result.params._repr_html_()
    assert(len(html_params) > 500)

    html_report = result._repr_html_()
    assert(len(html_report) > 1000)

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import make_params [as 别名]
def per_iteration(pars, iter, resid, *args, **kws):
    if iter < 3 or iter % 10 == 0:
        out = ['== %i ' % iter]
        for key, val in pars.valuesdict().items():
            out.append('%s=%.3f' % (key, val))
        print ', '.join(out)
        print args, kws


x = linspace(0., 20, 401)
y = gaussian(x, amplitude=24.56, center=7.6543, sigma=1.23)
y = y  - .20*x + 3.333 + random.normal(scale=0.23,  size=len(x))

mod = GaussianModel(prefix='peak_') + LinearModel(prefix='bkg_')

pars = mod.make_params()
pars['peak_amplitude'].value = 3.0
pars['peak_center'].value = 6.0
pars['peak_sigma'].value = 2.0
pars['bkg_intercept'].value = 0.0
pars['bkg_slope'].value = 0.0


out = mod.fit(y, pars, x=x, iter_cb=per_iteration)

pylab.plot(x, y, 'b--')

# print(' Nfev = ', out.nfev)
print( out.fit_report())

pylab.plot(x, out.best_fit, 'k-')