Clear Filters
Clear Filters

Drift removal of velocity obtained by integrating acceleration data

19 views (last 30 days)
Hi, can you help me?
I'm taking acceleration data from the IMU sensor and getting a velocity graph through integration and then looking at the velocity of motion by vector sum.
By the way, even though the drift has occurred and the movement is over, it remains at the speed of a point in time, not zero, is there any way to solve this?
I adjusted the offset because there was no movement before 0.47 seconds, and the blue line is the same, but based on the red line, there is no movement after 2.176 seconds (movement ends almost simultaneously with left and right)..
I used "detrend" to correct the trend, but as I did that, the graph tilted and transformed into a motionless section
Do you know how to solve this problem? Any help is appreciated
%% Data preprocessing
% --- left hand ---
T_left = readtable('ALT10.csv', 'VariableNamingRule', 'preserve');
timestamp_left = T_left.uni_time(2:end);
time_left = timestamp_left
time_left = timestamp_left-timestamp_left(1)
xacc_left = T_left.("watch_acc.y")(2:end);
yacc_left = T_left.("watch_acc.y")(2:end);
zacc_left = T_left.("watch_acc.y")(2:end);
% --- right hand ---
T_right = readtable('ART10.csv', 'VariableNamingRule', 'preserve');
timestamp_right = T_right.uni_time(2:end);
time_right = timestamp_right
time_right = timestamp_right-timestamp_right(1)
xacc_right = T_right.("watch_acc.y")(2:end);
yacc_right = T_right.("watch_acc.y")(2:end);
zacc_right = T_right.("watch_acc.y")(2:end);
%% Offset correction and Filtering
stationary_indices = 1:30; % Indices of samples where the device is assumed stationary
% Offset correction for all axes
offset_x_left = mean(xacc_left(stationary_indices));
xacc_left_corrected = xacc_left - offset_x_left;
offset_y_left = mean(yacc_left(stationary_indices));
yacc_left_corrected = yacc_left - offset_y_left;
offset_z_left = mean(zacc_left(stationary_indices));
zacc_left_corrected = zacc_left - offset_z_left;
% Repeat offset correction for the right hand
offset_x_right = mean(xacc_right(stationary_indices));
xacc_right_corrected = xacc_right - offset_x_right;
offset_y_right = mean(yacc_right(stationary_indices));
yacc_right_corrected = yacc_right - offset_y_right;
offset_z_right = mean(zacc_right(stationary_indices));
zacc_right_corrected = zacc_right - offset_z_right;
%% Velocity integration for vector magnitude
xvelocity_left_mag = getvelocity(xacc_left_corrected,time_left);
yvelocity_left_mag = getvelocity(yacc_left_corrected,time_left);
zvelocity_left_mag = getvelocity(zacc_left_corrected,time_left);
xvelocity_right_mag = getvelocity(xacc_right_corrected,time_right);
yvelocity_right_mag = getvelocity(yacc_right_corrected,time_right);
zvelocity_right_mag = getvelocity(zacc_right_corrected,time_right);
%for you information, getvelocity is this..
%function v = getvelocity(acc, time)
% v = zeros(1, length(acc));
% for i = 1:length(acc)-1
% dt = time(i+1) - time(i);
% numerical integration
% v(i+1) = v(i)+0.5*dt*(acc(i) + acc(i+1));
% end
fs = 1000
cutoff_freq = 5; % 보다 낮은 컷오프 주파수 설정
xlp_velocity_left = lowpass(xvelocity_left_mag, cutoff_freq, fs);
ylp_velocity_left = lowpass(yvelocity_left_mag, cutoff_freq, fs);
zlp_velocity_left = lowpass(zvelocity_left_mag, cutoff_freq, fs);
xvelocityR= lowpass(xvelocity_right_mag, cutoff_freq, fs);
yvelocityR = lowpass(yvelocity_right_mag, cutoff_freq, fs);
zvelocityR = lowpass(zvelocity_right_mag, cutoff_freq, fs);
vectorsumR = sqrt(xvelocityR.^2 + yvelocityR.^2 + zvelocityR.^2);
vectorsumL = sqrt(xlp_velocity_left.^2 + ylp_velocity_left.^2 + zlp_velocity_left.^2)
time_left = (timestamp_left - timestamp_left(1)) / 1000; % change to second(S) frome unix_time
time_right = (timestamp_right - timestamp_right(1)) / 1000;

Answers (1)

Walter Roberson
Walter Roberson on 21 Nov 2023
It turns out that even the best single-mode IMU systems encounter substantial drift
"Even the best accelerometers, with a standard error of 10 micro-g, would accumulate a 50-meter (164-ft) error within 17 minutes"
Your data only has about 3 significant figures, 4 at best, and your labels talk about mm; I figure that your error is probably on the order of 10 milli-g rather than 10 micro-g, so I estimate that you would probably be about 50 m off in... slightly over 1 second. If your actual error is 1 milli-g then 50 m in 10.2 seconds I estimate, and over your 2.176 seconds... you should still expect over 10 meters positioning error.


Community Treasure Hunt

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

Start Hunting!