%题目要求电场入射到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
