fdkaiser.m

matlabDigitalSigalProcess内有文件若干

        if any(passbandChange==0)
            f2.value = fm(2);
            i = [i 2];
        end
    case 2
        if any(passbandChange==1)
            f2.value = fm(2);
            i = 2;
        end
        if any(passbandChange==0)
            f1.value = fm(1);
            i = [1 i];
        end
    case 3 % bandpass - update f2 f3
        if any(passbandChange==1)
            f2.value = fm(2);
            f3.value = fm(3);
            i = [2 3];
        end
        if any(passbandChange==0)
            f1.value = fm(1);
            f4.value = fm(4);
            i = [1 i 4];
        end
    case 4 % bandstop - update f1 f4
        if any(passbandChange==1)
            f1.value = fm(1);
            f4.value = fm(4);
            i = [1 4];
        end
        if any(passbandChange==0)
            f2.value = fm(2);
            f3.value = fm(3);
            i = sort([i 2 3]);
        end
    end
    updateLines(passbandChange,type,i,f1,f2,f3,f4,Rp,Rs,Fs,Rp_range)
    fdutil('pokeFilterMeasurements','fdkaiser',type,f1,f2,f3,f4,Rp,Rs,Fs)
    
function updateLines(passbandChange,type,fchange_ind,f1,f2,f3,f4,Rp,Rs,Fs,Rp_range)
% assume values of f1,f2,f3,f4,Rp and Rs are correct now
% fchange_ind - vector of indices indicating which frequencies have
%               changed
    global L1 L2 ax
    
    f = getFrequencyValues(type,f1,f2,f3,f4,Fs)*Fs/2;
    
    if any(passbandChange==1)
        if nargin < 11
            Rp_range = passbandlimits(Rp.value);
        end
        minpass = Rp_range(1);
        maxpass = Rp_range(2);
        % update L1 xdata and ydata
        switch type
        case 1   % lowpass
            set(L1,'xdata',[f(1:2) NaN f(1:2)],...
                   'ydata',[maxpass maxpass NaN minpass minpass])
        case 2   % highpass
            set(L1,'xdata',[f(3:4) NaN f(3:4)],...
                   'ydata',[maxpass maxpass NaN minpass minpass ])
        case 3   % bandpass
            set(L1,'xdata',[f(3:4) NaN f(3:4)],...
                   'ydata',[maxpass maxpass NaN minpass minpass])
        case 4   % bandstop
            set(L1,'xdata',[f(1:2) NaN f(1:2) NaN f(5:6) NaN f(5:6)],...
                   'ydata',[ maxpass maxpass NaN minpass minpass NaN ...
                             maxpass maxpass NaN minpass minpass])
        end
        ylim = [minpass maxpass];
        dyl = (ylim(2)-ylim(1))*.15;
        ax.ylimPassband = ylim + [-dyl/2 dyl/2];
        if length(f)==4
            f(6)=0;  % zeropad for call to xlimpassband
        end
        ax.xlimPassband = fdutil('xlimpassband',type,...
                             Fs,f(2),f(3),f(4),f(5));
    end
    if any(passbandChange==0)
        minstop = -Rs.value;
        % update L2 xdata and ydata
        switch type
        case 1   % lowpass
            set(L2,'xdata',[f(3:4)],'ydata',[minstop minstop])
        case 2   % highpass
            set(L2,'xdata',[f(1:2)],'ydata',[ minstop minstop])
        case 3   % bandpass
            set(L2,'xdata',[f(1:2) NaN f(5:6)],...
                   'ydata',[ minstop minstop NaN minstop minstop])
        case 4   % bandstop
            set(L2,'xdata',[f(3:4)],'ydata',[minstop minstop])
        end
    end
%     fdutil('updateRanges',fchange_ind,f1,f2,f3,f4)
    set(f1,'range',[0 Fs/2])
    set(f2,'range',[0 Fs/2])
    set(f3,'range',[0 Fs/2])
    set(f4,'range',[0 Fs/2])
    
