高斯光束及其MATLAB仿真

高斯光束的定义

在光学中，高斯光束（英语：Gaussian beam）是横向电场以及辐照度分布近似满足高斯函数的电磁波光束。许多激光都近似满足高斯光束的条件，在这种情况中，激光在光谐振腔中以TEM00TEM_{00}TEM00波模（横向基模）传播。当它在满足近衍射极限的镜片中发生折射时，高斯光束会变换成另一种不同参数的高斯光束，因此，高斯光束是激光光学中一种方便、广泛应用的模型。
描述高斯光束的数学函数是亥姆霍兹方程的一个近轴近似解（属于小角近似的一种）。这个解具有高斯函数的形式，代表了光束中电场分量的复振幅。尽管电磁波的传播包括电场和磁场两部分，研究其中任一个场，就足以描述波在传播时的性质。
亥姆霍兹方程的近轴近似解可能不止一个。笛卡尔坐标系下求解可得一类称为厄米-高斯模的解，在柱坐标中求解则得到一类称为拉盖尔-高斯模的解。

基模高斯光束

一般表达式：
E(x,y,z)=w0wexp⁡(−x2+y2w2)exp⁡{i[k(x2+y2)2R−ψ]}exp⁡(ikz)E(x,y,z)=\frac{w_0}{w}\exp\mathrm{(}-\frac{x^2+y^2}{w^2}) \\ \exp \left\{ i\left[ \frac{k(x^2+y^2)}{2R}-\psi \right] \right\} \exp\mathrm{(}ikz) E(x,y,z)=ww0exp(−w2x2+y2)exp{i[2Rk(x2+y2)−ψ]}exp(ikz)
简化表达式为：
E(x,y,z)=11+iz/zRexp⁡[−(x2+y2)/w021+iz/zR]exp⁡(ikz)E(x,y,z)=\frac{1}{1+i{{z}\Bigg/{z_R}}}\exp \left[ -\frac{(x^2+y^2{{)}\Bigg/{w_{0}^{2}}}}{1+i{{z}\Bigg/{z_R}}} \right] \exp\mathrm{(}ikz) E(x,y,z)=1+iz/zR1exp⎣⎢⎢⎢⎢⎡−1+iz/zR(x2+y2)/w02⎦⎥⎥⎥⎥⎤exp(ikz)
w0w_{0}w0是中心束腰半径；k=2π/λk=2\pi/\lambdak=2π/λ是波数；λ\lambdaλ是光束的波长；zR=1/2kw0z_{R}=1/2kw_{0}zR=1/2kw0表示瑞利距离；w=w01+(z/zR)2w=w_{0}\sqrt{1+(z/z_{R})^2}w=w01+(z/zR)2表示光束传播到zzz处的束宽，R=zR(zzR+zRz)R=z_{R}(\frac{z}{z_{R}}+\frac{z_{R}}{z})R=zR(zRz+zzR)表示等相位面的曲率半径，Ψ\varPsiΨ表示相位因子。

拉盖尔-高斯光束

一般表达式：
Epl(x,y,z)=w0w(2rw)lLpl(−r2w2)exp⁡(−r2w2)exp⁡{i[kr22R−(2p+l+1)ψ]}exp⁡(ilφ)exp⁡(ikz)E_{p}^{l}(x,y,z)=\frac{w_0}{w}\left( \sqrt{2}\frac{r}{w} \right) ^lL_{p}^{l}\left( -\frac{r^2}{w^2} \right) \exp \left( -\frac{r^2}{w^2} \right) \\ \exp \left\{ i\left[ k\frac{r^2}{2R}-(2p+l+1)\psi \right] \right\} \\ \exp\mathrm{(}il\varphi \left) \exp \right( ikz) Epl(x,y,z)=ww0(2wr)lLpl(−w2r2)exp(−w2r2)exp{i[k2Rr2−(2p+l+1)ψ]}exp(ilφ)exp(ikz)
简化表达式：
E(x,y,z)=(2rw)lLpl(2r2w2)[1−iz/zR1+(z/zR)2]2p+l11+iz/zRexp⁡[−(x2+y2)/w021+iz/zR]exp⁡(ikz)E(x,y,z)=\left( \sqrt{2}\frac{r}{w} \right) ^lL_{p}^{l}\left( 2\frac{r^2}{w^2} \right) \left[ \frac{1-i{{z}\Bigg/{z_R}}}{\sqrt{1+({{z}\Bigg/{z_R}})^2}} \right] ^{2p+l} \\ \frac{1}{1+i{{z}\Bigg/{z_R}}}\exp \left[ -\frac{(x^2+y^2{{)}\Bigg/{w_{0}^{2}}}}{1+i{{z}\Bigg/{z_R}}} \right] \exp\mathrm{(}ikz) E(x,y,z)=(2wr)lLpl(2w2r2)⎣⎢⎢⎢⎢⎢⎡1+(z/zR)21−iz/zR⎦⎥⎥⎥⎥⎥⎤2p+l1+iz/zR1exp⎣⎢⎢⎢⎢⎡−1+iz/zR(x2+y2)/w02⎦⎥⎥⎥⎥⎤exp(ikz)

