高斯光束传播-matlab

%Gaussian_propagation.m%SimulationofdiffractionofGaussianBeamclear;%GaussianBeam%N:samplingnumberN=input('Numberofsamples(enterfrom100to500)=');L=10*10^-3;Ld=input('wavelengthoflightin[micrometers]=');Ld=Ld*10^-6;ko=(2*pi)/Ld;wo=input('WaistofGaussianBeamin[mm]=');wo=wo*10^-3;z_ray=(ko*wo^2)/2*10^3;sprintf('Rayleighrangeis%f[mm]',z_ray)z_ray=z_ray*10^-3;z=input('Propagationlength(z)in[mm]');z=z*10^-3;%dx:stepsizedx=L/N;forn=1:N+1form=1:N+1%Spaceaxisx(m)=(m-1)*dx-L/2;y(n)=(n-1)*dx-L/2;%GaussianBeaminspacedomainGau(n,m)=exp(-(x(m)^2+y(n)^2)/(wo^2));%FrequencyaxisKx(m)=(2*pi*(m-1))/(N*dx)-((2*pi*(N))/(N*dx))/2;Ky(n)=(2*pi*(n-1))/(N*dx)-((2*pi*(N))/(N*dx))/2;%FreespacetransferfunctionH(n,m)=exp(j/(2*ko)*z*(Kx(m)^2+Ky(n)^2));endend%GaussianBeaminFrequencydomainFGau=fft2(Gau);FGau=fftshift(FGau);%PropagatedGaussianbeaminFrequencydomainFGau_pro=FGau.*H;%PeakamplitudeoftheinitialGaussianbeamPeak_ini=max(max(abs(Gau)));sprintf('Initialpeakamplitudeis%f[mm]',Peak_ini)%PropagatedGaussianbeaminspacedomainGau_pro=ifft2(FGau_pro);Gau_pro=Gau_pro;%PeakamplitudeofthepropagatedGaussianbeamPeak_pro=max(max(abs(Gau_pro)));sprintf('Propagatedpeakamplitudeis%f[mm]',Peak_pro)%CalculatedBeamWidth[NM]=min(abs(x));Gau_pro1=Gau_pro(:,M);[N1M1]=min(abs(abs(Gau_pro1)-abs(exp(-1)*Peak_pro)));Bw=dx*abs(M1-M)*10^3;sprintf('Beamwidth(numerical)is%f[mm]',Bw)%TheoreticalBeamWidthW=(2*z_ray)/ko*(1+(z/z_ray)^2);W=(W^0.5)*10^3;sprintf('Beamwidth(theoretical)is%f[mm]',W)%axisinmmscalex=x*10^3;y=y*10^3;figure(1);mesh(x,y,abs(Gau))title('InitialGaussianBeam')xlabel('x[mm]')ylabel('y[mm]')axis([min(x)max(x)min(y)max(y)01])axissquarefigure(2);mesh(x,y,abs(Gau_pro))title('propagatedGaussianBeam')xlabel('x[mm]')ylabel('y[mm]')axis([min(x)max(x)min(y)max(y)01])axissquare