Filter löschen
Filter löschen

How to remove the ripple or artifacts caused by lowpass filter?

48 Ansichten (letzte 30 Tage)
Hi, I am in a complex situation, where the lowpass fitering with sharp cutoff freqeuncy is affecting a technique I was developing. The following first code use a fc = 80Hz, to extract the envelop. As shown in the second code, I want to reduce the fc to 15Hz and achieve a smooth envelop.i.e. gradient A' reduces in the second simulation compared to first due to reduction of fc from 80 Hz to 15 Hz.
But, reducing to 15Hz, creates additional ripples and sharp edges and changes the expected envelop, which is a problem for me. So,Is there any other way to escape from these artifacts? I tried reducing the filter order to 1, but butter filtering is not performing well. So I am looking if there is any other way to elminate those added artifacts either by changing the lowpass filteration or by post proceesing the extracted envelop?
clear all;
% close all; clc;
%-------------------------------------------------------------
% System definition %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%-------------------------------------------------------------
fs=8192; % sampling frequency
dt = 1/fs; % sample time
T=3; % duration of the signal
Nt = T*fs; % total number of samples
t = 0:dt:T-dt; % time vector
% Source definition
f0 = 900; % initial angular speed in Hz
fT = 2700; % final angular speed in Hz
finst = linspace(f0,fT,Nt);
phi = 2*pi*cumsum(finst)*dt;
a1 = 45;
b1 = 5000;
c1 = 2;
c2 = 1.5;
d1 = 2+3*1j;
% Definition of envelop
int_phi0 = 3*pi/4;
A = d1+(a1 * exp(-b1 * (t - T/c1).^2 + 1j * int_phi0));
% Definition of source signal
q = A.'.*exp(1j*phi).';
filtorder =6; fc =80;%----------------------------------% fc-cutoff
[bcoeff,acoeff] = butter(filtorder,fc/fs*2);
Afilt = filtfilt(bcoeff,acoeff,q.*exp(-1j*phi).');
% figure()
% freqz(bcoeff,acoeff,[],fs)
figure()
plot(t,q,'k')
Warning: Imaginary parts of complex X and/or Y arguments ignored.
hold on
plot(t,abs(Afilt),'r','LineWidth',2)
grid on; box on;
legend('q(t)','A(t)')
xlabel('s')
ylabel('Amp')
title('q(t) and A(t) - filtered')
set(gca,'FontSize',15)
xlim([1.2 1.8])
subplot()
clear all;
% close all; clc;
%-------------------------------------------------------------
% System definition %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%-------------------------------------------------------------
fs=8192; % sampling frequency
dt = 1/fs; % sample time
T=3; % duration of the signal
Nt = T*fs; % total number of samples
t = 0:dt:T-dt; % time vector
% Source definition
f0 = 900; % initial angular speed in Hz
fT = 2700; % final angular speed in Hz
finst = linspace(f0,fT,Nt);
phi = 2*pi*cumsum(finst)*dt;
a1 = 45;
b1 = 5000;
c1 = 2;
c2 = 1.5;
d1 = 2+3*1j;
% Definition of envelop
int_phi0 = 3*pi/4;
A = d1+(a1 * exp(-b1 * (t - T/c1).^2 + 1j * int_phi0));
% Definition of source signal
q = A.'.*exp(1j*phi).';
filtorder =6; fc =15;%----------------------------------% fc-cutoff
[bcoeff,acoeff] = butter(filtorder,fc/fs*2);
Afilt = filtfilt(bcoeff,acoeff,q.*exp(-1j*phi).');
% figure()
% freqz(bcoeff,acoeff,[],fs)
figure()
plot(t,q,'k')
Warning: Imaginary parts of complex X and/or Y arguments ignored.
hold on
plot(t,abs(Afilt),'r','LineWidth',2)
grid on; box on;
legend('q(t)','A(t)')
xlabel('s')
ylabel('Amp')
title('q(t) and A(t) - filtered')
set(gca,'FontSize',15)
xlim([1.2 1.8])

Akzeptierte Antwort

Mathieu NOE
Mathieu NOE am 6 Jul. 2023
hello
Butterworth filters do have ripples in their step response...
look at the step response of your Butterworth filter with fc = 80 Hz :
(x axis is samples not seconds as displayed !)
now that you have lowered the fc to 15 Hz, it's getting very long until the filter's output is steady again :
no surprise , it takes 5 time longer and with oscillations
so you should maybe use Bessel filters instead (step response do not exhibit any oscillation)
this is my suggestion below , you can change filtorder and fc as you wish, this will affect how close the envelop will be to the "real" one, but the ripples are gone
% clear all;
% close all; clc;
%-------------------------------------------------------------
% System definition %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%-------------------------------------------------------------
fs=8192; % sampling frequency
dt = 1/fs; % sample time
T=3; % duration of the signal
Nt = T*fs; % total number of samples
t = 0:dt:T-dt; % time vector
% Source definition
f0 = 900; % initial angular speed in Hz
fT = 2700; % final angular speed in Hz
finst = linspace(f0,fT,Nt);
phi = 2*pi*cumsum(finst)*dt;
a1 = 45;
b1 = 5000;
c1 = 2;
c2 = 1.5;
d1 = 2+3*1j;
% Definition of envelop
int_phi0 = 3*pi/4;
A = d1+(a1 * exp(-b1 * (t - T/c1).^2 + 1j * int_phi0));
% Definition of source signal
q = A.'.*exp(1j*phi).';
% filtorder =6; fc =80;%----------------------------------% fc-cutoff
% [bcoeff,acoeff] = butter(filtorder,fc/fs*2);
% alternative with Bessel filter
filtorder =3; fc =500;%----------------------------------% fc-cutoff
[b,a] = besself(filtorder,fc,'low'); % NB it's analog !
[bcoeff,acoeff] = c2dm(b,a,1/fs,'tustin');
Afilt = filtfilt(bcoeff,acoeff,q.*exp(-1j*phi).');
% figure()
% freqz(bcoeff,acoeff,[],fs)
figure()
dstep(bcoeff,acoeff)
figure()
plot(t,q,'k')
hold on
plot(t,abs(Afilt),'r','LineWidth',2)
grid on; box on;
legend('q(t)','A(t)')
xlabel('s')
ylabel('Amp')
title('q(t) and A(t) - filtered')
set(gca,'FontSize',15)
xlim([1.2 1.8])
%
% subplot()
clear all;
% close all; clc;
%-------------------------------------------------------------
% System definition %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%-------------------------------------------------------------
fs=8192; % sampling frequency
dt = 1/fs; % sample time
T=3; % duration of the signal
Nt = T*fs; % total number of samples
t = 0:dt:T-dt; % time vector
% Source definition
f0 = 900; % initial angular speed in Hz
fT = 2700; % final angular speed in Hz
finst = linspace(f0,fT,Nt);
phi = 2*pi*cumsum(finst)*dt;
a1 = 45;
b1 = 5000;
c1 = 2;
c2 = 1.5;
d1 = 2+3*1j;
% Definition of envelop
int_phi0 = 3*pi/4;
A = d1+(a1 * exp(-b1 * (t - T/c1).^2 + 1j * int_phi0));
% Definition of source signal
q = A.'.*exp(1j*phi).';
% filtorder =6; fc =15;%----------------------------------% fc-cutoff
% [bcoeff,acoeff] = butter(filtorder,fc/fs*2);
% alternative with Bessel filter
filtorder =3; fc =100;%----------------------------------% fc-cutoff
[b,a] = besself(filtorder,fc,'low'); % NB it's analog !
[bcoeff,acoeff] = c2dm(b,a,1/fs,'tustin');
Afilt = filtfilt(bcoeff,acoeff,q.*exp(-1j*phi).');
% figure()
% freqz(bcoeff,acoeff,[],fs)
figure()
dstep(bcoeff,acoeff)
figure()
plot(t,q,'k')
hold on
plot(t,abs(Afilt),'r','LineWidth',2)
grid on; box on;
legend('q(t)','A(t)')
xlabel('s')
ylabel('Amp')
title('q(t) and A(t) - filtered')
set(gca,'FontSize',15)
xlim([1.2 1.8])
  5 Kommentare
Star Strider
Star Strider am 6 Jul. 2023
Bessel filters cannot be represented as digital filters, and retain their phase-neutral characteristics.
Mathieu NOE
Mathieu NOE am 6 Jul. 2023
I don't see what is the trouble with the Bode plot of the bessel filter
(also done now in your code below)
at f = fc the filter has a -6 dB attenuation
also this is your code , now corrected . notice that with the 2pi corrction we finally are using the same fc as in your original code . It's just now bessel instead of butter and the filter order has been reduced from 6 to 3
clear all;
close all; clc;
%-------------------------------------------------------------
% System definition %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%-------------------------------------------------------------
fs=8192; % sampling frequency
dt = 1/fs; % sample time
T=3; % duration of the signal
Nt = T*fs; % total number of samples
t = 0:dt:T-dt; % time vector
% Source definition
f0 = 900; % initial angular speed in Hz
fT = 2700; % final angular speed in Hz
finst = linspace(f0,fT,Nt);
phi = 2*pi*cumsum(finst)*dt;
a1 = 45;
b1 = 5000;
c1 = 2;
c2 = 1.5;
d1 = 2+3*1j;
% Definition of envelop
int_phi0 = 3*pi/4;
A = d1+(a1 * exp(-b1 * (t - T/c1).^2 + 1j * int_phi0));
% Definition of source signal
q = A.'.*exp(1j*phi).';
% filtorder =6; fc =80;%----------------------------------% fc-cutoff
% [bcoeff,acoeff] = butter(filtorder,fc/fs*2);
% alternative with Bessel filter
filtorder =3; fc =80;%----------------------------------% fc-cutoff
[b,a] = besself(filtorder,2*pi*fc,'low'); % NB it's analog !
[bcoeff,acoeff] = c2dm(b,a,1/fs,'tustin');
Afilt = filtfilt(bcoeff,acoeff,q.*exp(-1j*phi).');
figure()
f = logspace(1,3,300);
freqz(bcoeff,acoeff,f,fs)
figure()
dstep(bcoeff,acoeff)
figure()
plot(t,q,'k')
hold on
plot(t,abs(Afilt),'r','LineWidth',2)
grid on; box on;
legend('q(t)','A(t)')
xlabel('s')
ylabel('Amp')
title('q(t) and A(t) - filtered')
set(gca,'FontSize',15)
xlim([1.2 1.8])
%
% subplot()
clear all;
% close all; clc;
%-------------------------------------------------------------
% System definition %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%-------------------------------------------------------------
fs=8192; % sampling frequency
dt = 1/fs; % sample time
T=3; % duration of the signal
Nt = T*fs; % total number of samples
t = 0:dt:T-dt; % time vector
% Source definition
f0 = 900; % initial angular speed in Hz
fT = 2700; % final angular speed in Hz
finst = linspace(f0,fT,Nt);
phi = 2*pi*cumsum(finst)*dt;
a1 = 45;
b1 = 5000;
c1 = 2;
c2 = 1.5;
d1 = 2+3*1j;
% Definition of envelop
int_phi0 = 3*pi/4;
A = d1+(a1 * exp(-b1 * (t - T/c1).^2 + 1j * int_phi0));
% Definition of source signal
q = A.'.*exp(1j*phi).';
% filtorder =6; fc =15;%----------------------------------% fc-cutoff
% [bcoeff,acoeff] = butter(filtorder,fc/fs*2);
% alternative with Bessel filter
filtorder =3; fc =15;%----------------------------------% fc-cutoff
[b,a] = besself(filtorder,2*pi*fc,'low'); % NB it's analog !
[bcoeff,acoeff] = c2dm(b,a,1/fs,'tustin');
Afilt = filtfilt(bcoeff,acoeff,q.*exp(-1j*phi).');
figure()
f = logspace(0,2,100);
freqz(bcoeff,acoeff,f,fs)
figure()
dstep(bcoeff,acoeff)
figure()
plot(t,q,'k')
hold on
plot(t,abs(Afilt),'r','LineWidth',2)
grid on; box on;
legend('q(t)','A(t)')
xlabel('s')
ylabel('Amp')
title('q(t) and A(t) - filtered')
set(gca,'FontSize',15)
xlim([1.2 1.8])

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (1)

Walter Roberson
Walter Roberson am 6 Jul. 2023
Imagine that you were doing a discrete fourier transform of a square wave. Maybe it looks not bad at a target frequency. But then you try to use it at a higher frequency, and you look carefully at the leading and trailing edges and see ripples. Would decreasing the number of fourier coefficients solve the problem? Obviously not. Would increasing the number of fourier coefficients solve the problem? No -- examine the result at a sufficiently high frequency and you would continue to see ripples at the leading and trailing edges.
Using butter filter has similar problems; it is a discrete reconstruction of an ideal filter shape, and if you use too few coefficients the result is going to look bad. You can increase the filter order to make the ripples less obvious, but you cannot eliminate them with any finite order. The best you can do is reduce the magnitude of the ripple to less than a level that is acceptable to you -- which might require increasing the order (which in turn implies increased delays). buttord might help in the design.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by