0

I want to get Velocity from Acceleration data(10 Hz) in MATLAB. I have the real Constant Velocity that I've calculated by (DeltaX/DeltaT). The problem is that when I use the code below, the graph has lot's of differences with the Real Velocity Graph. I've used both cumsum and omega arithmetic, could someone tell me if there is any problem in my codes of any of both methods? Thank you in advance. link to access txt files: https://spaces.hightail.com/space/19rl7

%%%%%%%%%%%%%%% cumsum METHOD %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fid1 = fopen( 'Acceleration_data.txt', 'r' );
fid2 = fopen( 'Actual_Velocity_data.txt', 'r' );
formatSpec = '%f';
vel_Fast_Running = fscanf(fid2,formatSpec);
acc_Fast_Running = fscanf(fid1,formatSpec);
Fs=1000;
Ns = 10;%% Sampling freq
T = 1/Fs; %% Period
L = length(acc_Fast_Running); %% signal length
acceleration_data=acc_Fast_Running-median(acc_Fast_Running);
[B,A] = butter(6,0.098,'high'); %% highpass filter
a_filtr = filtfilt(B,A,acc_Fast_Running); %% acceleratiof w/o DC
v_cumsum = cumsum(a_filtr)*0.1;%integration
v_cumsum = v_cumsum - mean(v_cumsum);
[B,A] = butter(3,0.12,'low'); 
v_cumsum_filtr = filtfilt(B,A,v_cumsum);
%subplot(1,2,2);plot([v_cumsum' v_cumsum_filtr']);
hold on;plot(vel_Fast_Running)
plot(v_cumsum_filtr)
legend('Actual Vel','Vel From integration')
xlabel('Time')
ylabel('Velocity m/s')

%%%%%%%%%%%%%%% OMEGA ARITHMETIC METHOD %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fid1 = fopen( 'Acceleration_data.txt', 'r' );
fid2 = fopen( 'Actual_Velocity_data.txt', 'r' );
formatSpec = '%f';
vel_Fast_Running = fscanf(fid2,formatSpec);
acc_Fast_Running = fscanf(fid1,formatSpec);
Fs=1000;
Ns = 10;%% Sampling freq
T = 1/Fs; %% Period
L = length(acc_Fast_Running ); %% signal length
%t=(1:L)/Fs;
f = Fs*(0:L-1)/L;
acc_Fast_Running =acc_Fast_Running -mean(acc_Fast_Running );
[B,A] = butter(6,0.098,'high'); %% highpass filter
a_filtr = filtfilt(B,A,acc_Fast_Running);
Y = (fft(a_filtr)); % Fast fourier Transform
for i = 1 : L 
   if f(i) ~= 0 && f(i) < Fs/2;
       Yf(i) = Y(i)/(2*pi*f(i)*sqrt(-1)) ;
   else  
       Yf(i) = 0;
   end
end
Yf(1:100)=0; %% 
d_oa = real(ifft(Yf));
d1 = designfilt('lowpassiir','FilterOrder',4, ...
 'HalfPowerFrequency',0.1,'DesignMethod','butter');
y = filtfilt(d1,d_oa);
figure(15)
%plot( d_oa)
subplot(2,1,1)
plot(y)
xlabel('Time')
ylabel('Vel be Omega Method')
subplot(2,1,2)
plot(vel_Fast_Running )
xlabel('Time')
ylabel('Actual Velocity m/s')

enter image description here

Shayan Shiri
  • 51
  • 1
  • 6

0 Answers0