FDTD.rar_FDTD反射系数_FDTD透射系数_反射系数_反射透射_计算反射系数

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% FDTD程序 --matlab %% %% 图1图2在计算时 动态显示有介质和无介质的情况下的波形 %% 图3第一行是图一中230处的时域波形 %% 第二行是图二中230处的时域波形 %% 第三行是前两行的差 即反射波的波形 %% 图4中显示的是反射波的频谱 %% 图5显示的是入射波的频谱 %% 图6显示的是反射系数.....结果貌似是0.6 %% 图7显示的是图1中315点的时域波形 即透射波形 %% 图8是透射波的频谱 %% 图9是透射波采样点无介质时的频谱 %% 图10是透射系数 %% 缺点:波形撞到后边界时，会形成透射，而不是反射， %% 但1200步的时间步反射波并未影响到抽样点 故此直接无视 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% 给常量赋初值及计算迭代初值E0,H0 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% A=1; % E和H的幅值 k0=175; % 给定的常数K0 W w=40; % a=0.5; % 阿尔法 E=[]; % 定义2个空数组；用于有介质时的循环迭代 H=[]; x=[-600:1:600]; % 高斯脉冲的范围，因为K0=175 所以峰值在0时刻是在175处 Ef=[]; % 定义2个空数组；用于无介质时的循环迭代 Hf=[]; %% 计算E0初值 第十三张幻灯片 E0=A.*exp(log(0.001)*((x-k0)/w).*((x-k0)/w)); %% 计算H0初值 第十三张幻灯片 H05=A.*exp(log(0.001)*((x-k0+0.25)/w).*((x-k0+0.25)/w)); E=E0; %放入E中 H=H05; %放入H中 EK0=E(1,1); %第一个值和最后一个值锁定作为边界 不参与循环 EK0V=E(1,1201); HK0V=H(1,1201); Ef=E0; %与上用处相同 Hf=H05; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% 迭代部分 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% for t=1:1200; %时间步长 1200步 E(t,1)=EK0; %边界值 不参与循环 for k=1:1200; if ((k>=850)&&(k<=909)) %介质空间 er=4; else er=1; end E(t+1,k+1)=E(t,k+1)-(a/er).*(H(t,k+1)-H(t,k)); %迭代公式 H(t+1,k)=H(t,k)-a.*(E(t+1,k+1)-E(t+1,k)); end E(t+1,1201)=EK0V; %边界值 H(t+1,1201)=HK0V; figure(1) %描图 clf plot([250 250],[-0.5 1]); %2条介质线 hold on; plot([309 309],[-0.5 1]); hold on; plot(x,E(t,:),'r'); %波形 hold on; Ef(t,1)=EK0; %无介质的图 for k=1:1200; if ((k>=850)&&(k<=909)) er=1; else er=1; end Ef(t+1,k+1)=Ef(t,k+1)-(a/er).*(Hf(t,k+1)-Hf(t,k)); Hf(t+1,k)=Hf(t,k)-a.*(Ef(t+1,k+1)-Ef(t+1,k)); end Ef(t+1,1201)=EK0V; Hf(t+1,1201)=HK0V; figure(2) %描图 clf plot([250 250],[-0.5 1]); hold on; plot([309 309],[-0.5 1]); hold on; plot(x,Ef(t,:),'r'); hold on; end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% 绘图部分 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% l=[1:1:1201]; fs1=E(:,830); %取有介质第830位置的值 figure(3); % 入射波和反射波的叠加 plot(l,fs1+1,'b'); hold on; fs2=Ef(:,830); %取无介质830位置的值 plot(l,fs2,'g'); %入射波 hold on; fs=fs1-fs2; %二者相减为反射波 plot(l,fs-1,'r'); % 反射波图形 hold on; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% for i=1:50000; %为FFT变换补足50000个数据 fs=[fs;0]; end l1=[1:1:51201]; rfs=fft(fs); ffs=abs(rfs); %fft变化部分 l2=[1:1:1400]; for i=1:1400 lfs(i,1)=ffs(i,1); end figure(4); plot(l2,lfs,'r'); %描频反射波谱图 hold on; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% 求入射波的频谱 描图 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% for i=1:50000; fs2=[fs2;0]; end l1=[1:1:51201]; rfs2=fft(fs2); ffs2=abs(rfs2); l2=[1:1:1400]; for i=1:1400 lfs2(i,1)=ffs2(i,1); end figure(5); plot(l2,lfs2,'g'); hold on; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% 频谱相除 反射系数的图 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% figure(6) last=lfs./lfs2; plot(l2,last,'c'); hold on; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% 透射部分 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ts=E(:,915); %取有介质915点的值观察它的频谱 透射波~ figure(7); plot(l,ts+1,'y'); hold on; ts2=Ef(:,915); %取无介质915点的值观察它的频谱 透射波~ figure(7); plot(l,ts2,'k'); hold on; for i=1:50000; %画有介质透射波频谱 ts=[ts;0]; end l1=[1:1:51201]; rts=fft(ts); fts=abs(rts); l2=[1:1:1400]; for i=1:1400 lts(i,1)=fts(i,1); end figure(8); plot(l2,lts,'y'); hold on; for i=1:50000; %无介质入射波频谱 ts2=[ts2;0]; end l1=[1:1:51201]; rts2=fft(ts2); fts2=abs(rts2); l2=[1:1:1400]; for i=1:1400 lts2(i,1)=fts2(i,1); end figure(9); plot(l2,lts2,'k'); hold on; figure(10) %透射系数图 last2=lts./lts2; plot(l2,last2,'c'); hold on; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%