Filtering FIRCLS1 and evaluate kernel and its power spectrum

Question

Jonas Bender 2022-3-15

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1672059-filtering-fircls1-and-evaluate-kernel-and-its-power-spectrum

回答： Jonas Bender 2022-3-15

Dear all,

I would like to filter my data using FIRCLS1. Unfortunately I did not get out filtered data. Filter properties are:

order = round (8*srate/cutoff); % length of the filter kernel in the time domain 
cutoff_cls1 = 0.3; % normalized cutoff value corresponding to nyquist
pb_ripple = 0.0002; % upper bound: bandpass (PB) ripple (Durchlassbereich: oberhalb der Grenzfrequenz von 3Hz) 
sb_ripple = 0.0002; % lower bound: bandstop (SB) ripple (Bandsperre unterhalb der Grenzfrequenz von 3Hz))

Ideal filter response:

shape_cls = [1 1 0 0]; % y-axis
frex_cls = [0, cutoff_cls1 * 10 , cutoff_cls1 * 10, srate/2] / (srate/2); % z-axis

Generate filter kernel

filtkern_cls1 = fircls1 (order,cutoff_cls1, pb_ripple, sb_ripple);

FFT to inspect the power spectrum of the kernel

hz_filtkern_cls1 = linspace (0, (srate/2), (pnts_first_last (1,:)/2) + 1); 
filtkernX_cls1 = abs (fft (filtkern_cls1, pnts_first_last (1,:))); 

Some plotting

% % evaluate filter kernel and its power spectrum
% 
figure (7), clf
subplot (221)
plot ((-order / 2 : order/2) /srate, filtkern_cls1, 'k', 'linew', .5);
xlabel ('time [sec]')
title ('filter kernel (fircls1)')
% 
% % plot the ideal filter
% 
figure (7), 
subplot (222); hold on;
plot (frex_cls * srate/2, shape_cls, 'r', 'linew', 1)
set (gca, 'xlim', [0 6])
set (gca, 'ylim', [0 2])
% 
% % plot the power spectrum of the filter kernel
% 
plot (hz_filtkern_cls1, filtkernX_cls1 (1:length (hz_filtkern_cls1)), 'k')
xlabel ('frequency (Hz)')
ylabel ('gain')
title ('filter kernel spectrum')
hold off;
% 

Apply the filter to the signal

for i=1:numel(data)
    data(i).filt_first_last_cls1 = filtfilt(filtkern_cls1, 1, data(i).acc_pelvis_first_last);
end

further plotting

subplot (2, 2, [3 4])
plot (data(1).timevec_first_last, data(1).filt_first_last_cls1)
set (gca, 'xlim', [200 215])
xlabel ('time (sec)')
ylabel ('amplitude')
title ('filtered acceleration signal pair1_FIRCLS1')

The result is not I would expect. The power spectrum of the kernel ends at 300 Hz instead of 3 Hz (see black line; ther red one is the the ideal filter response. Subsequently data is not filterd as prefered: