完整程序:

 %题目要求电场入射到xo平面
clear;
clc;
u0=4*pi*1e-7;           %自由空间中的磁导率
e0=1e-9/(36*pi);        %自由空间中的电介质常数
f=1e8;                  %电磁波的频率
w=2*pi*f;               
Ei=5*1.41;              %入射波幅度
R=0.052;                %反射系数
T=0.526;                %透射系数
Er=Ei*R;                %反射波幅值
Et=Ei*T;                %透射波幅值
k1=w*(u0*e0)^0.5;       %真空中（介质1）波数      
k2=k1*2;                %介质2波数            
gifname='E_view_40,0,0.gif';      %取图用
figure(1)         
for t=0:1:200            %为了消除波数与频率之间的数量级带来的影响，时间单位为ns   
    %入射波部分
    %===入射电场部分
    xei=-40:0.1:-1.5;              %入射电场采样范围及采样步长
    m1=zeros(size(xei));
    [yei,zei]=meshgrid(-40:0.1:-1.5,-40:0.1:-1.5);          %生成y,z取样数组网络
    Eyi=1.73*(1.73*yei+zei)+Ei*0.5*cos(w*t*1e-9-k1*0.5*(1.73*yei+zei));    %电场y方向上分量，来源教材公式6.8-22
    %时域分解来源于5.7-14，后面反射，折射亦同的相同
    %其中1.73*(1.73*yei+zei)相当于根号3y尖，为方向系数，是为了与画图函数plot3匹配而设置
    Ezi=-1*(1.73*yei+zei)+Ei*0.5*cos(w*t*1e-9-k1*0.5*(1.73*yei+zei));      %电场z方向上分量
    plot3(m1,Eyi,Ezi,'r','LineWidth',1);                    %绘制入射电场
    hold on
    %反射波部分
    %===反射电场部分
    xer=0:0.1:40;
    m2=zeros(size(xer));
    [yer,zer]=meshgrid(0:0.1:40,0:0.1:40);
    Eyr=1.73*(1.73*yer+zer)+Er*0.5*cos(k1*0.5*(1.73*yer+zer)-w*t*1e-9);
    Ezr=1*(1.73*yer+zer)-Er*0.5*cos(k1*0.5*(1.73*yer+zer)-w*t*1e-9);
    plot3(m2,Eyr,Ezr,'r','LineWidth',1);
    hold on
    %透射波部分
    %===透射电场部分
    xei=0.05:0.1:40;
    m3=zeros(size(xei));
    [yei,zei]=meshgrid(0.05:0.1:40,0.05:0.1:40);
    Eyt=1.73*(1.73*yei+3.61*zei)+Et*0.25*cos(k2*(1.73*yei+3.61*zei)-w*t*1e-9);
    Ezt=-3.61*(1.73*yei+3.61*zei)+Et*0.25*cos(k2*(1.73*yei+3.61*zei)-w*t*1e-9);
    plot3(m3,Eyt,Ezt,'r','LineWidth',1);
    hold on
    view([40,0,0]);               %从x轴方向上观察图形
    a=-40:1:40;                   %下列数行是为了生成xoy平面辅助观察
    b=a;
    [a,b]=meshgrid(a,b)
    c=a*t*0;
    surf(a,b,c);
    shading interp
    set(gca,'ZLim',[-40 40])
    alpha(0.6)
    h1=quiver3(0,0,0,40,0,0);           %下列数行生成矢量箭头指出各方向，辅助观察
    h2=quiver3(0,0,0,0,40,0);
    h3=quiver3(0,0,0,0,0,40);
    h4=quiver3(0,0,0,0,30,17.4);
    h5=quiver3(0,-30,17.4,0,33,-18.9);
    h6=quiver3(0,0,0,0,13.86,-28.84);
    text(36,0,0,'\leftarrow X');        %下列数行生成文字标签，便于观察和说明
    text(0,36,0,'\leftarrow Y');
    text(0,0,36,'\leftarrow Z');
    text(0,27,14,'\leftarrow   θr=60°');
    text(0,-2,2,'\rightarrow   θi=60°');
    text(0,13.36,-26.34,'\leftarrow   θt=30°');
    set(h1,'maxheadsize',0.5);            %控制箭头大小
    set(h2,'maxheadsize',0.5);
    set(h3,'maxheadsize',0.5);
    set(h4,'maxheadsize',0.5);
    set(h5,'maxheadsize',0.5);
    set(h6,'maxheadsize',0.5);
    hold on
    axis([-40,40,-40,40,-40,40]);         %划定图像显示范围
    grid on;
    hold off
    xlabel('x轴')
    ylabel('y轴')
    zlabel('z轴')
    title(['均匀平面波斜投射到不同介质分界面示意图','t=',num2str(t),'ns'],'fontsize',14)
    set(gca,'fontsize',12)
    drawnow
    frame=getframe(1);              %获取当前图像，生成gif图
    im=frame2im(frame);
    [imind,cm]=rgb2ind(im,500);
    if t==0
        imwrite(imind,cm,gifname,'gif');
    else
            imwrite(imind,cm,gifname,'gif','WriteMode','append','DelayTime',0.1);
    end
end