厄密-高斯光束

一般表达式：
Emn(x,y,z)=w0wHm(2xw)Hn(2yw)exp⁡(−r2w2)exp⁡{i[kr22R−(2p+l+1)ψ]}exp⁡(ikz)E_{mn}(x,y,z)=\frac{w_0}{w}H_m\left( \sqrt{2}\frac{x}{w} \right) H_n\left( \sqrt{2}\frac{y}{w} \right) \exp \left( -\frac{r^2}{w^2} \right) \\ \exp \left\{ i\left[ k\frac{r^2}{2R}-(2p+l+1)\psi \right] \right\} \exp\mathrm{(}ikz) Emn(x,y,z)=ww0Hm(2wx)Hn(2wy)exp(−w2r2)exp{i[k2Rr2−(2p+l+1)ψ]}exp(ikz)
简化表达式：
Emn(x,y,z)=Hm(2xw)Hn(2yw)[1−iz/zR1+(z/zR)2]2p+l11+iz/zRexp⁡[−(x2+y2)/w021+iz/zR]exp⁡(ikz)E_{mn}(x,y,z)=H_m\left( \sqrt{2}\frac{x}{w} \right) H_n\left( \sqrt{2}\frac{y}{w} \right) \left[ \frac{1-i{{z}\Bigg/{z_R}}}{\sqrt{1+({{z}\Bigg/{z_R}})^2}} \right] ^{2p+l} \\ \frac{1}{1+i{{z}\Bigg/{z_R}}}\exp \left[ -\frac{(x^2+y^2{{)}\Bigg/{w_{0}^{2}}}}{1+i{{z}\Bigg/{z_R}}} \right] \exp\mathrm{(}ikz) Emn(x,y,z)=Hm(2wx)Hn(2wy)⎣⎢⎢⎢⎢⎢⎡1+(z/zR)21−iz/zR⎦⎥⎥⎥⎥⎥⎤2p+l1+iz/zR1exp⎣⎢⎢⎢⎢⎡−1+iz/zR(x2+y2)/w02⎦⎥⎥⎥⎥⎤exp(ikz)

MATLAB仿真程序

采用MATLAB对上述高斯光束进行仿真模拟，部分子函数未编写，可以直接从其他途径搜索得到，这里只是展示主函数和最终的结果。

