我有这个代码,旨在将数据放入这里DriveFilelink
read_file函数用于以 struct 化方式提取信息.然而,我在高斯匹配方面遇到了挑战,这从高斯fit image中观察到的差异中可以看出.这个问题似乎源于对某些参数施加的限制,例如固定所有高斯的平均值和四个幅度中的三个.尽管存在这些限制,但它们是必要的,因为它们是基于现有信息的.
有人知道如何解决这个问题吗?为了更贴合.
代码如下:
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import find_peaks
from scipy.optimize import curve_fit
from google.colab import drive
import os
# Mount Google Drive
drive.mount('/content/drive')
# Function to read numeric data from a file
def read_numeric_data(file_path):
with open(file_path, 'r', encoding='latin1') as f:
data = []
live_time = None
real_time = None
for line in f:
if 'LIVE_TIME' in line:
live_time = line.strip()
elif 'REAL_TIME' in line:
real_time = line.strip()
elif all(char.isdigit() or char in ".-+e" for char in line.strip()):
row = [float(x) for x in line.split()]
data.extend(row)
return np.array(data), live_time, real_time
file_path = '/content/drive/MyDrive/ProjetoXRF_Manuel/April2024/20240402_In.mca'
data, _, _ = read_numeric_data(file_path)
a = -0.0188026396003431
b = 0.039549044037714
Data = data
# Función para convolucionar múltiples gaussianas
def multi_gaussian(x, *params):
eps = 1e-10
y = np.zeros_like(x)
amplitude_relations = [1,0.5538, 0.1673 , 0.1673*0.5185 ]
meanvalues= [24.210, 24.002 , 27.276 , 27.238]
amplitude = params[0]
for i in range(4):
sigma = params[i * 3 + 2] # Get the fitted standard deviation for this Gaussian
y += amplitude*amplitude_relations[i]* np.exp(-(x - meanvalues[i])**2 / (2 * sigma**2 + eps))
return y
sigma=[]
area=[]
# Función para graficar el espectro de energía convolucionado
def plot_convolved_spectrum(Data, a, b,i, ax=None):
maxim = np.max(Data)
workingarray = np.squeeze(Data)
# Definir puntos máximos
peaks, _ = find_peaks(workingarray / maxim, height=0.1)
peak_values = workingarray[peaks] / maxim
peak_indices = peaks
# Calcular valores de energía correspondientes a los picos
energy_spectrum = a + b * peak_indices
# Preparar datos para convolución
data = workingarray[:885] / maxim
data_y = data / data.max()
data_x = a + b * np.linspace(0, 885, num=len(data_y))
# Ajustar múltiples gaussianas al espectro de energía
peak_indices2, _ = find_peaks(data_y, height=0.1)
peak_amplitudes = [1,0.5538, 0.1673 , 0.1673*0.5185 ]
peak_means = [24.210, 24.002 , 27.276 , 27.238]
peak_sigmas = [0.1] * 4
params_init = list(zip(peak_amplitudes, peak_means, peak_sigmas))
params_init = np.concatenate(params_init)
# Ajustar curva
params, params_cov = curve_fit(multi_gaussian, data_x, data_y, p0=params_init)
# Obtener una interpolación de alta resolución del ajuste
x_fine = np.linspace(data_x.min(), data_x.max(), num=20000)
y_fit = multi_gaussian(x_fine, *params)
# Data Graphic
ax.scatter(data_x, data_y, color='black', marker='o', s=20, label = 'Data' )
y_data_error =np.sqrt(workingarray[:885])
ax.plot(data_x, data_y + y_data_error/maxim, color='black',linestyle='--')
ax.plot(data_x, data_y - y_data_error/maxim, color='black',linestyle='--')
# Fit Graphic
ax.plot(x_fine, y_fit, label="Fit", linewidth=1.5)
# Extraer desviaciones estándar y amplitudes
sigmas_array = params[2::3]
# Calcular sigma promedio
sigma.append(np.mean(sigmas_array))
# Configuración del gráfico
ax.set_xlabel('Energy (KeV)')
ax.set_ylabel('Normalized Data')
ax.legend()
ax.set_title('Convolved Energy Spectrum')
# Imprimir información
print("Standard deviations:", sigmas_array)
fig, ax = plt.subplots()
plot_convolved_spectrum(Data, a, b,1,ax=ax)
ax.set_xlim(22, 28)
plt.show()
100 我try 过设置界限和细化初始猜测以改善高斯匹配.然而,需要注意的是,使用边界并不等同于修复参数.我的目标是仅维持两个参数的自由度:第一高斯的幅度和所有高斯的标准差.在探索这些调整时,我正在努力在约束和灵活性之间取得平衡,以实现所需的配合准确性.
高斯函数的参数定义如下:A代表幅度,x_0代表平均值,西格玛代表标准差.在分析能谱的背景下,我拥有有关平均值和控制幅度关系的先验知识.因此,我的目标是根据该已知信息限制某些参数(例如A和平均值).具体来说,我打算仅适合第一个幅度,然后利用已知的关系,例如第二个峰值的A2 = 0.1673 * A1.并适应相应的标准偏差.
100
图中第一个峰值的明显奇异性可能表明单高斯.然而,事实并非如此.代表这个能量峰值的高斯实际上是由两个高斯加在一起组成的.这种复杂性之所以出现,是因为该水平的能量排放包括两个不同的值,24.002和24.210.由于实验的分辨率,我们无法仅通过看到情节来区分它们.
出于所有这些原因,我正在试图修复一些参数.
如果有人对这里的数据有问题,那就是: