Finite Difference Time Domain with multi materials

This is the finite element time domain for one barrier

I want to add a barrier , the loop was doubled, But I can not see the effect of the second barrier on figure.

I have no error but the results was not resonable

Any advice please

clear;
clc;
%%%%% for the first material
%epsr is the dielectric constant
apsilonr=41.04;
%epsz0 is the permitivity of space
apsilonz0=8.854223e-12;
%sigma is the conductivity
sigma=0.86674;
%KE is the number of cells to be used
KE=400;
dx=zeros(KE,1);
Ex=zeros(KE,1);
Hy=zeros(KE,1);
ix=zeros(KE,1);
gax=ones(KE,1);% amplitude of wave in +ve y axis direction
gbx=zeros(KE,1); %amplitude of wave in -ve y axis direction
%kc is the number of cell in which asinusoidal wave originates
p = zeros(KE,1);
sar = zeros(KE,1);

kc=20;
%f0 the corresponding frequancy of the wave
f0=1800e6;
%c0 the speed of light in avacuum
c0=3e8;
lambdaM1=c0/(sqrt(apsilonr)*f0);
ddx=lambdaM1/10;
dt=ddx./(2*c0);
Ex_low_m2=0;
Ex_low_m1=0;
T =0;
c=0;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%   for the second material
%epsr is the dielectric constant
apsilonr1=52.725;
%epsz0 is the permitivity of space
apsilonz0=8.854223e-12;
%sigma is the conductivity
sigma1=0.942;
%KE is the number of cells to be used

KE=400;
dx1=zeros(KE,1);
Ex1=zeros(KE,1);
Hy1=zeros(KE,1);
ix1=zeros(KE,1);
gax1=ones(KE,1);% amplitude of wave in +ve y axis direction
gbx1=zeros(KE,1); %amplitude of wave in -ve y axis direction
%kc is the number of cell in which asinusoidal wave originates
p1 = zeros(KE,1);
sar1 = zeros(KE,1);

kc=20;
%f0 the corresponding frequancy of the wave
f0=1800e6;
%c0 the speed of light in avacuum
c0=3e8;
lambdaM1=c0/(sqrt(apsilonr1)*f0);
ddx1=lambdaM1/10;
dt=ddx./(2*c0);
Ex1_low_m2=0;
Ex1_low_m1=0;
T =0;
c=0;
%mu=     % permeability 
apsilonrr1=  sigma1/(2*pi*f0) ;   %imaginery factor of apsilonr
%   attenuation factor  alpha(Np/m),  phase factor beta(rad/m) 
%alpha= (2*pi* f0)*sqr((mu*apsilonr/2)*(sqr(1+(apsilonr/apsilonrr)^2)-1))
%beta= (2*pi* f0)*sqr((mu*apsilonr/2)*(sqr(1+(apsilonr/apsilonrr)^2)+1))

%FDTD loop  

%%iteration = 3000; %%% here is the division to different layers



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%mu=     % permeability 
apsilonrr1=  sigma1/(2*pi*f0) ;   %imaginery factor of apsilonr
%   attenuation factor  alpha(Np/m),  phase factor beta(rad/m) 
%alpha= (2*pi* f0)*sqr((mu*apsilonr/2)*(sqr(1+(apsilonr/apsilonrr)^2)-1))
%beta= (2*pi* f0)*sqr((mu*apsilonr/2)*(sqr(1+(apsilonr/apsilonrr)^2)+1))

%FDTD loop  

%%iteration = 3000; %%% here is the division to different layers


for k=120:KE
 gax(k,1)=1/(apsilonr1+(sigma1*dt/apsilonz0));
 gbx(k,1)=sigma1*dt/apsilonz0;
end

for k=100:KE
 gax(k,1)=1/(apsilonr+(sigma*dt/apsilonz0));
 gbx(k,1)=sigma*dt/apsilonz0;
end
for m=100:500:1000
 c=c + 1;
for n=1:m
 T = T + 1;
 for k=2:KE
 dx(k)=dx(k)+.5*(Hy(k-1)-Hy(k));
 end
 source=sin(2*pi*f0*n*dt);
 dx(kc)=dx(kc)+source;
 %calculate the E field of the FDTD
 for k=2:KE
 Ex(k)=gax(k,1)*(dx(k)-ix(k));
 ix(k)=ix(k)+gbx(k,1)*Ex(k);
 end
 %absorbing boundary condition added
 Ex1(1)=Ex1_low_m2;
 Ex1_low_m2=Ex1_low_m1;
 Ex1_low_m1=Ex1(2);
 %calulate the H field of the FDTD
 for k=1:(KE-1)
 Hy(k)=Hy(k)+.5*(Ex(k)-Ex(k+1));
 end
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5
 for k=2:KE
 dx1(k)=dx1(k)+.5*(Hy1(k-1)-Hy(k));
 end
 source1=dx(kc);
 dx1(kc)=dx(kc)+source1;
 %calculate the E field of the FDTD
 for k=2:KE
 Ex1(k)=gax1(k,1)*(dx1(k)-ix(k));
 ix1(k)=ix(k)+gbx(k,1)*Ex1(k);
 end
 %absorbing boundary condition added
 Ex(1)=Ex_low_m2;
 Ex_low_m2=Ex_low_m1;
 Ex_low_m1=Ex(2);
 %calulate the H field of the FDTD
 for k=1:(KE-1)
 Hy1(k)=Hy1(k)+.5*(Ex1(k)-Ex1(k+1));
 end
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 figure(c); plot(Ex,'b'); hold on;
 xa=[120 120];
 ya=[-2 2];
xb=[100 100];
yb=[-2 2];
 plot(xa,ya,xb,yb,'r','LineWidth',2,'LineStyle','--')
 title(['T = ',num2str(n)]);
 ylabel('Propagation of E')
 xlabel('FDTD ')
 axis([1 KE -2 2])
 text(250,-1.15,'Skin')
 text(250,-1.3,'in head')
 text(30,-1.3,'Epsilon=1','color','k')
 hold off
end
end

 for k = 1:19
 p (k,1) = (sigma * (Ex(k,1)))^2 / 2;
 end
 for k = 1:20
 sar(k,1) = p(k,1) / k;
 end
 figure;
 plot(Ex);
 ylabel(' Ex')
 xlabel('x ')
 figure;
 plot(Hy);
ylabel('Hx')
 xlabel('y ')
 figure;
 plot(p);
 ylabel('P')
 xlabel('x ')
 figure;
 plot(sar);
 ylabel('SAR')
 xlabel('x ')