微波EDA网,见证研发工程师的成长!
首页 > 研发问答 > 微波和射频技术 > 电磁仿真讨论 > fdtd s parameter

fdtd s parameter

时间:03-31 整理:3721RD 点击:
can anyone please help me how to extract the s-parameter from the fdtd simulation result.if any one have matlab coding for that please upload it.

this may be helpful, this is 3d
s parameters are found via fourier transform
------------------------
mr=1,er=1
m=mr*4*pi*10^-7
e=er*8.854*10^-12
h=0.075*10^-2
dx=h;dz=h;
c=1/(m*e)^0.5
dt=0.5*h/c
T=20*dt
a=1.5*10^-2
Ex=zeros(370);
Ez=zeros(370);
Hy=zeros(370);
c1=350;c2=350;
AB1=(c*dt-dz)/(c*dt+dz);
AB2=(c*c*e*dz*dt)/(2*dx*(c*dt+dz));
volt_history1=zeros(1,3000);
volt_history2=zeros(1,3000);
Hy_eskilk1=zeros(20,3000);
Hy_eskilk2=zeros(20,3000);
Hy_eskison1=zeros(20,3000);
Hy_eskison2=zeros(20,3000);
son=zeros(1,511);
son1=zeros(1,511);
fftvoltson=zeros(1,4096);
fftvolthis=zeros(1,4096);
bolum=zeros(1,4096);
bolum1=zeros(1,4096);
bolumson=zeros(1,103);
for j=1:3000

if j<=6*T/dt;
Ex(c2:c2+19,1)=exp(-(j*dt-3*T)^2/T^2);
%elseif j>6*T/dt && j <=6*T/dt + 10
%Ex(c2:c2+19,1)=0;
%elseif j >6*T/dt + 10
%for b=c2:c2+19
% Hy(b,1)=Hy_eskilk2(b-c2+1)-AB1*(Hy_eskilk1(b-c2+1)-Hy(b,2))+AB2*(Ez(b+1,2)+Ez(b+1,1)-Ez(b,2)-Ez(b,1));
%end
end
%abc buraya
for x=1:20
Hy_eskilk1(x,j)=Hy(c2+x-1,1);
Hy_eskilk2(x,j)=Hy(c2+x-1,2);
Hy_eskison1(x,j)=Hy(1,c2+x-1);
Hy_eskison2(x,j)=Hy(2,c2+x-1);
end
for k=1:370

for i=1:369

if k>1
if k==350 && i<350
Ex(i,k)=0;
elseif k==370
Ex(i,k)=0;
else
Ex(i,k)=Ex(i,k)-(dt/(dz*e))*(Hy(i,k)-Hy(i,k-1));
end

else
Ex(i,k)=Ex(i,k)-(dt/(dz*e))*Hy(i,k);
end

if k<370
if k==1 && i>349 && j>6*T/dt + 10
Hy(i,1)=Hy_eskilk2(i-c2+1,j-1)-AB1*(Hy(i,1)-Hy(i,2))+AB2*(Ez(i+1,2)+Ez(i+1,1)-Ez(i,2)-Ez(i,1));
%elseif i==1 && k>349
% Hy(1,i)=Hy_eskison2(k+1-c1)-AB1*(Hy_eskison1(k+1-c1)-Hy(2,i))-AB2*(Ex(2,i+1)+Ex(1,i+1)-Ex(2,i)-Ex(1,i));
elseif i==1 && k>349 && j>3
Hy(1,k)=Hy_eskison2(k+1-c1,j-1)-AB1*(Hy(1,k)-Hy(2,k))-AB2*(Ex(2,k+1)+Ex(1,k+1)-Ex(2,k)-Ex(1,k));
elseif i~=1
Hy(i,k)=Hy(i,k)-(dt/(dz*m))*(Ex(i,k+1)-Ex(i,k))+(dt/(m*dx))*(Ez(i+1,k)-Ez(i,k));
end
if i==349 && k<350
Ez(i+1,k)=0;
elseif i==369
Ez(i+1,k)=0;
elseif i>1
Ez(i,k)=Ez(i,k)+(dt/(dx*e))*(Hy(i,k)-Hy(i-1,k));
else
Ez(i,k)=Ez(i,k)+(dt/(dx*e))*Hy(i,k);
end
end
end
end

