Designed filter shows half the cut-off frequency that it should

hi everyone,

I tried to implement the designed filter (with fdesign), using its coefficients and test it with a chirp signal. here is my code:

clear screen

clearvars

format long

SamplingFrequency = 2000;

SamplingPeriod = 1/SamplingFrequency;

CenterFrequency = 450;

TransitionWidth = 100;

PassbandRipple = 1;

StopbandAttenuation = 80;

DesignSpec = fdesign.lowpass('Fp,Fst,Ap,Ast',...

CenterFrequency-TransitionWidth/2, ...

CenterFrequency+TransitionWidth/2, ...

PassbandRipple,StopbandAttenuation, ...

SamplingFrequency);

LPF = design(DesignSpec,'equiripple','SystemObject',true);

fvtool(LPF)

Coefficients = LPF.Numerator;

%%% here we define different input scenarios to be tested on the FIR filter

%%% Impulse Response

FIR_In = cat(2, 1, zeros(1,size(Coefficients, 2)));

FIR_Out = myFIR(FIR_In, Coefficients);

%%% Chirp Response

clear FIR_Out FIR_In

TotalTime = 100; %seconds

StartFrequency = 5; %Herz

StopFrequency = 1000; %Herz

t = 0:SamplingPeriod:TotalTime;

InstantanousFrequency = linspace(StartFrequency,StopFrequency,size(t,2));

disp (size(InstantanousFrequency))

x = zeros(1, size(t,2));

for i = 1:size(t,2)

x(i) = cos (2*pi*InstantanousFrequency(i)*t(i));

end

FIR_In = x;

FIR_In = cat(2, FIR_In, zeros(1,size(Coefficients, 2)));

FIR_Out = myFIR(FIR_In, Coefficients);

figure

plot (InstantanousFrequency,FIR_In(1,1:size(InstantanousFrequency,2)),InstantanousFrequency,FIR_Out(1,1:size(InstantanousFrequency,2)))

function FIR_Out = myFIR(FIR_In, Coefficients)

FIR_Out = zeros(1, size(FIR_In, 2));

x = zeros(1, size(Coefficients, 2));

for CurrentInputIndex = 1:1:size(FIR_In, 2)

%%% delay line

for DelayLineIndex = size(Coefficients, 2):-1:2

x(DelayLineIndex) = x(DelayLineIndex-1);

end

x(1) = FIR_In(CurrentInputIndex);

%%%

Acc = 0;

for TapIndex = 1:1:size(Coefficients, 2)

Acc = Acc + x(TapIndex) * Coefficients(TapIndex);

end

FIR_Out(CurrentInputIndex) = Acc;

end

end

as it is shown in the result, I get the transistion band on half the value that it should be (should be 400-500) and also the double passband. I think I am missing something about the nyquist criteria but I don't know what. Fs is 2000 and I am only chirping to 1000Hz. what am I doing wrong here?

Answer 1

[Abhimenyu]

Hi Keyvan,

From the information shared by you, I could infer that you are trying to implement the designed filter using the "fdesign" function of MATLAB. The "Nyquist criterion" is correctly defined. The generation of the "chirp" signal is affecting the performance of the filter because the "InstantaneousFrequency" is not integrated to get the phase of the signal. The "cumtrapz" function of MATLAB can be used to integrate the "InstantaneousFrequency" to get the phase of the "chirp" signal, which is then used to generate the signal itself.

Please follow the example MATLAB R2023b code below to understand more:

%% NO CHANGE HAS BEEN MADE TO THE CODE IN THIS SECTION
clear screen
clearvars
format long

% Filter specifications
SamplingFrequency = 2000;
CenterFrequency = 450;
TransitionWidth = 100;
PassbandRipple = 1;
StopbandAttenuation = 80;

% Design the filter
DesignSpec = fdesign.lowpass('Fp,Fst,Ap,Ast',...
    CenterFrequency - TransitionWidth/2, ...
    CenterFrequency + TransitionWidth/2, ...
    PassbandRipple, StopbandAttenuation, ...
    SamplingFrequency);
LPF = design(DesignSpec, 'equiripple', 'SystemObject', true);
Coefficients = LPF.Numerator;

% Generate a chirp signal
TotalTime = 5; % seconds
StartFrequency = 5; % Hz
StopFrequency = 1000; % Hz

t = 0:1/SamplingFrequency:TotalTime;
InstantaneousFrequency = linspace(StartFrequency, StopFrequency, numel(t));

The "x" variable which is the "chirp" signal has been computed below using the MATLAB function "cumtrapz".

%% CHANGE MADE : GENERATION USING "cumtrapz" function
x = cos(2 * pi * cumtrapz(t, InstantaneousFrequency));

FIR_In = x;
% Apply the filter using the custom myFIR function
FIR_Out = myFIR(FIR_In, Coefficients);

figure
plot (InstantaneousFrequency,FIR_In(1,1:size(InstantaneousFrequency,2)),InstantaneousFrequency,FIR_Out(1,1:size(InstantaneousFrequency,2)))

%% NO CHANGE HAS BEEN MADE TO THE CODE IN THIS SECTION

function FIR_Out = myFIR(FIR_In, Coefficients)
    FIR_Out = zeros(1, size(FIR_In, 2)); % Initialize the output signal
    x = zeros(1, size(Coefficients, 2)); % Buffer for the filter's delay line

    % Process each sample of the input signal
    for CurrentInputIndex = 1:size(FIR_In, 2)
        % Shift the delay line to make room for the new sample
        for DelayLineIndex = size(Coefficients, 2):-1:2
            x(DelayLineIndex) = x(DelayLineIndex-1);
        end
        x(1) = FIR_In(CurrentInputIndex); % Insert the new sample at the beginning

        % Perform the dot product of the delay line and the filter coefficients
        Acc = 0;
        for TapIndex = 1:size(Coefficients, 2)
            Acc = Acc + x(TapIndex) * Coefficients(TapIndex);
        end
        FIR_Out(CurrentInputIndex) = Acc; % Store the result in the output signal
    end
end

For more information on the "cumtrapz" function, follow this MATLAB R2023b documentation link: https://www.mathworks.com/help/matlab/ref/cumtrapz.html

Regards,

Abhimenyu