fdfirls.m

matlabDigitalSigalProcess内有文件若干

    order.value = N;
    wt1 = pbspecs(3);
    wt2 = sbspecs(3);
    set(wt1,'value',wt(1+(type == 2 | type == 3)),...
       'label','Weight','callback','fdfirls(''dirty'')')
    set(wt2,'value',wt(2-(type == 2 | type == 3)),...
       'label','Weight','callback','fdfirls(''dirty'')')
    L1.color = co(min(2,size(co,1)),:);
    L2.color = co(min(2,size(co,1)),:);

    switch type
    case 1 % lowpass
        f1 = pbspecs(1);
        f2 = sbspecs(1);
        pbspecs(2).visible = 'off';
        sbspecs(2).visible = 'off';
        set(f1,'value',f(2)*Fs/2,'label','Fp')
        set(f2,'value',f(3)*Fs/2,'label','Fs')
    case 2 % highpass
        f1 = sbspecs(1);
        f2 = pbspecs(1);
        pbspecs(2).visible = 'off';
        sbspecs(2).visible = 'off';
        set(f1,'value',f(2)*Fs/2,'label','Fs')
        set(f2,'value',f(3)*Fs/2,'label','Fp')
    case 3 % bandpass
        f2 = pbspecs(1);
        f3 = pbspecs(2);
        f1 = sbspecs(1);
        f4 = sbspecs(2);
        pbspecs(2).visible = 'on';
        sbspecs(2).visible = 'on';
        set(f1,'value',f(2)*Fs/2,'label','Fs1')
        set(f2,'value',f(3)*Fs/2,'label','Fp1')
        set(f3,'value',f(4)*Fs/2,'label','Fp2')
        set(f4,'value',f(5)*Fs/2,'label','Fs2')
    case 4
        f1 = pbspecs(1);
        f4 = pbspecs(2);
        f2 = sbspecs(1);
        f3 = sbspecs(2);
        pbspecs(2).visible = 'on';
        sbspecs(2).visible = 'on';
        set(f1,'value',f(2)*Fs/2,'label','Fp1')
        set(f2,'value',f(3)*Fs/2,'label','Fs1')
        set(f3,'value',f(4)*Fs/2,'label','Fs2')
        set(f4,'value',f(5)*Fs/2,'label','Fp2')
    end
    ax.xlimPassband = fdutil('xlimpassband',type,...
                         Fs,f1.value,f2.value,f3.value,f4.value);

    if ax.xlimbound(2) ~= Fs/2
        set(ax,'xlimbound',[0 Fs/2],'xlim',[0 Fs/2])
    end
    
    f = f(:)';
    minpass = L1.ydata(end);
    maxpass = L1.ydata(1);
    minstop = L2.ydata(1);
    fdutil('setLines','fdfirls',L1,L2,0,type,f,Fs,minpass,maxpass,minstop,1)

%---------------------------------------------------------------------
% fdfirls('newtype')
%  Button down fcn band configuration popup
case 'newtype'
    filt = filtdes('filt');
    newtype = bandpop.value;
    oldtype = passframe.userdata;
    if isempty(oldtype)
        oldtype = filt.specs.fdfirls.type;
    end
    passframe.userdata = newtype;
    if (newtype ~= oldtype) | strcmp(filtdes('getenable'),'on')
        if oldtype < 3  % low or high pass
            edges = [f1.value f2.value]'*2/Fs;
        else
            edges = [f1.value f2.value f3.value f4.value]'*2/Fs;
        end
        edges = fdutil('changeFilterType',newtype,oldtype,edges);
        f = [0; edges(:); 1];
        wp = get(wt1,'value');
        ws = get(wt2,'value');
        switch newtype
        case 1
            wt = [wp ws];
        case 2
            wt = [ws wp];
        case 3
            wt = [ws wp ws];
        case 4
            wt = [wp ws wp];
        end
        fdfirls('newfilt',1,...
                     newtype,f,Rp.value,...
                     Rs.value,wt,...
                     order.value)
        fdfirls('dirty')
    else
       % disp('no change of type')
    end
    
