%multi-meida perpendicular input%電場對理想介質的垂直入射；
%copyright 屬于研學論壇bbs.matwav.com和plasma
clear all;
k1=1;
k2=2;
k3=4;

ei10=1;
mu0=4*pi*1.e-7;
epsilon0=1.e-9/(36*pi);
epsilonr1=2;%二區的相對介電常數
epsilonr2=8;%三區的相對介電常數
%epsilonr2大則inta2小，inta2=sqrt(mu0/epsilon2);
%并且導致intad變化，同時改變了3區的分布情形。

d=-10;
%注意！！！！！！！！！！！
%當d＝npi/k2并且inta3=inta0的時候反射為零。


inta0=sqrt(mu0/epsilon0);
inta1=inta0;
inta2=inta0/sqrt(epsilonr1);%1-2-3區的波阻抗
inta3=inta0/sqrt(epsilonr2);%三區的相對介電常數

%如果inta3=inta0那么k3=k1;
%！！！！！！！！！！！！！！！！！！！！！！！！！！！
%lambla/2介質演示；d=-10*pi/k2;inta3=inta1;k3=k1;
%這里做調節就可以了。
%如果inta3!=inta0，那么也會有發射發生。
%！！！！！！！！！！！！！！！！！！！！！！！！！！！

gama2=(inta3-inta2)/(inta3+inta2);%1區反射系數
intad=inta2*((exp(i*k2*d)+gama2*exp(-i*k2*d))/(exp(i*k2*d)-gama2*exp(-i*k2*d)));%等效波阻抗

%！！！！！！！！！！！！！！！！！！！！！！！！！！！
%做lambla/4介質的演示。d=-4*pi-pi/(2*k2);inta3=inta0/sqrt(5);inta2=sqrt(inta1*inta3);intad=inta1;
%！！！！！！！！！！！！！！！！！！！！！！！！！！！
framemax=500;
omega=20*pi/framemax;

delz=.2;
z1min=-30;
z1max=d;%the first zone;
z1=[z1min:delz:z1max];

z2min=d;
z2max=0;%the second zone;
z2=[z2min:delz:z2max];

y2area=[-1:2/length(z2):1-2/length(z2)];
y2area=y2area';
area2=y2area*z2;%畫出一個介質區（2區）；
area2(:,:)=1.;
area2(:,1)=0.;area2(:,length(z2))=0.;
area2(1,:)=0.;area2(length(z2),:)=0.;


z3min=0;
z3max=20;%the third zone;
z3=[z3min:delz:z3max];

jr1=0;ji1=0;ji2=0;jr2=0;ji3=0;

for t=1:framemax;

contourf(z2,y2area,area2);
colormap winter;%畫出一個介質區（2區）；
hold on;            grid on;

    
    ei1=ei10*exp(i*omega*t-i*k1*(d+z1));
    %1區入射波；
    er1=((intad-inta1)/(intad+inta1))*ei10*exp(i*omega*t+i*k1*(d+z1));
    %1區反射波；
    %公式：er1=((intad-inta1)/(intad+inta1))*ei10*exp項
    
    ei2=(1-abs(gama2))*ei10*exp(i*omega*t-i*k2*(d+z2));
    %2區入射波；
    %公式：ei2=(1-gama2)*ei10*exp項;
    er2=abs(gama2)*(1-abs(gama2))*ei10*exp(i*omega*t+i*k2*(d+z2));
    %2區反射波；
    %公式：gama2*ei20*exp項;
    
    ei3=inta3*((1-abs(gama2))-abs(gama2)*(1-abs(gama2)))*ei10*exp(i*omega*t-i*k3*(d+z3))/inta3;
    %公式：ei3=(ei20+er20)*exp項;
                if t<=length(z1)
                    if ji1<length(z1);
                        ji1=ji1+1;
                        for jij=t:1:length(z1);
                            ei1(jij)=0;
                        end
                    end
                end

           plot(z1,real(ei1),                'MarkerFaceColor','r',...
                    'MarkerEdgeColor','k',    'Color','blue',...
                        'LineWidth',2,...
                'MarkerSize',50);  %1區入射波
            hold on;

            if t>length(z1)
                if jr1<length(z1);
                    jr1=jr1+1;
                    for jr1j=1:length(z1)-jr1;
                        er1(jr1j)=0;
                    end
                end
               plot(z1,real(er1),                'MarkerFaceColor','r',...
                    'MarkerEdgeColor','k',    'Color','red',...
                        'LineWidth',2,...
                'MarkerSize',50); %1區反射波
          plot(z1,ei1+er1,                'MarkerFaceColor','r',...
                    'MarkerEdgeColor','k',    'Color','black',...
                        'LineWidth',2,...
                'MarkerSize',50); %1區反射波和入射波疊加。

            
            end
            hold on;
            if t>=length(z1);
                if ji2<length(z2);
                    ji2=ji2+1;
                    for ji2j=t+1-length(z1):1:length(z2);
                        ei2(ji2j)=0;
                    end
                end

                plot(z2,real(ei2),                'MarkerFaceColor','r',...
                    'MarkerEdgeColor','k',    'Color','magenta',...
                        'LineWidth',2,...
                'MarkerSize',50); %2區入射波。

            end
            hold on;
            if t>=length(z1)+length(z2);
                if jr2<length(z2);
                    jr2=jr2+1;
                    for jr2j=1:length(z2)-jr2;
                        er2(jr2j)=0;
                    end
                end
            plot(z2,real(er2),                'MarkerFaceColor','r',...
                    'MarkerEdgeColor','k',    'Color','yellow',...
                        'LineWidth',2,...
                'MarkerSize',50); %2區反射波。

         end
            hold on;
            if t>=length(z1)+length(z2);
                if ji3<length(z3);
                    ji3=ji3+1;
                    for ji3j=t+1-length(z1)-length(z2):1:length(z3);
                        ei3(ji3j)=0;
                    end
                end
                

            plot(z3,real(ei3),                'MarkerFaceColor','r',...
                    'MarkerEdgeColor','k',    'Color','green',...
                        'LineWidth',2,...
                'MarkerSize',50); %3區入射波。
             end

    title('藍色1區正向行波，紅色為1區負向行波（反射波），黑色為1區形成的駐波。');
    text(z1min+3,1.8*ei10,'粉紅為2區正向行波，黃色為2區反向行波。','FontSize',18,'EdgeColor','red');
    text(z1min+3,1.3*ei10,'綠色為3區正向行波。','FontSize',18,'EdgeColor','red');

            axis([z1min z3max -2 2]);
grid on;            
           hold off
    getframe;
  %  t=t/2.;
end