您如何在Python中使用高斯函数拟合测得的发射线数据？ （原子光谱）

Question

对于一个物理实验室项目，我正在测量来自各种元素的各种发射线。高强度峰出现在某些波长处。我的目标是在Python中拟合高斯函数，以便找到强度在哪个波长处达到峰值。

我已经尝试过使用scipy.stats库中的norm函数。下面是代码和生成的图形。

import numpy as np
from scipy.stats import norm
import matplotlib.pyplot as plt

mean, std = norm.fit(he3888_1[:,0])
plt.plot(he3888_1[:,0], he3888_1[:,1], color='r')
x = np.linspace(min(he3888_1[:,0]), max(he3888_1[:,0]), 100)
y = norm.pdf(x, mean, std)
plt.plot(x, y)
plt.xlabel("Wavelength (Angstroms)")
plt.ylabel("Intensity")
plt.show()

这可能是因为在此之前相对较长的时间内强度较低吗？

Answer 1

Lmfit对于您的情况似乎是一个不错的选择。下面的代码模拟添加了线性背景的高斯峰，并说明如何使用lmfit提取参数。后者具有许多其他易于集成的内置模型（Lorentzian，Voight等）。

import numpy as np
from lmfit.models import Model, LinearModel
from lmfit.models import GaussianModel, LorentzianModel
import matplotlib.pyplot as plt

def generate_gaussian(amp, mu, sigma_sq, slope=0, const=0):
    x = np.linspace(mu-10*sigma_sq, mu+10*sigma_sq, num=200)
    y_gauss = (amp/np.sqrt(2*np.pi*sigma_sq))*np.exp(-0.5*(x-mu)**2/sigma_sq)
    y_linear = slope*x + const
    y = y_gauss + y_linear
    return x, y

# Gaussiand peak generation
amplitude = 6
center = 3884
variance = 4
slope = 0
intercept = 0.05
x, y = generate_gaussian(amplitude, center, variance, slope, intercept)


#Create a lmfit model: Gaussian peak + linear background
gaussian = GaussianModel()
background = LinearModel()
model = gaussian + background

#Find what model parameters you need to specify
print('parameter names: {}'.format(model.param_names))
print('independent variables: {}'.format(model.independent_vars))

#Model fit
result = model.fit(y, x=x, amplitude=3, center=3880,
                   sigma=3, slope=0, intercept=0.1)
y_fit = result.best_fit #the simulated intensity

result.best_values #the extracted peak parameters

# Comparison of the fitted spectrum with the original one
plt.plot(x, y, label='model spectrum')
plt.plot(x, y_fit, label='fitted spectrum')
plt.xlabel('wavelength, (angstroms')
plt.ylabel('intensity')
plt.legend()

输出：

parameter names: ['amplitude', 'center', 'sigma', 'slope', 'intercept']

independent variables: ['x']

result.best_values
Out[139]: 
{'slope': 2.261379140543626e-13,
 'intercept': 0.04999999912168238,
 'amplitude': 6.000000000000174,
 'center': 3883.9999999999977,
 'sigma': 2.0000000000013993}