import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

modelová funkce pro difrakci na dvojštěrbině

def model_fce(alpha,d,b,alpha0):
    wavelength=550 #vlnová délka světla (nm)
    return (np.sinc(d*np.pi/wavelength*np.sin(alpha-alpha0))*np.cos(np.pi*b/wavelength*np.sin(alpha-alpha0)))**2
    

simulace spektra difrakce na šterbině

d=500 #šířka štěrbiny (nm)
b=6000 #vzdálenost štěrbin od sebe (nm)
alpha0=-5*np.pi/180 #poloha maxima (deg)
alpha=np.linspace(-50,50,1000)
alpha=alpha*np.pi/180 #konverze deg. na rad. 
y=model_fce(alpha,d,b,alpha0)+np.random.normal(0,0.005,np.size(alpha))

nakreslení simulovaného spektra

plt.scatter(alpha*180/np.pi,y)
plt.xlabel(r"$\alpha $ (deg.)")
plt.ylabel('intenzita')
plt.show()

fit nelineární metodou nejmenších čtverců

#počáteční odhady
d_initial=500
b_initial=5000
alpha0_initial=-2*np.pi/180
par,cov=curve_fit(model_fce,alpha,y,p0=[d_initial,b_initial,alpha0_initial],bounds=([0,0,-15*np.pi/180],[10000,10000,15*np.pi/180]))
d_fit=par[0]
b_fit=par[1]
alpha0_fit=par[2]

Nakreslení výsledků

plt.scatter(alpha*180/np.pi,y,label='simulated data')
plt.plot(alpha*180/np.pi,model_fce(alpha,d_fit,b_fit,alpha0_fit),c='red',label='fit')
plt.xlabel(r"$\alpha$ (deg.)")
plt.ylabel('intenzita')
plt.legend()
plt.show()

print(f"šířka šterbiny     d = ({d_fit:6.1f} +/-  {np.sqrt(cov[0,0]):3.1f})  nm")
print(f"vzdálenost šterbin b = ({b_fit:6.1f} +/-  {np.sqrt(cov[1,1]):3.1f})  nm")
print(f"poloha maxima alpha0 = ({alpha0_fit*180/np.pi:6.3f} +/-  {np.sqrt(cov[2,2])*180/np.pi:4.3f}) stupňů")

šířka šterbiny     d = ( 499.8 +/-  0.3)  nm
vzdálenost šterbin b = (6000.4 +/-  1.1)  nm
poloha maxima alpha0 = (-4.998 +/-  0.001) stupňů