Use the correct filter on a signal

Hello everyone,

I am looking to filter the following signal in such a way that all vibrations after 1.2 seconds disappear. I'm not sure what type of filter to apply. To filter out the noise, I had used a Butterworth low-pass filter. In the first signal, I actually want to obtain only 9 peaks using the `findpeaks` function.

I would like to get something like that with no vibration after 1.2 seconds :

You can find the text file with the signal data attached (VIA = f(t)).

Below is my current code :

(make sure to change the pathfile, if you want to help me)

% Programme "Banc_Essai" Determine si le signal est Amorti ou Non Amorti
clear all;
clc ;
close all;

%% Choix de l'essai
% Nom des cases
prompt = {'\fontsize{9}Série :','\fontsize{9}Numéro :','\fontsize{9}Essai :'}; 
% Nom de la boîte de dialogue
dlgtitle = 'Choix du Dossier';   
 % Dimension des zones de texte
dims = [1 40;1 40;1 40];   
% Entrées prédéfinies
definput = {'5_23_019','05031610','ACQ001.txt'};  
opts.Interpreter = 'tex'; 
% Réponse sous forme matrice
boite = inputdlg(prompt,dlgtitle,dims,definput,opts); 

% 27:02271611  33: 04261043 34: 04271505 & 04281125 35: 05021455 36: 05031137

Serie = boite{1}; % 1ère donnée Dialogue Box
Numero = boite{2}; % 2ème donnée Dialogue Box
Essai = boite{3}; % 3ème donnée Dialogue Box

% Création chemin d'accès
txt_nom = sprintf('U:\\theo\\%s\\%s',Serie,Numero);
fname =fullfile(txt_nom,Essai);

txt_titre = sprintf('%s\\%s\\%s',Serie,Numero,Essai);

%% Lecture du fichier
A=dlmread(fname,'',1); % lit le fichier et supprime la 1r ligne.

% Récupération des colonnes temps et VIA
temps = A(:,1);
VIA = A(:,2);
VIAm = -VIA; % Signal opposé car findpeaks trouve les pics "max" avec cette méthode on détermnie les pics "min"

%% FIltrage du signal

order = 4;  % Ordre du filtre %4
cutoffFreq =30;  % Fréquence de coupure du filtre (en Hz) % 50
fs = 1000;
[b, a] = butter(order, cutoffFreq/(fs/1), 'low');  % Conception du filtre passe-bas Butterworth %/2

y = filter(b, a, VIA);  % Application du filtre au signal
y_2 = -y; % Signal inversé

%% Localisation des pics 

Fs = 1000 ;%Hz car pas = 0.001s
balayage = 0.075 ;%s

% Pics au dessus de 5 en prenant le signal inverse
[Y0 , X0] =findpeaks(VIAm,Fs,'MinPeakDistance',balayage,'MinPeakHeight',5); 
Y0 = -Y0;
sY0= length(Y0);

% Pics au dessus de 5 en prenant le signal filtré
[Y , X] =findpeaks(y_2,Fs,'MinPeakDistance',balayage,'MinPeakHeight',5); 
Y = -Y; % Vecteur inverse 

% Pics au dessus de 2 en prenant le signal inverse
[Y2, X2]=findpeaks(VIAm,Fs,'MinPeakDistance',balayage,'MinPeakHeight',2); 

lX = length(X); % Nb de pics pour les calculs de periode, fréquence et A.

txt_nb = sprintf('Nb pics : %d',lX);

if lX <= 10

        % Calcul Average Frequency
        T = 0;
        
        for i = 1:lX-1
        
        T = (X(i+1) - X(i)) + T;
        
        end
        
        Tot = T/(lX-1);
        F =  1/Tot;
        
        % Calcul de A
        A = 0;
        
        for j =1:lX-1
        
            A = A + (Y(j)/Y(j+1))/(lX-1);
        end
        
        J = 0.028;%kg.m²
        Rp = 0.017; % Amortissement banc
        R = 2*J*F*log(A) - Rp;

        txt_A = "AMORTI";

 elseif lX > 10

     condition = 0;
     for k = 1:5

        P = Y(k)/Y(lX+k-5);
        if P <= 2 
            txt_A = "NON AMORTI";

             T = 0;
        
        for i = 1:lX-1
        
        T = (X(i+1) - X(i)) + T;
        
        end
        
        Tot = T/(lX-1);
        F =  1/Tot;
        
        % Calcul de A
        A = 0;
        
        for j =1:lX-1
        
            A = A + (Y(j)/Y(j+1))/(lX-1);
        end
        
        J = 0.028;%kg.m²
        Rp = 0.017; % Amortissement banc
        R = 2*J*F*log(A) - Rp;

            condition = 1;
            break
        end 
     end

    if condition == 0

        T = 0;
        
        for i = 1:lX-1
        
        T = (X(i+1) - X(i)) + T;
        
        end
        
        Tot = T/(lX-1);
        F =  1/Tot;
        
        % Calcul de A
        A = 0;
        
        for j =1:lX-1
        
            A = A + (Y(j)/Y(j+1))/(lX-1);
        end
        
        J = 0.028;%kg.m²
        Rp = 0.017; % Amortissement banc
        R = 2*J*F*log(A) - Rp;

        txt_A = "AMORTI";

    end

