matlab_对光纤中光的传播特性进行仿真（包括单模和多模光纤）

clc,clear,close all V = 5.0000; %设置LP11模的归一化工作频率 V 的取值 U = 3.1527; %设置LP11模的归一化横向相位参数 U 的取值 W = 3.8808; %设置LP11模的归一化横向衰减参数 W 的取值 %网格设定及计算 Npoint = 201; R1 = linspace(0,1,Npoint); R2 = linspace(1,5,Npoint); Theta1 = linspace(0,2*pi,Npoint); Theta2 = linspace(0,2*pi,Npoint); E1 = zeros(Npoint,Npoint); E2 = zeros(Npoint,Npoint); I1 = E1; I2= E2; %光纤纤芯中的光场分布 for i = 1:Npoint for j = 1:Npoint E1(i,j) = besselj(1,U*R1(i))*cos(Theta1(j)); I1(i,j) = E1(i,j).^2; end end %光纤包层中的光场分布 for i = 1:Npoint for j = 1:Npoint E2(i,j) = besselj(1,U).*besselk(1,W.*R2(i))./ ... besselk(1,W)*cos(Theta1(j)); I2(i,j) = E2(i,j).^2; end end %归一化处理 E2 = E2/max(max(E1)); E1 = E1/max(max(E1)); I2 = I2/max(max(I1)); I1 = I1/max(max(I1)); %极坐标变换为柱坐标，以便于作图 [Theta1 R1] = meshgrid(Theta1, R1); [Theta2 R2] = meshgrid(Theta2, R2); [X1 Y1] = pol2cart(Theta1,R1); [X2 Y2] = pol2cart(Theta2,R2); %三维作图并标注 figure,surf(X1,Y1,E1) shading interp hold on surf(X2,Y2,E2) shading interp xlabel('x'),ylabel('y'),zlabel('E') title(['LP_{11} Mode, V = ' num2str(V)]) figure,surf(X1,Y1,I1) shading interp hold on surf(X2,Y2,I2) shading interp xlabel('x'),ylabel('y'),zlabel('I') title(['LP_{11} Mode, V = ' num2str(V)]) %生成动画 numFrames = 100; mov = moviein(numFrames); figure('Renderer','zbuffer'); %打开新图形窗口 surf(X1,Y1,E1) shading interp xlabel('x'),ylabel('y'),zlabel('I') title(['LP_{11} Mode, V = ' num2str(V)]) % axis tight set(gca,'NextPlot','replaceChildren'); for k = 1:numFrames scale = cos(2*pi*k/numFrames); surf(X1,Y1,E1*scale) shading interp xlabel('x'),ylabel('y'),zlabel('I') title(['LP_{11} Mode, V = ' num2str(V)]) axis([-inf,inf,-inf,inf,-1,1]) view([-38,12]) hold on surf(X2,Y2,E2*scale) shading interp axis([-inf,inf,-inf,inf,-1,1]) mov(k) = getframe; hold off end movie(mov) %播放动画 %将动画存储为 GIF 文件并播放 animated(1,1,1,numFrames) = 0; for k=1:numFrames if k == 1 [animated, cmap] = rgb2ind(mov(k).cdata,256,'nodither'); else animated(:,:,1,k) = rgb2ind(mov(k).cdata,cmap,'nodither'); end end filename = 'LP11.gif'; imwrite(animated,cmap,filename,'DelayTime',0.1,'LoopCount', inf); web(filename)