Подгонка свертки ферми-функции с гауссовой функцией к реальным данным с питоном - PullRequest
Новая
       55

Подгонка свертки ферми-функции с гауссовой функцией к реальным данным с питоном

0 голосов
/ 06 февраля 2019

Я недавно имел дело с данными фотоэмиссии и хочу оценить инструментальное разрешение путем подгонки свертки функции ферми с гауссианой к реальным данным.На данный момент я пренебрегаю тем, что я также должен включить фон.Я пытался работать над душой, основываясь на том, что нашел в Link1

К сожалению, я застрял и прошу вас о помощи.Я приведу небольшой пример, в котором были смоделированы данные, которые (судя по глазу) на самом деле являются хорошим представлением фактического экспета данных, в котором отсутствует фон: 

import array as ar
import math
from math import e, sqrt, pi
import pylab as plt
from numpy import array, linspace, arange, zeros, ceil, amax, amin, argmax, argmin, abs
from numpy import polyfit, polyval, seterr, trunc, mean
from numpy.linalg import norm
from scipy.interpolate import interp1d
import scipy.optimize
import numpy as np
import matplotlib.pyplot as plt

kb = 8.6173e-5 # Boltzmann k in eV/K


def f(e,scale,T=300,muf=0):
    return scale*(1.0 / (np.exp((-e - muf)/(kb*T)) + 1)) # -e due to binding energy scale

def gaussian(x, mu, sig,scale):
    return scale*(np.exp(-np.power(-x - mu, 2.) / (2 * np.power(sig, 2.))))  # -x due to binding energy scale




Points=10000 # If i change the points the "modelt data" doesn not look good anymore?
espan = np.linspace(-10, 10, Points)


# simulate the fermi edge data such that is close to the read real data
# note that it actually is a convoluiton of a fermi and a gaussian function.
Sim_data = np.random.normal((8.07+np.convolve(f(espan,scale=1.05), gaussian(espan, mu=0, sig=0.3,scale=1.05),mode="same")),10)

# min max norm not ( not used for now )
# Sim_data_norm = ((Sim_data - Sim_data.min()) / (Sim_data.max() - Sim_data.min()))


# convolution
# I know the actual data was meassured at T=300K so i preset this value in the fermi function
# The position of the fermi function and the gaussian should be 0 so mu and muf are set to 0
# i guess x_values would have to be replaced for the x axis of my real data later?
def Gaus_fermi_convolution(size,sigC,scaleC,xC=espan):
    return size* np.convolve(gaussian(x=xC,sig=sigC,scale=scaleC,mu=0), f(e=espan,scale=scaleC,T=300,muf=0), mode="same") 

# try to fit
lowerBounds = [0.0, 0.0,0.0] # A, B, C ,D lower bounds
upperBounds = [440.0, 0.6,440.0] # A, B, C, D upper bounds

popt, pcov = scipy.optimize.curve_fit(Gaus_fermi_convolution, espan, Sim_data,bounds=[lowerBounds, upperBounds])
Opt_size,Opt_sigC,Opt_scaleC = popt

# Print out the coefficients determined by the optimization
print(Opt_size,Opt_sigC,Opt_scaleC)

# Plotting 
plt.plot(espan,Sim_data)
#The result is clearly wrong
plt.plot(espan,Gaus_fermi_convolution(size=Opt_size,sigC=Opt_sigC,scaleC=Opt_scaleC,xC=espan))

plt.xlim([2, -2])
plt.xlabel( r'$Binding \: Energy  \: [eV]$')
plt.ylabel(r'$Intensity\: [arb. u.] $')
plt.show()
...