voltage1=0;
for x=1:20
voltage1=voltage1+Ex(c2+x-1,c1-15);
end
volt_history1(j)=voltage1*dx;
voltage2=0;
for x=1:20
voltage2=voltage2+Ez(c1-15,c2+x-1);
end
volt_history2(j)=voltage2*dx;
%for a=c1:c1+19
% Hy(1,a)=Hy_eskison2(a+1-c1)-AB1*(Hy_eskison1(a+1-c1)-Hy(2,a))-AB2*(Ex(2,a+1)+Ex(1,a+1)-Ex(2,a)-Ex(1,a));
%end
end
%set(gcf,'renderer','zbuffer');
%mesh(Ex);
%for i=1700:2100
% volt_history1(i)=volt_history1(i);
% volt_history2(i)=volt_history2(i);
%end
%subplot(2,1,1);
%plot(volt_history1)
%subplot(2,1,2);
%plot(volt_history2)
%volt1=fft(volt_history1,4096);
%plot(abs(volt1))
%-----------------------------------------------------------------------
Ex1=zeros(20,2900);
Ez1=zeros(20,2900);
Hy1=zeros(20,2900);
c1=350;c2=350;
AB1=(c*dt-dz)/(c*dt+dz);
AB2=(c*c*e*dz*dt)/(2*dx*(c*dt+dz));
volt_history11=zeros(1,3000);

for j=1:3000

if j<=6*T/dt;
Ex1(1:19,1)=exp(-(j*dt-3*T)^2/T^2);
end

for k=1:2000

for i=1:19
if k>1
Ex1(i,k)=Ex1(i,k)-(dt/(dz*e))*(Hy1(i,k)-Hy1(i,k-1));
else
Ex1(i,k)=Ex1(i,k)-(dt/(dz*e))*Hy1(i,k);
end

Hy1(i,k)=Hy1(i,k)-(dt/(dz*m))*(Ex1(i,k+1)-Ex1(i,k))+(dt/(m*dx))*(Ez1(i+1,k)-Ez1(i,k));
if i==1 || i==20
Ez1(i,k)=0;
else
Ez1(i,k)=Ez1(i,k)+(dt/(dx*e))*(Hy1(i,k)-Hy1(i-1,k));
end

end
end
voltage11=0;
for x=1:19
voltage11=voltage11+Ex1(x,c1-15);
end
volt_history11(j)=voltage11*dx;
end
%plot(volt_history11)
%set(gcf,'renderer','zbuffer');
%mesh(Ez);
for i=1:3000
voltson(i)=volt_history1(i)-volt_history11(i);
end
%plot(volt_history1);
fftvoltson=fft(voltson,4096);
fftvolthis=fft(volt_history11,4096);
fftvolthis1=fft(volt_history2,4096);
for i=1:4096
bolum(i)=fftvoltson(i)/fftvolthis(i);
end
%plot(abs(bolum));
%plot(fftvolthis)
for i=1:4096
bolum1(i)=fftvolthis1(i)/fftvolthis(i);
end


plot(abs(bolumson));legend('|S11|^2+|S12|^2');
xlabel('frequency(GHz)');
set(gca,'XTick',[0 52 103])
set(gca,'XTickLabel',{'0';'10';'20'})
set(gca,'YLim',[0.3 1.1])
axis tight

Fairly easy :
with oAB being the FDTD reponse at the port A to the extication coming from B.
S11=FFT(o11,N)./FFT(i1,N);
S21=FFT(o21,N)./FFT(i1,N);
..
S34=FFT(o34,N)./FFT(i4,N);

and so on.
choose N with respect to the frequency resolution you want and the time sampling you have, p.ex N=T/timestep

Regards,
David

It's easy, Supaswing is right!

hallo supaswing,
can you please again explain how this statement is apply in the code posted by cakalhunter,to determine the s-parameters.

thanks friendds these codes are very useful

Copyright © 2017-2020 微波EDA网 版权所有

网站地图

Top