function [f1,f2,f3,f4] = setupFrequencyObjects(type,module,pbobjects,...
                   sbobjects,f,Fs,ax,setValueFlag)

    if nargin<8
        setValueFlag = 1;
    end
    
    fdutil('changeToEdit',[pbobjects(1:2); sbobjects(1:2)])
    
    switch type
    case 1 % lowpass
        f1 = pbobjects(1);
        f2 = sbobjects(1);
        f3 = sbobjects(2);
        f4 = pbobjects(2);
        pbobjects(2).visible = 'off';
        sbobjects(2).visible = 'off';
        set(f1,'label','Fp')
        set(f2,'label','Fs')
    case 2 % highpass
        f1 = sbobjects(1);
        f2 = pbobjects(1);
        f3 = pbobjects(2);
        f4 = sbobjects(2);
        pbobjects(2).visible = 'off';
        sbobjects(2).visible = 'off';
        set(f1,'label','Fs')
        set(f2,'label','Fp')
    case 3 % bandpass
        f1 = sbobjects(1);
        f2 = pbobjects(1);
        f3 = pbobjects(2);
        f4 = sbobjects(2);
        pbobjects(2).visible = 'on';
        sbobjects(2).visible = 'on';
        set(f1,'label','Fs1')
        set(f2,'label','Fp1')
        set(f3,'label','Fp2')
        set(f4,'label','Fs2')
    case 4
        f1 = pbobjects(1);
        f2 = sbobjects(1);
        f3 = sbobjects(2);
        f4 = pbobjects(2);
        pbobjects(2).visible = 'on';
        sbobjects(2).visible = 'on';
        set(f1,'label','Fp1')
        set(f2,'label','Fs1')
        set(f3,'label','Fs2')
        set(f4,'label','Fp2')
    end
    if setValueFlag
        if type < 3
            set(f1,'value',f(2)*Fs/2,'range',[0 f(3)]*Fs/2)
            set(f2,'value',f(3)*Fs/2,'range',[f(2) 1]*Fs/2)
        else
            set(f1,'value',f(2)*Fs/2,'range',[0 f(3)*Fs/2])
            set(f2,'value',f(3)*Fs/2,'range',[f(2) f(4)]*Fs/2)
            set(f3,'value',f(4)*Fs/2,'range',[f(3) f(5)]*Fs/2)
            set(f4,'value',f(5)*Fs/2,'range',[f(4) 1]*Fs/2)
        end
        ax.xlimPassband = fdutil('xlimpassband',type,...
                             Fs,f1.value,f2.value,f3.value,f4.value);
    end

    f1.callback = [module '(''fchange'',1)'];
    f2.callback = [module '(''fchange'',2)'];
    f3.callback = [module '(''fchange'',3)'];
    f4.callback = [module '(''fchange'',4)'];

function [Rp,Rs,f1,f2,f3,f4] = ...
      setupMeasurementObjects(specStruc,Fs,order,order1,passframe1,stopframe1,ax,...
                     pbmeasures,sbmeasures,Rp,Rs,...
                     Wn1m,Wn2m,Betam,Rpm,Rsm,f1,f2,f3,f4)
%  set values of MEASUREMENTS objects ... assumes specStruct is current

    setOrderFlag = specStruc.setOrderFlag;
    type = specStruc.type;
    f = specStruc.f;
    Rpass = specStruc.Rp;
    Rstop = specStruc.Rs;
    Wn = specStruc.Wn;
    N = specStruc.order;
    beta_val = specStruc.Beta;
    
    if ~setOrderFlag  % estimate order

        set(order1,'visible','on')

        set(passframe1,'visible','off')
        set(stopframe1,'visible','off')
        
        fdutil('changeToText',[Rpm Rsm Wn1m Wn2m Betam])
        set(sbmeasures(3),'visible','off')
        
        set(Rpm,'visible','on','label','Actual Rp')
        set(Rsm,'visible','on','label','Actual Rs')
        set(Betam,'visible','on','label','Beta','value',beta_val)
        
        set(Wn1m,'visible','on','value',Wn(1)*Fs/2,...
            'format','%1.4g')
        if type > 2
            set(Wn1m,'label','Fc1')
            set(Wn2m,'visible','on','value',Wn(2)*Fs/2,...
                'format','%1.4g','label','Fc2')
        else
            set(Wn1m,'label','Fc')
            set(Wn2m,'visible','off')
        end

        order.value = N;
        
    else
        set(order1,'visible','off')
    
        set(passframe1,'visible','on')
        set(stopframe1,'visible','on')
        
        fdutil('changeToEdit',[pbmeasures(1:3) sbmeasures(1:3)])

        pbmeasures(1).visible = 'on';
        sbmeasures(1).visible = 'on';

        [f1,f2,f3,f4] = setupFrequencyObjects(type,'fdkaiser',pbmeasures,...
                                     sbmeasures,f,Fs,ax,0);
                                             
        Rp = pbmeasures(3);        
        Rs = sbmeasures(3);        
        set(Rp,'value',Rpass,'label','Rp',...
            'visible','on','callback','fdkaiser(''Rpchange'')')
        set(Rs,'value',Rstop,'label','Rs',...
            'visible','on','callback','fdkaiser(''Rschange'')')
    
    end


