Python GaussianModel.fit方法代码示例

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [as 别名]
def fit_profile(profile, guess):
    "Fit a profile to a Gaussian + Contant"
    x = np.arange(len(profile))

    model = GaussianModel(missing='drop') + ConstantModel(missing='drop')
    result = model.fit(profile, x=x, verbose=False, **guess)
    return result

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [as 别名]
 def peakFit(self, x, y, model = 'gau', pi = None, di = None):
     ti = time.time()
     if pi:
         NumPeaks = len(pi)
         center = []
         fwhm = []
         amp = []
         numVal = len(x)
         for i in range(NumPeaks):
             pImin = pi[i]-di
             if pImin < 0:
                 pImin = 0
             pImax = pi[i] + di
             if pImax > (numVal-1):
                 pImax = numVal-1
             __y = y[pImin:pImax]
             __x = x[pImin:pImax]
             
             __y = np.power(10,__y/10) #np.array(y)- np.min(y)
             
             mod = GaussianModel()
             pars = mod.guess(__y, x=__x)
             out  = mod.fit(__y, pars, x=__x)
             center.append(out.best_values['center'])
             fwhm.append(out.best_values['sigma']*2.3548)
             amp.append(out.best_values['amplitude'])
         #print 'fit:', time.time()-ti
         return center, fwhm ,amp

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

    area = 1
    cen = 0
    std = 0.2
    x = np.arange(-3, 3, 0.01)
    y = gaussian(x, area, cen, std)

    g = GaussianModel()

    fit_option1 = {'maxfev': 5000, 'xtol': 1e-2}
    result1 = g.fit(y, x=x, amplitude=1, center=0, sigma=0.5, fit_kws=fit_option1)

    fit_option2 = {'maxfev': 5000, 'xtol': 1e-6}
    result2 = g.fit(y, x=x, amplitude=1, center=0, sigma=0.5, fit_kws=fit_option2)

    assert_true(result1.values!=result2.values)
    return

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [as 别名]
def gaussian_fit(x, y, title_name):
	mod = GaussianModel()
	pars = mod.guess(y, x=x)
	out = mod.fit(y, pars, x=x)
	plt.figure()
	plt.plot(x, y)
	plt.plot(x, out.best_fit, 'r-')
	plt.title(title_name)
	print(out.fit_report(min_correl = 0.25))
	print('Center at ' + str(out.best_values['center']) + ' Angstrom')
	plt.show()

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [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 fit [as 别名]
def fluxError(counts, wavelength, error, continuum):
	flux_vector = []
	E_W_vector = []
	cont_avg = np.mean(continuum)
	#plt.close('all')

	for i in range(100):
		#plt.errorbar(wavelength, counts, yerr=error)
		new_counts=[]
		j = 0
		for point in counts:
			new_counts.append(np.random.normal(point, error[j]))
			j = j + 1
		new_counts = np.array(new_counts)
		#So for each N in 1000 a new counts vector is generated randomly
		#Take this data against the wavelength values and fit a gaussian 
		#each time to compute the flux. Append this to a vector and then 
		#find the standard deviation to get the flux error for that emission line
		#Note this is to be encorporated in the fitLines module so that each of the emission 
		#lines is fit in turn. Next step here is to construct the model with lmfit, 
		#guess the initial parameters and then fit the gaussian and take 
		#out.best_values('amplitude') as the flux and store in flux_vector	

		#Now use the lmfit package to perform gaussian fits to the data	
		#Construct the gaussian model
		mod = GaussianModel()

	#Take an initial guess at what the model parameters are 
	#In this case the gaussian model has three parameters, 
	#Which are amplitude, center and sigma
		pars = mod.guess(new_counts, x=wavelength)

	#We know from the redshift what the center of the gaussian is, set this
	#And choose the option not to vary this parameter 
	#Leave the guessed values of the other parameters
		pars['center'].set(value = np.mean(wavelength))


	#Now perform the fit to the data using the set and guessed parameters 
	#And the inverse variance weights form the fits file 
		out  = mod.fit(new_counts, pars, x=wavelength)
		flux = out.best_values['amplitude']
		E_W = out.best_values['amplitude'] / cont_avg
		flux_vector.append(flux)
		E_W_vector.append(E_W)
		#plt.scatter(wavelength, new_counts)
		#plt.plot(wavelength, out.best_fit)
		

	print 'Hello', flux_vector
	#Now return the standard deviation of the flux_vector as the flux error 
	return {'flux_error' : np.std(flux_vector), 'E_W_error' : np.std(E_W_vector)}	

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [as 别名]
def fit_gaussian(y,x):    
    x=array(x)
    y=array(y)
    mod=GaussianModel()
    pars=mod.guess(y,x=x)
    result=mod.fit(y,pars,x=x)
    a=result.params['amplitude'].value
    b=result.params['center'].value
    c=result.params['sigma'].value
    best=result.best_fit
    chsqred=result.redchi
    chisq=result.chisqr
    fwhm=result.params['fwhm'].value
    
    return a,b,c,best,fwhm,chisq,chsqred

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [as 别名]
 def peakFit(self, x, y):
     if len(x) == 0:
         y = self.getdBmSpec()
         y = y[self.__scalePos]
         x = self.__scaledWavelength
     y = np.power(10,y/10)
     mod = GaussianModel()
     pars = mod.guess(y, x=x)
     out  = mod.fit(y, pars, x=x)
     
     print(out.fit_report(min_correl=0.25))
     center = out.best_values['center']
     fwhm = out.best_values['sigma']*2.3548
     
     return center, fwhm#, amp

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [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 fit [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 fit [as 别名]
    dmu = (2. * A * (xx-mu)*s)/(np.pi * fac**2)
    return np.array([ds, dmu, da])


if __name__ == '__main__':
    xs = np.linspace(-4, 4, 100)

    print('**********************************')
    print('***** Test Gaussian **************')
    print('**********************************')
    ys = gaussian(xs, 2.5, 0, 0.5)
    yn = ys + 0.1*np.random.normal(size=len(xs))

    mod = GaussianModel()
    pars = mod.guess(yn, xs)
    out = mod.fit(yn, pars, x=xs)
    out2 = mod.fit(yn, pars, x=xs, fit_kws={'Dfun': dfunc_gaussian,
                                            'col_deriv': 1})
    print('lmfit without dfunc **************')
    print('number of function calls: ', out.nfev)
    print('params', out.best_values)
    print('lmfit with dfunc *****************')
    print('number of function calls: ', out2.nfev)
    print('params', out2.best_values)
    print('\n \n')
    out2.plot(datafmt='.')

    print('**********************************')
    print('***** Test Lorentzian ************')
    print('**********************************')
    ys = lorentzian(xs, 2.5, 0, 0.5)

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [as 别名]
# MODEL = 'loren'
# MODEL = 'voigt'
# gamma_free = True

if MODEL.lower().startswith('g'):
    mod = GaussianModel()
    gamma_free = False
    figname = '../doc/_images/models_peak1.png'
elif MODEL.lower().startswith('l'):
    mod = LorentzianModel()
    gamma_free = False
    figname = '../doc/_images/models_peak2.png'
elif MODEL.lower().startswith('v'):
    mod = VoigtModel()
    figname = '../doc/_images/models_peak3.png'

pars = mod.guess(y, x=x)

if gamma_free:
    pars['gamma'].set(value=0.7, vary=True, expr='')
    figname = '../doc/_images/models_peak4.png'

out = mod.fit(y, pars, x=x)
print(out.fit_report(min_correl=0.25))

plt.plot(x, y, 'b-')
plt.plot(x, out.best_fit, 'r-')
# plt.savefig(figname)
plt.show()
# <end examples/doc_builtinmodels_peakmodels.py>

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [as 别名]
import numpy as np

from lmfit.models import GaussianModel

data = np.loadtxt('model1d_gauss.dat')
x = data[:, 0]
y = data[:, 1]

y[44] = np.nan
y[65] = np.nan

# nan_policy = 'raise'
# nan_policy = 'propagate'
nan_policy = 'omit'

gmodel = GaussianModel()
result = gmodel.fit(y, x=x, amplitude=5, center=6, sigma=1,
                    nan_policy=nan_policy)

print(result.fit_report())

# make sure nans are removed for plotting:
x_ = x[np.where(np.isfinite(y))]
y_ = y[np.where(np.isfinite(y))]

plt.plot(x_, y_, 'bo')
plt.plot(x_, result.init_fit, 'k--')
plt.plot(x_, result.best_fit, 'r-')
plt.show()
# <end examples/doc_model_with_nan_policy.py>

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [as 别名]
    def _gen(self):
        # We checked in the __init__ that this import works.
        from lmfit.models import GaussianModel, LinearModel
        # For thread safety (paranoia) make copies of stuff
        dets = self.detectors
        target_field = self.target_field
        motor = self.motor
        initial_center = self.initial_center
        initial_width = self.initial_width
        tol = self.tolerance
        min_cen = self.min_cen
        max_cen = self.max_cen
        seen_x = deque()
        seen_y = deque()
        for x in np.linspace(initial_center - self.RANGE * initial_width,
                             initial_center + self.RANGE * initial_width,
                             self.NUM_SAMPLES, endpoint=True):
            yield Msg('set', motor, x)
            yield Msg('create')
            ret_mot = yield Msg('read', motor)
            key, = ret_mot.keys()
            seen_x.append(ret_mot[key]['value'])
            for det in dets:
                yield Msg('trigger', det, block_group='B')
            for det in dets:
                yield Msg('wait', None, 'B')
            for det in dets:
                ret_det = yield Msg('read', det)
                if target_field in ret_det:
                    seen_y.append(ret_det[target_field]['value'])
            yield Msg('save')

        model = GaussianModel() + LinearModel()
        guesses = {'amplitude': np.max(seen_y),
                'center': initial_center,
                'sigma': initial_width,
                'slope': 0, 'intercept': 0}
        while True:
            x = np.asarray(seen_x)
            y = np.asarray(seen_y)
            res = model.fit(y, x=x, **guesses)
            old_guess = guesses
            guesses = res.values
            if np.abs(old_guess['center'] - guesses['center']) < tol:
                break
            next_cen = np.clip(guesses['center'] +
                            np.random.randn(1) * guesses['sigma'],
                            min_cen, max_cen)
            yield Msg('set', motor, next_cen)
            yield Msg('create')
            ret_mot = yield Msg('read', motor)
            key, = ret_mot.keys()
            seen_x.append(ret_mot[key]['value'])
            for det in dets:
                yield Msg('trigger', det, block_group='B')
            for det in dets:
                yield Msg('wait', None, 'B')
            for det in dets:
                ret_det = yield Msg('read', det)
                if target_field in ret_det:
                    seen_y.append(ret_det[target_field]['value'])
            yield Msg('save')

        yield Msg('set', motor, np.clip(guesses['center'], min_cen, max_cen))

        if self.output_mutable is not None:
            self.output_mutable.update(guesses)
            self.output_mutable['x'] = np.array(seen_x)
            self.output_mutable['y'] = np.array(seen_y)
            self.output_mutable['model'] = res

# 需要导入模块: from lmfit.models import GaussianModel [as 别名]
# 或者: from lmfit.models.GaussianModel import fit [as 别名]
            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-')
pylab.show()

#<end examples/doc_with_itercb.py>