%---------------------------------------------------------------------
% fdfirls('Lrespdown')
%  Button down fcn of Lresp (response line) - pan
case 'LrespDown'
    bounds.xlim = [0 Fs/2];
    bounds.ylim = [-500 30];
    h = ax.h;
    panfcn('Ax',h,...
           'Bounds',bounds,...
           'UserHand',get(h,'zlabel'))
           
%---------------------------------------------------------------------
% fdfirls('L1down')
%  Button down fcn of L1
case 'L1down'
    L1.erasemode = 'xor';
    L2.erasemode = 'xor';
    
%---------------------------------------------------------------------
% fdfirls('L2down')
%  Button down fcn of L2
case 'L2down'
    L1.erasemode = 'xor';
    L2.erasemode = 'xor';
    
%---------------------------------------------------------------------
% fdfirls('L1up')
%  Button up fcn of L1
case 'L1up'
    L1.erasemode = 'normal';
    L2.erasemode = 'normal';
    Lresp.erasemode = 'normal';
    
%---------------------------------------------------------------------
% fdfirls('L2up')
%  Button up fcn of L2
case 'L2up'
    L1.erasemode = 'normal';
    L2.erasemode = 'normal';
    Lresp.erasemode = 'normal';

%---------------------------------------------------------------------
% fdfirls('L1drag',type,ind)
%  drag callback of L1 - passband line
%  Inputs:
%     type - band configuration 1==low, 2=high, 3=pass, 4=stop
%     ind - index of vertex being dragged
case 'L1drag'
    type = varargin{2};
    ind = varargin{3};
    xd = L1.xdata;
    minspacing = Fs/500;
    switch type
    case 1  % lowpass
        newf1 = inbounds(xd(ind),[minspacing f2.value-minspacing]);
        xd([2 5]) = newf1;
        L1.xdata = xd;
        f1.value = newf1;
        i = 1;  % index of changed frequency
    case 2  % highpass
        newf2 = inbounds(xd(ind),[f1.value+minspacing Fs/2-minspacing]);
        xd([1 4]) = newf2;
        L1.xdata = xd;
        f2.value = newf2;
        i = 2;
    case 3  % bandpass
    % L1 xdata = [f2 f3 NaN f2 f3]
        if any(ind == [1 4])   % dragging f2
            newf2 = inbounds(xd(ind),[f1.value+minspacing f3.value-minspacing]);
            xd([1 4]) = newf2;
            L1.xdata = xd;
            f2.value = newf2;
            i = 2;
        else % dragging f3
            newf3 = inbounds(xd(ind),[f2.value+minspacing f4.value-minspacing]);
            xd([2 5]) = newf3;
            L1.xdata = xd;
            f3.value = newf3;
            i = 3;
        end
    case 4   % bandstop
    % L1 xdata = [0 f1 NaN 0 f1 NaN f4 Fs/2 NaN f4 Fs/2]
        if any(ind == [2 5])   % dragging f1
            newf1 = inbounds(xd(ind),[minspacing f2.value-minspacing]);
            xd([2 5]) = newf1;
            L1.xdata = xd;
            f1.value = newf1;
            i = 1;
        else % dragging f4
            newf4 = inbounds(xd(ind),[f3.value+minspacing Fs/2-minspacing]);
            xd([7 10]) = newf4;
            L1.xdata = xd;
            f4.value = newf4;
            i = 4;
        end
    end
    if type < 3
        f = {f1.value, f2.value};
    else
        f = {f1.value, f2.value, f3.value, f4.value};
    end
    ax.xlimPassband = fdutil('xlimpassband',type,Fs,f{:});
    
