写论文用到了,找了好久,在电磁波理论(葛德彪,魏兵著)找到了相应的代码,发到网上吧,希望大家方便点。
第一个:分层介质反射和透射系数传播矩阵法
clear
epsle_0=8.854187817e-12; %air epsion
mu0=4*pi*1e-7; %air mu0
epsion1=1;
epsion2=2;
epsion3=3;
% % epsion=epsle_0*[1,4.2,epsion1,4.2, epsion2,4.2,epsion3,4.2,1];%9
% % mur=mu0*[1,1,1,1,1,1,1,1,1];
% theta=90*pi/180;
% sigma=[0,0,0,0,0,0,0,0,0];
% sigmam=[0,0,0,0,0,0,0,0,0];
% d=[0,-0.001,-0.001006,-0.001012,-0.002012,-0.003012,-0.003018,-0.004018,-0.005018];
epsion=epsle_0*[1,4,1,4,1,4,1,4,1,4,1,4,1];%9
mur=mu0*[1,1.,1,1.,1,1.,1,1.,1,1.,1,1.,1];
d=[0,-0.01,-0.02,-0.03,-0.04,-0.05,-0.06,-0.07,-0.08,-0.09,-0.1,-0.11,-0.12];
sigma=[0,0,0,0,0,0,0,0,0,0,0,0,0];
sigmam=[0,0,0,0,0,0,0,0,0,0,0,0,0];
theta=30*pi/180;
M=size(epsion);
N=M(1,2)-1;
disp('Please choose the incident wave:TM=-1 or TE=1');
wavetype=input('wavetype is:');
fid1=fopen('RT.dat','wt');
for f=10:10:10000
omega=2*pi*f*1e6;%入射波圆频率(实际频率为循环变量的1e6倍)
k = omega*sqrt(epsion .* mur);
kx=k(:,1)*sin(theta);%入射波矢量在X方向的投影(连续)
kz=sqrt(k.^2-kx.^2);%存放各介质区域的波矢量z分量大小矩阵
if wavetype==1
for n=1:N
p(n)=mur(:,n+1)*kz(:,n)/(mur(:,n)*kz(:,n+1));
end
elseif wavetype==-1
for n=1:N
p(n)=epsion(:,n+1)*kz(:,n)/(epsion(:,n)*kz(:,n+1));
end
end
for n=1:N
r(n)=(1-p(:,n))/(1+p(:,n));
end
VV=eye(2);
for n=1:N-1
v=0.5*(1+p(:,n));
V11=v*exp(-j*kz(:,n+1)*(d(:,n+1)-d(:,n)));
V12=v*r(:,n)*exp(-j*kz(:,n+1)*(d(:,n+1)-d(:,n)));
V21=v*r(:,n)*exp(j*kz(:,n+1)*(d(:,n+1)-d(:,n)));
V22=v*exp(j*kz(:,n+1)*(d(:,n+1)-d(:,n)));
V=[V11,V12;V21,V22];
VV=VV*V;
end
v1=0.5*(1+p(:,N));
Vtn11=v1*exp(j*kz(:,N+1)*d(:,N));
Vtn12=v1*r(:,N)*exp(j*kz(:,N+1)*d(:,N));
Vtn21=v1*r(:,N)*exp(-j*kz(:,N+1)*d(:,N));
Vtn22=v1*exp(-j*kz(:,N+1)*d(:,N));
Vtn=[Vtn11,Vtn12;Vtn21,Vtn22];
R=(-VV(1,2)/VV(1,1))*exp(j*2*kz(:,1)*d(:,1));
T=(-VV(2,1)*VV(1,2)/VV(1,1)+VV(2,2))*exp(-j*(kz(:,N+1)*d(:,N)-kz(:,1)*d(:,1)));
plot(f,abs(R),'*',f,abs(T),'o');
hold on
fprintf(fid1,'\n%8d%+12.5E%+12.5E',f,abs(R),abs(T));
end
legend('反射系数','透射系数');
xlabel('频率 GHZ');
ylabel('反射透射系数');
fclose(fid1);
TE波:
TM波:
下面是课本给出的
还有个分层介质反射系数的连分法
下面展示一些 内联代码片
。
// A code block
var foo = 'bar';
%%%%%%分层反射连数法
clear
epsle_0=8.854187817e-12; %air epsion
mu0=4*pi*1e-7; %air mu0
epsion1=1;
epsion2=2;
epsion3=3;
% % epsion=epsle_0*[1,4.2,epsion1,4.2, epsion2,4.2,epsion3,4.2,1];%9
% % mur=mu0*[1,1,1,1,1,1,1,1,1];
% theta=90*pi/180;
% sigma=[0,0,0,0,0,0,0,0,0];
% sigmam=[0,0,0,0,0,0,0,0,0];
% d=[0,-0.001,-0.001006,-0.001012,-0.002012,-0.003012,-0.003018,-0.004018,-0.005018];
epsion=epsle_0*[1,4,1,4,1,4,1,4,1,4,1,4,1];%9
mur=mu0*[1,1.,1,1.,1,1.,1,1.,1,1.,1,1.,1];
d=[0,-0.01,-0.02,-0.03,-0.04,-0.05,-0.06,-0.07,-0.08,-0.09,-0.1,-0.11,-0.12];
sigma=[0,0,0,0,0,0,0,0,0,0,0,0,0];
sigmam=[0,0,0,0,0,0,0,0,0,0,0,0,0];
theta=30*pi/180;
disp('Please choose the incident wave:TM=-1 or TE=1');
wavetype=input('wavetype is:');
fid1=fopen('RT.dat','wt');
for f=100:50:10000
omega=2*pi*f*1e6;%入射波圆频率(实际频率为循环变量的1e6倍)
k = omega*sqrt(epsion .* mur);
kx=k(:,1)*sin(theta);%入射波矢量在X方向的投影(连续)
kz=sqrt(k.^2-kx.^2);%存放各介质区域的波矢量z分量大小矩阵
M=size(epsion);
N=M(1,2)-1;
if wavetype==1
for n=1:N
p(n)=mur(:,n+1)*kz(:,n)/(mur(:,n)*kz(:,n+1));
end
elseif wavetype==-1
for n=1:N
p(n)=epsion(:,n+1)*kz(:,n)/(epsion(:,n)*kz(:,n+1));
end
end
for n=1:N
r(n)=(1-p(:,n))/(1+p(:,n));
end
R=0;
for n=N:-1:1
A=1/r(:,n);
B=exp(i*2*kz(:,n+1)*d(:,n));
C=exp(i*2*kz(:,n)*d(:,n));
R=A*C+(1-A^2)*B*C/(A*B+R);
end
plot(f,abs(R),'*');
hold on;
VV=eye(2);
for n=1:N-1
v=0.5*(1+p(:,n));
V11=v*exp(-j*kz(:,n+1)*(d(:,n+1)-d(:,n)));
V12=v*r(:,n)*exp(-j*kz(:,n+1)*(d(:,n+1)-d(:,n)));
V21=v*r(:,n)*exp(j*kz(:,n+1)*(d(:,n+1)-d(:,n)));
V22=v*exp(j*kz(:,n+1)*(d(:,n+1)-d(:,n)));
V=[V11,V12;V21,V22];
VV=VV*V;
end
end
% v1=0.5*(1+p(:,N));
% Vtn11=v1*exp(j*kz(:,N+1)*d(:,N));
% Vtn12=v1*r(:,N)*exp(j*kz(:,N+1)*d(:,N));
% Vtn21=v1*r(:,N)*exp(-j*kz(:,N+1)*d(:,N));
% Vtn22=v1*exp(-j*kz(:,N+1)*d(:,N));
% Vtn=[Vtn11,Vtn12;Vtn21,Vtn22];
%
% R=(-VV(1,2)/VV(1,1))*exp(j*2*kz(:,1)*d(:,1));
%
% T=(-VV(2,1)*VV(1,2)/VV(1,1)+VV(2,2))*exp(-j*(kz(:,N+1)*d(:,N)-kz(:,1)*d(:,1)));
%
% plot(f,abs(R),'*',f,abs(T),'o');
% hold on
%
% fprintf(fid1,'\n%8d%+12.5E%+12.5E',f,abs(R),abs(T));
% end
%
% legend('反射系数','透射系数');
% xlabel('频率 GHZ');
% ylabel('反射透射系数');
% fclose(fid1);
TE波
TM波:
课本给的:
课本链接:
通过百度网盘分享的文件:电磁波理论_葛德彪,魏兵著_978-7-03-032006-3.pdf
链接:https://pan.baidu.com/s/1YmZ3uOLFt9j9WiMW2a2_Dg
提取码:g8js