ISO 2631 filter not working

Hi everybody,

I'm trying to implement an ISO 2631 filter so as to weight a mechanical vibration, but it is returning NaNs and huge numbers. I'm no vibration expert, so I have no idea of what could be going on. The input is an acceleration vector 'a' and the output should be the weighted acceleration vector 'aw'.

if true
a = rand(100,1);
f1 = 0.4;  % Wk
f2 = 100;  % Wk
f3 = 12.5; % Wk 
f4 = 12.5; % Wk 
f5 = 2.37; % Wk
f6 = 3.35; % Wk
Q4 = 0.63; % Wk  
Q5 = 0.91; % Wk  
Q6 = 0.91; % Wk
w1 = 2*pi*f1; 
w2 = 2*pi*f2; 
w3 = 2*pi*f3; 
w4 = 2*pi*f4; 
w5 = 2*pi*f5; 
w6 = 2*pi*f6;
s = tf('s');
Hh = 1 / (1 + sqrt(2)*w1/s + (w1/s)^2);
Hl = 1 / (1 + sqrt(2)*s/w2 + (s/w2)^2);
Ht = (1 + s/w3) / (1 + s/(Q4*w4) + (s/w4)^2);
Hs = (1 + s/(Q5*w5) + (s/w5)^2) * (w5/w6)^2 / (1 + s/(Q6*w6) + (s/w6)^2);
[num, den] = tfdata(Hh*Hl*Ht*Hs, 'v'); % return as vectors
aw = filter(num, den, a); % filtered acceleration
end

Any hint would be much appreciated.

Thanks, Alfonso

[Romina Cornistein]

Try this code:

function y = filterwk(x, fs)                                    

f1 = 0.4;                                                       
f2 = 100;                                                       
f3 = 12.5;                                                      
f4 = 12.5;                                                      
Q4 = 0.63;                                                      
f5 = 2.37;                                                      
Q5 = 0.91;                                                      
f6 = 3.35;                                                      
Q6 = 0.91;                                                      

w1 = 2*pi*f1; 
w2 = 2*pi*f2; 
w3 = 2*pi*f3;                                                   
w4 = 2*pi*f4;                                                   
w5 = 2*pi*f5;                                                   
w6 = 2*pi*f6;                                                   
Fs=fs;                                                                

A1=[1 0 0];
    B1=[1 2^(0.5)*w1 w1^2];

      [b1, a1] = bilinear(A1, B1, Fs);
   num1=[1];
      den1=[ (1/w2)^2 2^(0.5)/w2 1];

[b2, a2] = bilinear(num1, den1, Fs);

% a-v transition                                                
%                                                               
B3 = [1/w3 1];                                                  
A3 = [1/(w4^2) 1/(Q4*w4) 1];                                       
[b3,a3] = bilinear(B3,A3,fs);                                   
%                                                               
% upward step                                                   
%                                                               
B4 = [1/(w5^2) 1/(Q5*w5) 1]*((w5^2)/(w6^2));                           
A4 = [1/(w6^2) 1/(Q6*w6) 1];                                       
[b4,a4] = bilinear(B4,A4,fs);                                   
%                                                               
                                                     %                                                               
y = filter(b2,a2,x);                                            
y = filter(b1,a1,y);                                            
y = filter(b3,a3,y);                                            
y = filter(b4,a4,y);

Answer 1

[David Camacho Alcocer]

Hi Alfonso,

Did the solution provided work?

Also, was your data obtained trough measurements or was simulated? I asked this because I do not get why the analog to digital conversion suggested by Romina. In my case the data is already (if I am not mistaking) in digital form. However, when I apply the filters, I also get huge numbers.

So, please let us know if the solution suggested worked for you and if not why. This would be extremely helpful

Thanks