qiangdu.m

这是一个计算一维光子晶体滤波器传输特性的软件

clear all
wavelenth=600*10^(-9); %中心波长
N_w=1000;
n_a=4.6;n_b=1.6;n_base=1;   %介质的折射率参数
d_a=0.25*wavelenth/n_a;  d_ai=(0.25-0.089)*wavelenth/n_a; d_a0=(0.25+0.089)*wavelenth/n_a;%A介质的厚度
d_b=0.25*wavelenth/n_b;  d_bi=(0.25-0.089)*wavelenth/n_b; d_b0=(0.25+0.089)*wavelenth/n_b;%B介质的厚度
for u=1:2          %TE和TM模的选择
    for p=1:91         %角度值范围
    theta(p)=0+(p-1).*1;    %入射角
    theta_0=theta(p)/180*pi;  %入射角弧度
        for t=1:(N_w+1)        % 频率的变化范围
            c=3.*10^(8);         %光速
            w_0=2.*pi.*c./wavelenth;  %中心波长的角频率
            w_r=0+(t-1).*2./N_w;      %相对频率的份数从0：2倍
            w=w_0.*w_r;                          %相对频率的角频率
            theta_a=asin(sin(theta_0)./n_a);    % A材料的折射角  
            theta_b=asin(sin(theta_0)./n_b);   %B材料折射角
            itate_a=n_a.*cos(theta_a);          % TE A导纳
            dirta_a=w.*n_a.*d_a.*cos(theta_a)./c;dirta_ai=w.*n_a.*d_ai.*cos(theta_a)./c;dirta_a0=w.*n_a.*d_a0.*cos(theta_a)./c;%TE A位相角 
            itatm_a=n_a./cos(theta_a);          %TM  A导纳
            itate_b=n_b.*cos(theta_b);           %TE  B导纳
            dirta_b=w.*abs(n_b).*d_b.*cos(theta_b)./c;dirta_bi=w.*n_b.*d_bi.*cos(theta_b)./c;dirta_b0=w.*n_b.*d_b0.*cos(theta_b)./c;%TE  B位相角  
            itatm_b=n_b./cos(theta_b);            %TM B导纳  
            ita0=[cos(theta_0),1/cos(theta_0)];%空气介质TE和TM模的选择的导纳计算方法
            ita_0=ita0(u);   %空气介质的导纳
            theta_base=asin(sin(theta_0)./n_base);%入射角空气
            itabase=[n_base.*cos(theta_base),n_base./cos(theta_base)];  %空气介质的导纳
            ita_base=itabase(u);      %空气介质的U导纳
            itaa=[itate_a,itatm_a];   %A介质的U导纳
            itab=[itate_b,itatm_b];   %B介质的U导纳             
            ita_a=itaa(u);ita_b=itab(u);%介质的U导纳选择
            M_A=[cos(dirta_a),-i.*sin(dirta_a)./ita_a;-i.*sin(dirta_a).*ita_a,cos(dirta_a)];%A介质的矩阵
            M_B=[cos(dirta_b),-i.*sin(dirta_b)./ita_b;-i.*sin(dirta_b).*ita_b,cos(dirta_b)];%B介质的矩阵
            M_Ai=[cos(dirta_ai),-i.*sin(dirta_ai)./ita_a;-i.*sin(dirta_ai).*ita_a,cos(dirta_ai)];%Ai介质的矩阵
            M_Bi=[cos(dirta_bi),-i.*sin(dirta_bi)./ita_b;-i.*sin(dirta_bi).*ita_b,cos(dirta_bi)];%Bi介质的矩阵
            M_A0=[cos(dirta_a0),-i.*sin(dirta_a0)./ita_a;-i.*sin(dirta_a0).*ita_a,cos(dirta_a0)];%Ai介质的矩阵
            M_B0=[cos(dirta_b0),-i.*sin(dirta_b0)./ita_b;-i.*sin(dirta_b0).*ita_b,cos(dirta_b0)];%Bi介质的矩阵
            M_total=(M_Bi*M_Ai)^8*(M_B*M_A)^10*(M_A*M_B)^10*(M_A0*M_B0)^8;
            %M_total=(M_A0*M_B0)^8;
            A=M_total(1,1);B=M_total(1,2);C=M_total(2,1);D=M_total(2,2);                    %矩阵个元素的值
            R(t,p,u)=abs((A.*ita_0+B.*ita_0.*ita_base-C-D.*ita_base)./(A.*ita_0+B.*ita_0.*ita_base+C+D.*ita_base)).^2;%反射率
            T(t,p,u)=abs(2.*ita_0./(A.*ita_0+B.*ita_0.*ita_base+C+D.*ita_base)).^2;  %透射率
        end
    end
end
T_te=T(:,:,1);T_tm=T(:,:,2);T_total=[fliplr(T_tm(:,2:91)),T_te];
R_te=R(:,:,1);R_tm=R(:,:,2);R_total=[fliplr(R_tm(:,2:91)),R_te];
clear j;clear t;clear k; clear p; clear u;clear s;