%============================================
clc;
close all
clear;%% =========基本参数的输入===================
lambda0=6.328e-6;                 %基本光束的输入波长
k0=2*pi/lambda0;                   %波数
w0=1.5*lambda0;                    %高斯光束的中心束腰半径
n_i=1.0;                                   %背景空间折射率
k=k0*n_i;                                 %背景空间的波数
ZR=k*w0^2/2;                         %瑞利距离
z=0;                                          %光束传播距离
w=w0*sqrt(1+(z/ZR)^2);           %传播到z处的束宽
R=ZR*(z/ZR+ZR/z);                 %等相位面曲率半径
Phi=atan(z/ZR);                       %相位因子
m=1;
n=2;                                        %高阶高斯光束的阶数
p=1;                                         %拉盖尔高斯光束的阶数
l=1;                                          %拉盖尔高斯光束的拓扑荷数%% =========坐标参数的设置=====================
Dx=10*lambda0;                        %x轴方向距离
Dy=10*lambda0;                        %y轴方向距离
step=0.01*lambda0;                   %步长
Numx=Dx/step+1;
Numy=Dy/step+1;                      %分割数
x=linspace(-Dx/2,Dx/2,Numx);  %x方向网格
y=linspace(-Dy/2,Dy/2,Numy);   %y方向网格
[x1,y1]=meshgrid(x,y);                 %行成网格数据
[phi,r]=cart2pol(x1,y1);                %转化成极坐标%% ========不同光束的表达式观察平面===============
% 此处采用一般表达式进行仿真计算
%基模高斯光束表达式
E_Gauss=(1/(1+1i*z./ZR)).*exp(-(r.^2./w0^2)./(1+1i*z./ZR)).*exp(1i*k.*z);
I_Gauss=E_Gauss.*conj(E_Gauss);           %计算基模高斯光束的光强
normalizationprocessI_Gauss=I_Gauss/max(max(I_Gauss));   %光强的归一化处理%厄米高斯光束
E_Hermite=Hermite(m,sqrt(2)*x1/w).*Hermite(n,sqrt(2)*y1/w).*((1-1i.*z/ZR)./...sqrt(1+(z/ZR)^2)).^(m+n).*E_Gauss;
I_Hermite=E_Hermite.*conj(E_Hermite);     %计算厄密高斯光束的光强
normalizationprocessI_Hermite=I_Hermite/max(max(I_Hermite)); %光强的归一化处理%拉盖尔高斯光束
E_Lague=power(sqrt(2)*r./w,l).*laguerre(p,l,2*r.^2/w.^2).*power((1-1i.*z/ZR)/...sqrt(1+(z/ZR).^2),2*p+l).*E_Gauss.*exp(1i*l*phi);
I_Lague=E_Lague.*conj(E_Lague);           %计算拉盖尔高斯光束的光强
normalizationprocessI_Lague=I_Lague/max(max(I_Lague));      %光强的归一化处理%% =========观察平面处画图设置==========================
%高斯光束图形设置
figure()
plot(x,normalizationprocessI_Gauss((size(x1,1))/2,:),'linewidth',1.5); %一维图像
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('Normalization Process Density','fontname','times new roman','fontsize',20); %y方向标注figure()
pcolor(x1,y1,normalizationprocessI_Gauss);                             %二维图像
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('\ity(\it\lambda)','fontname','times new roman','fontsize',20); %y方向标注
caxis([0,1]);
shading interp;colormap jet;colorbar; figure()
surf(x1,y1,normalizationprocessI_Gauss);                               %三维图像
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('\ity(\it\lambda)','fontname','times new roman','fontsize',20); %y方向标注
caxis([0,1]);
shading interp;colormap jet;colorbar; %厄密高斯光束图形设置
figure()
plot(x,normalizationprocessI_Hermite((size(x1,1))/2,:),'linewidth',1.5); %一维图像
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('Normalization Process Density','fontname','times new roman','fontsize',20); %y方向标注figure()
pcolor(x1,y1,normalizationprocessI_Hermite);                             %二维图像
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('\ity(\it\lambda)','fontname','times new roman','fontsize',20); %y方向标注
caxis([0,1]);
shading interp;colormap jet;  colorbar; figure()
pcolor(x1,y1,angle(E_Hermite));
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('\ity(\it\lambda)','fontname','times new roman','fontsize',20); %y方向标注
caxis([0,1]);
shading interp;  colormap jet;  colorbar;figure()
surf(x1,y1,normalizationprocessI_Hermite);                               %三维图像
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('\ity(\it\lambda)','fontname','times new roman','fontsize',20); %y方向标注
caxis([0,1]);
shading interp;colormap jet;colorbar; %拉盖尔高斯光束图形设置
figure()
plot(x,normalizationprocessI_Lague((size(x1,1))/2,:),'linewidth',1.5); %一维图像
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('Normalization Process Density','fontname','times new roman','fontsize',20); %y方向标注figure()
pcolor(x1,y1,normalizationprocessI_Lague);                             %二维图像
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('\ity(\it\lambda)','fontname','times new roman','fontsize',20); %y方向标注
caxis([0,1]);
shading interp;colormap jet;colorbar; figure()
pcolor(x1,y1,angle(E_Lague));
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('\ity(\it\lambda)','fontname','times new roman','fontsize',20); %y方向标注
shading interp;colormap jet;figure()
surf(x1,y1,normalizationprocessI_Lague);                               %三维图像
set(gca,'fontname','times new roman','fontsize',15);     %坐标轴设置
xlabel('\itx(\it\lambda)','fontname','times new roman','fontsize',20); %x方向标注
ylabel('\ity(\it\lambda)','fontname','times new roman','fontsize',20); %y方向标注
caxis([0,1]);
shading flat;colormap jet;colorbar;
theAxes=axis;                                        %保存坐标变量%% =======绘制动态图像=========================
%保存11组数据,用来表征传播距离变化时光强的变化
z=0:w0:10*w0;                                     %将传播距离进行离散化
for i=1:size(z,2)%基模高斯光束表达式E_Gauss=(1/(1+1i*z(i)./ZR)).*exp(-(r.^2./w0^2)./(1+1i*z(i)./ZR)).*exp(1i*k.*z(i));I_Gauss=E_Gauss.*conj(E_Gauss);           %计算基模高斯光束的光强normalizationprocessI_Gauss(:,:,i)=I_Gauss/max(max(I_Gauss));   %光强的归一化处理%厄米高斯光束E_Hermite=Hermite(m,sqrt(2)*x1/w).*Hermite(n,sqrt(2)*y1/w).*((1-1i.*z(i)/ZR)./...sqrt(1+(z(i)/ZR)^2)).^(m+n).*E_Gauss;I_Hermite=E_Hermite.*conj(E_Hermite);     %计算厄密高斯光束的光强normalizationprocessI_Hermite(:,:,i)=I_Hermite/max(max(I_Hermite)); %光强的归一化处理%拉盖尔高斯光束E_Lague=power(sqrt(2)*r./w,l).*laguerre(p,l,2*r.^2/w.^2).*power((1-1i.*z(i)/ZR)/...sqrt(1+(z(i)/ZR).^2),2*p+l).*E_Gauss.*exp(1i*l*phi);I_Lague=E_Lague.*conj(E_Lague);           %计算拉盖尔高斯光束的光强normalizationprocessI_Lague(:,:,i)=I_Lague/max(max(I_Lague));      %光强的归一化处理
end