%---------------------------------------------------------------------
% fdfirls('L2drag',type,ind)
%  drag callback of L2 - stopband line
%  Inputs:
%     type - band configuration 1==low, 2=high, 3=pass, 4=stop
%     ind - index of vertex being dragged
case 'L2drag'
    type = varargin{2};
    ind = varargin{3};
    xd = L2.xdata;
    minspacing = Fs/500;
    switch type
    case 1  % lowpass
        newf2 = inbounds(xd(ind),[f1.value+minspacing Fs/2-minspacing]);
        xd(1) = newf2;
        L2.xdata = xd;
        f2.value = newf2;
        i = 2;
    case 2  % highpass
        newf1 = inbounds(xd(ind),[minspacing f2.value-minspacing]);
        xd(2) = newf1;
        L2.xdata = xd;
        f1.value = newf1;
        i = 1;
    case 3  % bandpass
    % L2 xdata = [0 f1 NaN f4 Fs/2]
        if ind == 2   % dragging f1
            newf1 = inbounds(xd(ind),[minspacing f2.value-minspacing]);
            xd(2) = newf1;
            L2.xdata = xd;
            f1.value = newf1;
            i = 1;
        else % dragging f4
            newf4 = inbounds(xd(ind),[f3.value+minspacing Fs/2-minspacing]);
            xd(4) = newf4;
            L2.xdata = xd;
            f4.value = newf4;
            i = 4;
        end
    case 4   % bandstop
    % L2 xdata = [f2 f3]
        if ind == 1   % dragging f2
            newf2 = inbounds(xd(ind),[f1.value+minspacing f3.value-minspacing]);
            xd(1) = newf2;
            L2.xdata = xd;
            f2.value = newf2;
            i = 2;
        else % dragging f3
            newf3 = inbounds(xd(ind),[f2.value+minspacing f4.value-minspacing]);
            xd(2) = newf3;
            L2.xdata = xd;
            f3.value = newf3;
            i = 3;
        end
    end
    
%---------------------------------------------------------------------
% fdfirls('Rpdrag',type,ind)
%  drag callback of L1 - passband line
%  Inputs:
%     type - band configuration 1==low, 2=high, 3=pass, 4=stop
%     ind - index of segment being dragged
case 'Rpdrag'
    type = varargin{2};
    ind = varargin{3};
    yd = L1.ydata;
    
    switch ind
    case {4,10}    % dragging lower line
        below = yd(ind);
        if below >= 0 
            below = -.00001;
        elseif below < -10
            below = -10;
        end
        dev = 1-10^(below/20);
        above = 20*log10(1+dev);
    case {1,7}    % dragging upper line
        above = yd(ind);
        if above > 10 
            above = 10;
        elseif above <= 0
            above = .00001;
        end
        dev = 10^(above/20)-1;
        below = 20*log10(1-dev);
    end
    
    newRp = above-below;
    if type ~= 4
    % 'ydata':[maxpass maxpass NaN minpass minpass]
        yd = [above above NaN below below];
    else
    % 'ydata': [ maxpass maxpass NaN minpass minpass NaN ...
    %                    maxpass maxpass NaN minpass minpass])
        yd = [above above NaN below below ...
               NaN above above NaN below below ];
    end
    L1.ydata = yd;
    ylim = [below above];
    dyl = (ylim(2)-ylim(1))*.15;
    ax.ylimPassband = ylim + [-dyl/2 dyl/2];
    
    devs = [(10^(newRp/20)-1)/(10^(newRp/20)+1) 10^(-Rsm.value/20)].^2;
    bw = bandwidths(bandpop.value,f1.value,f2.value,f3.value,f4.value,Fs);
    wt = 1./(devs.*bw);
    wt = wt/min(wt);
    set(wt1,'value',wt(1))
    set(wt2,'value',wt(2))
 
%---------------------------------------------------------------------
% fdfirls('Rsdrag',type,ind)
%  drag callback of L2 - stopband line
%  Inputs:
%     type - band configuration 1==low, 2=high, 3=pass, 4=stop
%     ind - index of segment being dragged
case 'Rsdrag'
    type = varargin{2};
    ind = varargin{3};
  
    yd = L2.ydata;
    newRs = -yd(ind);
    if newRs < 0
        newRs = 0;
    end
    switch type
    case {1,2,4}
        L2.ydata = [-newRs -newRs];
    case 3
        L2.ydata = [-newRs -newRs NaN -newRs -newRs];
    end
    devs = [(10^(Rpm.value/20)-1)/(10^(Rpm.value/20)+1) 10^(-newRs/20)].^2;
    bw = bandwidths(bandpop.value,f1.value,f2.value,f3.value,f4.value,Fs);
    wt = 1./(devs.*bw);
    wt = wt/min(wt);
    set(wt1,'value',wt(1))
    set(wt2,'value',wt(2))
        
