% 《數字信號處理教程——MATLAB釋義與實現》第八章模擬濾波器頻帶變換演示程序的子程序
% 電子工業出版社出版  陳懷琛編著 2004年9月
% call_l2bs6.m 繪制帶阻濾波器的圖形 坐標的對應關系已改  同時繪圖  data:12.12
%
function  []=call_l2bs6(h1,h2);
axes(h1);hold off;
axes(h2);hold off;

%這里設低通濾波器的通帶衰減頻率為fc=1，阻帶衰減頻率為fs=2;帶通或帶阻濾波器的中心頻率f0 可以更改這些值
%給定的參數如下：（根據給定的設計其他濾波器的轉化公式得出fc ,fs,f0
%如果有必要的話，這些參數可以從主函數傳遞過來，這樣子函數的原形變為: call_l2bs6(h1,h2,dt,T,fc,fs,f0)
%還可以在圖形上加edit控件，在調用本函數時得到其值 ,如: hedit1=uicontrol('style','edit',position',[ ]);  dt=get(hedit1,'values'); 

dt=0.1;T=3; fc=1; fs=1.5; f0=5;omegac=2;

ws1=(T-fs)/dt+1; wc1=(T-fc)/dt+1; wc2=(T+fc)/dt+2; ws2=(T+fs)/dt+2; w0=T/dt+1;
pl2=[0:dt:T];  %原始低通濾波器所選的坐標范圍，可以改變！
pl1=[-T:dt:T];
p1=pl1*j;
plp5=transform(pl2,f0,omegac,'lp2lp');
p5 =plp5*j;
php2=transform(pl2,f0,omegac,'lp2hp');  %從低通濾波器坐標范圍變換來的高通濾波器的范圍,一一對應的關系
p2=php2*j;
pbp3=transform(pl2,f0,omegac,'lp2bp');
p3=j*pbp3;
pbs4=transform(pl2,f0,omegac,'lp2bs');  %在作為帶阻濾波器設計時可以用到
p4=j*pbs4;
for i=1:length(p1)
    hl1(i)=1/abs((p1(i)^2 + 0.6180*p1(i) + 1)*(p1(i)^2 + 1.6180*p1(i) + 1)*(p1(i)+1));
end
hh1=hl1;
hbp1=hl1;
hbs1=hl1;

p1=[-T:dt:T];  %在定坐標時要把虛號去掉
mark1=1;mark2=1;

for j=1:2
    for i=2:length(p2)
        %低通濾波器
        axes(h1);
        if (j==2 & mark1==1)
            mark1=0;
            hold off;
        end
        if (i<=ws1 | i>=ws2)  %這里10表示 1/0.1;
            plot([p1(i-1),p1(i)],[hl1(i-1),hl1(i)],'color','r');axis([-T,T,0,1.2]);hold on
        elseif ((i>ws1 & i<=wc1)|(i>=wc2 & i<ws2)) %這里24表示2.4/0.1;
            plot([p1(i-1),p1(i)],[hl1(i-1),hl1(i)],'color','b');axis([-T,T,0,1.2]);hold on
        else
            plot([p1(i-1),p1(i)],[hl1(i-1),hl1(i)],'color','c');axis([-T,T,0,1.2]);hold on
        end
        
        
        if (j==1)
            
            %帶阻濾波器的繪制
            axes(h2);
            if (i<=ws1 | i>=ws2)  %這里10表示 0.5/0.1;帶通濾波器關于中心頻率處對稱，具體計算參見p166
                plot([-pbs4(i-1),-pbs4(i)],[hbs1(i-1),hbs1(i)],'color','r');axis([-f0-T,f0+T,0,1.2]);hold on
            elseif ((i>ws1 & i<=wc1) | (i>=wc2 & i<ws2))  %這里13表示(6.3-5)/0.1;
                plot([-pbs4(i-1),-pbs4(i)],[hbs1(i-1),hbs1(i)],'color','b');axis([-f0-T,f0+T,0,1.2]);hold on
            elseif (i==w0)
                plot([-pbs4(i),-pbs4(i)],[hbs1(i),hbs1(i)],'color','c');axis([-f0-T,f0+T,0,1.2]);hold on
            else
                plot([-pbs4(i-1),-pbs4(i)],[hbs1(i-1),hbs1(i)],'color','c');axis([-f0-T,f0+T,0,1.2]);hold on
            end
        else
            
            %帶阻濾波器的繪制
            axes(h2);
            if (i<=ws1| i>=ws2)  %這里10表示 0.5/0.1;帶通濾波器關于中心頻率處對稱，具體計算參見p166
                plot([pbs4(i-1),pbs4(i)],[hbs1(i-1),hbs1(i)],'color','r');axis([-f0-T,f0+T,0,1.2]);hold on
            elseif ((i>ws1 & i<=wc1) | (i>=wc2 & i<ws2))  %這里13表示(6.3-5)/0.1;
                plot([pbs4(i-1),pbs4(i)],[hbs1(i-1),hbs1(i)],'color','b');axis([-f0-T,f0+T,0,1.2]);hold on
            elseif (i==w0)
                plot([pbs4(i),pbs4(i)],[hbs1(i),hbs1(i)],'color','c');axis([-f0-T,f0+T,0,1.2]);hold on
            else
                plot([pbs4(i-1),pbs4(i)],[hbs1(i-1),hbs1(i)],'color','c');axis([-f0-T,f0+T,0,1.2]);hold on
            end
        end
        pause(0.1);
    end
    pause(1);
end