end

limx = max(temps); %limites en fonction des valeurs max des données
limy = max(VIA) + 6;

txt_R = sprintf("Coeff amortissement : %.3f",R);

%% Affichage de la courbe 
Fig = figure( 'Position', [300 100 800 800]);
t = tiledlayout(2,1,'TileSpacing','Compact','Padding','Compact');
title(t,txt_titre,'FontSize',13,'FontWeight','normal');

nexttile (1)
hold on
plot(temps, VIA,'LineWidth',1.5);
plot(X0,Y0,'og',X0,Y0,'g')
xlabel('Temps (s)');
ylabel('Vitesse angulaire (rad/s)');
title('Signal d''origine');
xlim([0 limx])
ylim([-limy limy])
grid on;
NE = [max(xlim) max(ylim)]-[diff(xlim) -diff(ylim)]*0.0001;
SE = [max(xlim) min(ylim)]+[-diff(xlim) diff(ylim)]*0.05;
text(NE(1), NE(2), txt_A, 'VerticalAlignment','top', 'HorizontalAlignment','right',"FontSize",15,'Color','r')



nexttile (2)
hold on
plot(temps, y,'LineWidth',1.5);
plot(X,Y,'color',[0.9, 0.5, 0.1],'LineWidth',2) 
plot(X,Y,'rx', 'MarkerSize', 10, 'LineWidth', 2)
hold off
xlabel('Temps (s)');
ylabel('Vitesse angulaire (rad/s)'); 
title('Signal filtré');
xlim([0 limx])
ylim([-limy limy])
grid on;
text(SE(1), SE(2), txt_nb, 'VerticalAlignment','bottom', 'HorizontalAlignment','right','FontSize',12)
text(NE(1), NE(2), txt_R, 'VerticalAlignment','top', 'HorizontalAlignment','right',"FontSize",12)

Thank you!

Answer 1

[Star Strider]

Frequency-selective filters do not work that way. (Intrigued by this, I attenpted to design one by taking the z-transform of the attenuation function and failed to get the desired result when I applied it.) They only woirk in the frequency domain, not time domain, even though digital (discrete) filters work in the time domain.

It is possible to do what you want by multiplying an exponential function by your signal. I designed one that will deliver the desired attenuation in decibells at 1.2 seconds.

Try this —

T1 = readtable('ACQ001.txt')

T1 = 4008×2 table
    temps     VIA  
    _____    ______

        0     1.329
    0.001     0.307
    0.002     0.102
    0.003     0.102
    0.004     0.102
    0.005     0.102
    0.006     0.102
    0.007     0.613
    0.008    -0.103
    0.009         0
     0.01         0
    0.011    -0.103
    0.012    -0.103
    0.013         0
    0.014    -0.103
    0.015         0

VN = T1.Properties.VariableNames;
t = T1.temps;
VIA = T1.VIA;

Fs = 1/mean(diff(t));

figure
plot(t, VIA)
grid
xlabel(VN{1})
ylabel(VN{2})
title('Original Signal')

dB = 15;                                                                    % Desired Attenuation (dB) At 1.2 s
tau = @(dB) log(db2mag(-dB))/1.2;                                           % Calculate Parameter
efcn = @(dB,t) exp(tau(dB)*t);                                              % Exponential Attenuation Function

figure
plot(t, efcn(dB,t))
grid
xlabel(VN{1})
ylabel('Attenuation')
title(sprintf('Attenuation %2d dB at 1.2 s (\\tau = %.3f s^{-1})', dB, tau(dB)))

figure
plot(t, VIA.*efcn(dB,t))
grid
xlabel(VN{1})
ylabel(VN{2})
title('Attenuated Signal')

This is of course not a frequency-domain filter.

.

[Théo]

Thanks a lot, I will try this immediatly. ;)

[Star Strider]

My pleasure!