Распространение оптического поля с помощью метода углового спектра и проверки ⇐ Python

[Anonymous]

Я написал код для распространения гауссовского луча в свободном пространстве методом углового спектра. Ширина луча поле с теоретическими формулами. Ниже приведен мой код, я не понимаю, что в нем не так, ширина луча не соответствует теоретической ширине луча, я бы очень признателен за какое -либо руководство. Спасибо < /p>
Вот мой код: < /p>
wvl = 1550e-9 # wavelength
wz = 0.05 # beam waist
gdim = 1 # spatial extent of the grid;
resol = 512 # grid resolution
I = 1 # intensity amplitude
zr = (np.pi * wz**2)/wvl # rayleigh range
k = 2*np.pi/wvl

# Theoretical beam width at reciever (Lz km) # Gaussian
wzt = lambda Lz, wg, lmda: wg * np.sqrt(1 + ((Lz * lmda)/(np.pi * wg**2))**2)

wz_ = wzt(0, wz, wvl)
z = zr

# define grid
dx = np.sqrt((wvl*z)/resol) #(wvl * np.sqrt(z**2 + gdim**2))/gdim
gdim = dx * resol
x, y = np.meshgrid(np.arange(-gdim/2, gdim/2, dx),
np.arange(-gdim/2, gdim/2, dx))

r = np.sqrt(x**2 + y**2)

def gaussian(wz, r, I):
"""
wz: Beam width at z=0;
r : Radial coordinates
I : Intensity distribution
"""
Fin = I*np.exp(-2*(r/wz)**2)
return Fin
gbeam = gaussian(wz, r, I)
def ft2(g, delta=None):
"""
Computes the 2D Fourier transform of g with scaling.

Parameters:
g : numpy.ndarray
Input 2D array to transform.
delta : float
Spacing in the spatial domain.

Returns:
numpy.ndarray
The scaled 2D Fourier transform of g.
"""
G = np.fft.fftshift(np.fft.fft2(np.fft.ifftshift(g)))
return G

def ang_spec_prop(e, z, resol, dx, wvl):
k = (2*np.pi/wvl)

fxx = np.fft.fftshift(np.fft.fftfreq(resol, dx))
fxx, fyy = np.meshgrid(fxx, fxx)

alfa = k ** 2 - (4 * (np.pi ** 2)) * (fxx ** 2 + fyy ** 2)

tmp = np.sqrt(np.abs(alfa))
kz = np.where(alfa >= 0, tmp, 1j*tmp)

h1 = np.exp(z * kz * 1j)
e_ = ft2(ft2(e)*h1)
return e_
prop_as = ang_spec_prop(gbeam, z, gbeam.shape[0], dx, wvl)

# Calculate the beam width
def anacal_gbw(ub, r):
# Find the intensity where it falls to 1/e^2 of its maximum
ub = np.abs(ub)**2
r_flat = r.flatten()
imax = ub.max()/(np.exp(1)**2)
in_flat = ub.flatten()
wg_calc = r_flat[np.abs(in_flat-imax).argmin()]
return wg_calc
ana_w = anacal_gbw(prop_as, r)
the_w = wzt(z, wz, wvl)
< /code>
С этими параметрами, the_w (теоретическая ширина луча) = 0,0707; но ana_w (аналитическая ширина луча) = 0,0792


Подробнее здесь: https://stackoverflow.com/questions/793 ... rification