基于Fraunhofer衍射近似的衍射光栅仿真（MATLAB）

I'm trying to generate a diffraction grating pattern using the following amplitude transmission function:

tA=1/2(1-cos(2*pi*x/P))*rect(x/D1)*rect(y/D1);

This method comes from David G.Voelz's Computational Fourier Optics: a MATLAB Tutorial. The book states that when the incident light is a unit-amplitude plane wave, the source field (U1(x1,y1) = tA(x1,y1)). I want to ask: if the incident light is a Gaussian beam, can it be treated as a plane wave? Or do I need to modify the source field (U1)?

Below is the code I learned from the book:

Far-field approximation function #

function[u2,L2]=propFF(u1,L1,lambda,z)
%propagation - Fraunhofer pattern
%assumes uniform sampling
%u1 - source plane field
%L1 - source plane side lenght
%lambda - wavelength
%z - propagation distance
%L2 - observation plane side length
%u2 - observation plane field
[M,N]=size(u1); %get input filed array size
dx1 = L1/M; %source sample interval
k=2*pi/lambda;
L2=lambda*z/dx1;    %obs sidelength
dx2=lambda*z/L1;    %obs sample interval
x2=-L2/2:dx2:L2/2-dx2;  %obs coords
[X2,Y2]=meshgrid(x2,x2);
c=1/(1i*lambda*z)*exp(1i*k/(2*z)*(X2.^2+Y2.^2));
u2=c.*ifftshift(fft2(fftshift(u1)))*dx1^2;
end

Code to generate diffraction grating pattern #

lambda= 0.5e-6; %wavelength
f=0.5;  %propagation distance
P=1e-4; %grating period
D1=1.02e-3; %grating side length
L1=1e-2;    %array side length
M=500;  %number of samples
dx1=L1/M;
x1=-L1/2:dx1:L1/2-dx1;  %source coordinates
[X1,Y1]=meshgrid(x1,x1);
% Grating field and irradiance
u1=1/2*(1-cos(2*pi*X1/P)).*rect(X1/D1).*rect(Y1/D1);
% Fraunhofer pattern
[u2,L2]=propFF(u1,L1,lambda,f);
dx2=L2/M;
x2=-L2/2:dx2:L2/2-dx2; 
y2=x2;
I2=abs(u2).^2;
figure(1)
plot(x2,I2(M/2+1,:));
xlabel('x(m)');
ylabel('Irradiance');

• Can a Gaussian beam be treated as a plane wave?
  No, not generally. Gaussian beams have a non-uniform amplitude profile and curved wavefronts (except at the waist), which are fundamentally different from plane waves. Only if the beam waist is extremely large compared to the grating size and wavelength (making divergence negligible over the propagation distance) could you approximate it as a plane wave—and this is a very specific edge case.

• How to modify the source field (U1)?
  You need to multiply the grating's transmission function (tA) by the Gaussian beam's field distribution at the grating plane. The field of a Gaussian beam at position (z) (relative to its waist) is:
  [
  U_{\text{Gaussian}}(x,y) = \frac{w_0}{w(z)} \exp\left(-\frac{x^2 + y^2}{w(z)2}\right) \exp\left(-ik\left(z + \frac{x^2 + y^2}{2R(z)}\right) + i\zeta(z)\right)
  ]
  Where:

  • (w_0): Beam waist radius
  • (w(z) = w_0\sqrt{1 + \left(\frac{\lambda z}{\pi w_0^2}\right)2}): Beam radius at distance (z) from waist
  • (R(z) = z\left(1 + \left(\frac{\pi w_0^2}{\lambda z}\right)^2\right)): Wavefront curvature radius at (z)
  • (\zeta(z) = \arctan\left(\frac{\lambda z}{\pi w_0^2}\right)): Gouy phase shift

  Modify your source field calculation to include the Gaussian beam (example code):

  % Define Gaussian beam parameters
  w0 = 1e-3; % Adjust to your actual waist radius
  z_waist_to_grating = 0; % Assume grating is at the beam waist plane; adjust if not
  w_z = w0 * sqrt(1 + (lambda*z_waist_to_grating/(pi*w0^2))^2);
  R_z = z_waist_to_grating * (1 + (pi*w0^2/(lambda*z_waist_to_grating))^2);
  zeta_z = atan(lambda*z_waist_to_grating/(pi*w0^2));
  
  % Gaussian beam field at grating plane
  U_gaussian = (w0/w_z) * exp(-(X1.^2 + Y1.^2)/w_z^2) ...
               .* exp(-1i*k*(z_waist_to_grating + (X1.^2 + Y1.^2)/(2*R_z)) + 1i*zeta_z);
  
  % Source field = Gaussian beam * grating transmission
  u1 = U_gaussian .* (1/2*(1-cos(2*pi*X1/P)) .* rect(X1/D1) .* rect(Y1/D1));
  

  Note: If rect isn't available in your environment, implement it as rect(x) = 1 where |x| <= 0.5, else 0.

• Impact on the diffraction pattern
  The Gaussian's amplitude profile will modulate the diffraction orders—orders closer to the beam center will have higher irradiance, while edge orders are suppressed. The curved wavefront may introduce small shifts or shape changes to the pattern, especially if the beam's divergence is significant over the propagation distance (f).

内容的提问来源于stack exchange，提问作者Keith Duong