function [Rp,Rs,f1,f2,f3,f4,specStruc] = ...
      measureFilt(objSetupFlag,specStruc,Fs,order,order1,passframe1,stopframe1,ax,...
                  pbmeasures,sbmeasures,Rp,Rs,Wn1m,Wn2m,Betam,Rpm,Rsm,f1,f2,f3,f4,...
                  ff,Hlog,L1,L2)

    if objSetupFlag
        [Rp,Rs,f1,f2,f3,f4] = ...
          setupMeasurementObjects(specStruc,Fs,order,order1,passframe1,stopframe1,ax,...
                     pbmeasures,sbmeasures,Rp,Rs,...
                     Wn1m,Wn2m,Betam,Rpm,Rsm,f1,f2,f3,f4);
    end

    n = specStruc.order;
    f = specStruc.f;
    setOrderFlag = specStruc.setOrderFlag;
    Wn = specStruc.Wn;
    beta_val = specStruc.Beta;
    
    for i=1:3
        set(pbmeasures(i),'enable','on')
        set(sbmeasures(i),'enable','on')
    end
    set(order1,'enable','on')
    
    set(f1,'range',[0 Fs/2])
    set(f2,'range',[0 Fs/2])
    set(f3,'range',[0 Fs/2])
    set(f4,'range',[0 Fs/2])

    if ~setOrderFlag
        order.value = n;
        Wn1m.value = Wn(1) * Fs/2;
        if length(Wn)>1
            Wn2m.value = Wn(2) * Fs/2;
        end
        [maxpass,minpass,minstop] = getMagMeasurements(ff,Hlog,specStruc.type,...
                                          f1.value,f2.value,f3.value,f4.value);
        set(Betam,'value',beta_val)
        set(Rpm,'value',max(maxpass - minpass))
        set(Rsm,'value',-max(minstop))
        yd = passbandlimits(Rp.value);
        maxpass = yd(2);
        minpass = yd(1);
        minstop = -Rs.value;
    else
        [fm,Rp_range,Rs_val] = ...
             getFreqMeasurements(1,ff,Hlog,specStruc.type,Rp.value,Rs.value);
        f1.value = fm(1);  f2.value = fm(2);
        if length(fm)>2
            f3.value = fm(3);  f4.value = fm(4);
        end
        Rp.value = diff(Rp_range);
        maxpass = Rp_range(2);   minpass = Rp_range(1);   minstop = -Rs_val;
        Rs.value = Rs_val;
        f = [0 fm(:)' Fs/2 0];  % for updating L1 and L2 later;
                                % extra zero is for call to xlimpassband when
                                % type == 1 or 2
        ax.xlimPassband = fdutil('xlimpassband',specStruc.type,...
                             Fs,f(2),f(3),f(4),f(5));
        f = f*2/Fs;
        specStruc.f = f;    
    end
    ylim = [minpass maxpass];
    dyl = (ylim(2)-ylim(1))*.15;
    ax.ylimPassband = ylim + [-dyl/2 dyl/2];
    
    % update L1 and L2 xdata and ydata, pointers, dragmodes
    % note - can't drag upper passband measurement for bandpass
    fdutil('setLines','fdkaiser',L1,L2,0,...
                      specStruc.type,f(:)',Fs,minpass,maxpass,minstop,...
                        ~(setOrderFlag & (specStruc.type==3)))


function f = getFrequencyValues(type,f1,f2,f3,f4,Fs);            
if type < 3  % low or high pass
    f = [0 f1.value f2.value Fs/2]*2/Fs;
else
    f = [0 f1.value f2.value f3.value f4.value Fs/2]*2/Fs;
end


function [setOrderFlag_init, type_init, f_init, Rp_init, Rs_init, ...
         order_init, Wn_init, Beta_init] = initSpecs(filt)
%initSpecs   Initial specifications for kaiser filter, from
%            filt input
%Switches off of filt.currentModule and if it finds any of
% fdcheby1, fdcheby2, fdellip, fdremez, fdfirls, or fdbutter
% retains the type, order, band edge, and any other relevant
% parameters

    % first define default values
    setOrderFlag_init = 0;   % by default, estimate order 
    type_init = 1;  % 1=lowpass, 2=highpass, 3=bandpass, 4=bandstop
    f_init = [0 .2 .25 1];
    Rp_init = 3;
    Rs_init = 20;
    order_init = 30;
    Wn_init = .225;
    Beta_init = 5;
    switch filt.specs.currentModule
    case {'fdcheby1','fdcheby2','fdellip','fdremez','fdbutter','fdfirls'}
        s = eval(['filt.specs.' filt.specs.currentModule]);
        setOrderFlag_init = s.setOrderFlag;
        type_init = s.type;
        f_init = s.f;
        Rp_init = s.Rp;
        Rs_init = s.Rs;
        order_init = s.order;
        switch filt.specs.currentModule
        case 'fdcheby1'
            Wn_init = s.Fpass;
        case 'fdcheby2'
            Wn_init = s.Fstop;
        case 'fdellip'
            Wn_init = s.Fpass;
        case {'fdremez','fdfirls'}
            switch s.type
            case {1,2}
                Wn_init = s.f(2);
            case 3
                Wn_init = s.f(3:4);
            case 4
                Wn_init = s.f([2 5]);
            end            
        case 'fdbutter'
            Wn_init = s.w3db;
        end
        if any(strcmp(filt.specs.currentModule,...
                 {'fdbutter','fdcheby1','fdcheby2','fdellip'}))
            order_init = ceil(order_init*10); 
             %  FIR filters are much higher order than IIR
        end
        if isfield(filt.specs,'fdkaiser')
            Beta_init = filt.specs.fdkaiser.Beta;
        end
    case 'fdkaiser'
        if isfield(filt.specs,'fdkaiser')
            setOrderFlag_init = filt.specs.fdkaiser.setOrderFlag;
            type_init = filt.specs.fdkaiser.type;
            f_init = filt.specs.fdkaiser.f;
            Rp_ini