end  % of function switch-yard

%---------------------------------------------------------------------
% -------- LOCAL FUNCTIONS START HERE  --------- 
%---------------------------------------------------------------------

function sethelp
global minOrdCheckbox bandpop 
global passframe stopframe passframe1 stopframe1
global f1 f2 f3 f4 wt1 wt2 order Rp Rm pbspecs sbspecs
global Rpm Rsm wt1m wt2m
global ax L1 L2 Lresp
global Fs

% disp('setting help ... stub')

 
function yd = passbandlimits(Rp)
% return ydata = [minpass maxpass] of passband 
%   given Rp decibels of ripple
%   in passband (centered at 1 in linear scale)
    dev = (10^(Rp/20)-1)/(10^(Rp/20)+1);
    above = 20*log10(1+dev); below = 20*log10(1-dev);
    yd = [below above];


function [maxpass,minpass,minstop] = getMeasurements(ff,Hlog,type,...
                                          f1,f2,f3,f4);
    switch type
    case 1   % lowpass
        passInd = find(ff<=f1);
        stopInd = find(ff>=f2);
    case 2   % highpass
        stopInd = find(ff<=f1);
        passInd = find(ff>=f2);
    case 3   % bandpass
        stopInd = find((ff<=f1)|(ff>=f4));
        passInd = find((ff>=f2)&(ff<=f3));
    case 4   % bandstop
        passInd = find((ff<=f1)|(ff>=f4));
        stopInd = find((ff>=f2)&(ff<=f3));
    end
    maxpass = max(Hlog(passInd));
    minpass = min(Hlog(passInd));
    minstop = max(Hlog(stopInd));


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

    % first define default values
    setOrderFlag_init = 1;   % ALWAYS 1 for firls, since we can't
                             % estimate the order...
                             % this field is a placeholder used by
                             % the other modules 
    type_init = 1;  % 1=lowpass, 2=highpass, 3=bandpass, 4=bandstop
    f_init = [0 .1 .15 .5];
    Rp_init = 3;   % for measurements only
    Rs_init = 20;  % for measurements only
    order_init = 30;
    wt_init = [1 1];
    switch filt.specs.currentModule
    case {'fdcheby1','fdcheby2','fdellip','fdbutter','fdkaiser','fdremez'}
        s = eval(['filt.specs.' filt.specs.currentModule]);
        type_init = s.type;
        f_init = s.f;
        order_init = s.order;
        
        if strcmp(filt.specs.currentModule,'fdremez')
            wt_init = s.wt;
        else
            devs = [(10^(s.Rp/20)-1)/(10^(s.Rp/20)+1) 10^(-s.Rs/20)];
            switch s.type
            case 3
                devs = devs([2 1 2]);
            case 4
                devs = devs([1 2 1]);
            end
            wt_init = ones(size(devs))*max(devs)./devs;
        end
        if any(strcmp(filt.specs.currentModule,...
                 {'fdcheby1','fdcheby2','fdellip','fdbutter'}))
            order_init = order_init*10;  %  FIR filters are much higher order
                                         %  than IIR
        end
    
    case 'fdfirls'
        if isfield(filt.specs,'fdfirls')
            type_init = filt.specs.fdfirls.type;
            f_init = filt.specs.fdfirls.f;
            Rp_init = filt.specs.fdfirls.Rp;
            Rs_init = filt.specs.fdfirls.Rs;
            order_init = filt.specs.fdfirls.order;
            wt_init = filt.specs.fdfirls.wt;
        end
    end

function bw = bandwidths(type,f1,f2,f3,f4,Fs)
% Returns bandwidths of passband and stopband in a two element vector
% The larger of the two bandwidths will be 1, and the smaller
% is expressed as a fraction of the larger

    switch type
    case 1 % lowpass
        bw = [f1 Fs/2-f2];
    case 2 % high
        bw = [Fs/2-f2 f1];
    case 3 % bandpass
        bw = [f3-f2 f1+(Fs/2-f4)];
    case 4 % bandstop
        bw = [f1+(Fs/2-f4) f3-f2];
    end
    bw = bw/max(bw);