Integrate Acceleration to velocity, wrong plot

I try to use my Smartphonesensor to check if the theoretical speed of the Roboter ist the same like the real speed.
So I installed the Matlab App and used the acceleration sensor to check it. Because i didn´t mounted the smartphone totally zero in every direction (x,y,z) i needed to subtract an offset. I took the offset from the resting time bevore the movement, here i took the values from 300-500 (100hz) and calculated the average. Because i had sensor noice i also added a gauß filter. At the End i had two plots one in x-direction the over in y-direction.
So the plot for the Acceleration in the x-Direction looks very good and seems to be right.
And the Velocity of the x-direction, which i calculated with this code, doesnt look right:
velocity_y = cumtrapz(t_LP,y_ac');
plot(t_LP,velocity_y,'g')
How you can see there is an rising chart, which in the end has an velocity of 0.5 m/s but i the velocity should be zero.
I think the offset is not totally exact and after the integration i get an big mistake so the chart ist frequently rising. Does anyone knows how to erease this mistake?

5 个评论

hello
I would suggest you high pass filter the acceleration before cumtrapz
of course you will have to tweak the filter order (start low like N= 2) and cut off frequency to not remove to much of the useful acceleration data.
At first thanks for your fast anwser. I have an further question, do you have an example how to design the high pass filter on my values and does i have to erase the gauß filter? Im would be realy thankfull for an anweser!
I have now made an passfilter with the order 2, but I dont understand what i need to change in the passband frequency, ripple and samplerate.
hpFilt = designfilt('highpassiir','FilterOrder',2, ...
'PassbandFrequency'10,'PassbandRipple',0.1, ...
'SampleRate',100);
dataOut = filter(hpFilt,x_ac);
plot(t_LP, dataOut)
After it i do cumtrapz, the Velocity Diagramm looks like the Acceleration one. Why ist it so?
velocity_x = cumtrapz(t_LP,dataOut);
plot(t_LP,velocity_x,'g')
Thanks for help
hello
I will post my suggestion in the answer section - I already has someone asking me for help
Your velocity graph looks ok to me - what did you expect ?
you can do the following test : derivation of the velocity (using diff or gradient) and you'll fall back on your initial acceleration data.

请先登录,再进行评论。

回答(2 个)

so finally... try this / adapt it to your own needs .
I can further help you if needed
FYI : this code does two integrations, to get velocity and displacement (this may not really interest you, but I leave it just in case)
clc
clearvars
load('Acceleration.mat');
% Acceleration data
acc = acc*9.81; % Add gravity effect
acc = detrend(acc);
% Time
tStep = 0.00048828; % Length of each time step
N = length(acc)*tStep;
t = 0:tStep:N;
t(end) = [];
N = length(t);
dt = mean(diff(t)); % Average dt
fs = 1/dt; % Frequency [Hz] or sampling rate
% some additionnal high pass filtering
N = 2;
fc = 0.5; % Hz
[B,A] = butter(N,2*fc/fs,'high');
acc2 = filter(B,A,acc);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
velocity = cumtrapz(dt,acc2);
velocity = detrend(velocity);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
disp = cumtrapz(dt,velocity);
%
figure(1)
plot(t,acc);
figure(2)
plot(t,velocity);
figure(3)
plot(t,disp);

6 个评论

Thank you very much for your help! Im sorry but I have one more question. I know that in the timeline between 11 seconds and 30 seconds i have a constant speed of ~50mm/s but the diagramm shows sometimes some peaks but not a constant drive.
So it looks a bit different which i had bevor, there you can see nicely the nearly constant speed:
Can you give me advise how i can accomplishe something like this with your skript.
well , can you share the new code with your data , so I can better help you ?
Here i just added the test.1 data and i noticed, that in the y-direction i need has to have a speed from constant ~50mm/s.
But how you can see in the start i reach the 50 mm/s but at ~15 seconds the chart falls and i the value stays nearly 0.
clc
clearvars
load('test1.mat');
x_ac = Acceleration.X; %Beschl. x-Achse
y_ac = Acceleration.Y; %Beschl. x-Achse
z_ac = Acceleration.Z; %Beschl. x-Achse
% Acceleration data
acc = y_ac*9.81; % Add gravity effect
acc = detrend(acc);
% Time
tStep = 0.02; % Length of each time step
N = length(acc)*tStep;
t = 0:tStep:N;
t(end) = [];
N = length(t);
dt = mean(diff(t)); % Average dt
fs = 1/dt; % Frequency [Hz] or sampling rate
% some additionnal high pass filtering
N = 4;
fc = 0.2; % Hz
[B,A] = butter(N,2*fc/fs,'high');
acc2 = filter(B,A,acc);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
velocity = cumtrapz(dt,acc2);
velocity = detrend(velocity);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
disp = cumtrapz(dt,velocity);
%
figure(1)
plot(t,acc);
figure(2)
plot(t,velocity);
figure(3)
plot(t,disp);
I have a second Measurement where i can cleary see that i have a constant drive since second ~4 with 50mm/s till second 10. The Measurement ist attached, but after ~10 Sekunds the chart is rising repeadly. I calculated the absolute of the x and y vector after filtering with gauß, i want nearly the same graph like this one but without the rapid rising of the chart at second 10 because 1 100% know its wrong. I know that i have a constant speed of ~50 mm/s since second 4 till 30. But the chart doesnt show me right.
I now added to the skript, that all values between -0.03 and 0.03 should be 0 so i could manage to get an slightly better result.
In this diagramm you can see after the adding the small "fix" that i have a constant drive of ~0.05 m/s. but it doesnt help in the end. It should be till second ~31 50 mm/s after it I have a 1 second delay and when i very fast drive.
hello again
I just come back on test1 data and did a quick check ... with still a couple of questions to you :
  • your acceleration data are expressed in g or in m/s² ?
  • in the beginning of your post, you mention that your "rested" signal (sensor at rest) correspond to samples 300 to 500 but in the code below , for the y direction I could not really see that (would correspond to time between 6 and 10 s)
  • so in my case, as "reference" zero acceleration I took the average accel value from t = 0 to 1.7 s and that give me this result. As we can see the velocity remains zero until t = 1.7 s , which is ok - as we consider this as "zero" acceleration - once we have corrected the accel data
  • but after we have a trend and I cannot say if it's real of an numerical artifact (which shoud not be the case if we are sure if the calibration)
  • if we change the static acceleration correction (calibration), this of course will create a different linear trend but that is not "scientific" by any means , as the only truly reliable method is to a have a record (long enough) with garanteed zero mechanical acceleration; after that there can be arguments about how accurate and reliable is the sensor and electronic system used in this experiment
code for my quick test
clc
clearvars
load('test1.mat');
x_ac = Acceleration.X; %Beschl. x-Achse
y_ac = Acceleration.Y; %Beschl. x-Achse
z_ac = Acceleration.Z; %Beschl. x-Achse
% Time
tStep = 0.02; % Length of each time step
N = length(x_ac);
t_ac = (0:N-1)*tStep;
fs = 1/tStep; % Frequency [Hz] or sampling rate
% acc = y_ac*9.81; % Add No gravity factor (assumes data in g)
acc = y_ac; % No gravity factor (assumes data in m/s²)
% remove "rest" acceleration (must start at zero)
ind = find(t_ac>0 & t_ac< 1.7);
acc_ref = mean(acc(ind));
acc = acc - acc_ref; % ??
% acc = acc - mean(acc); % ??
velocity = cumtrapz(tStep,acc);
figure(1)
subplot(211),plot(t_ac,acc);
title('acceleration (unit ?)');
subplot(212),plot(t_ac,velocity);
title('velocity (unit ?)');
xlabel('Time (s)')
plot :

请先登录,再进行评论。

I would suggest your problem may be a simple one.
When you integrate a function, there is a constant of integration, which is unknown. It can be any value. Effectively, it comes down to the idea that you need to know the INITIAL velocity. Or the final velocity. Or the velocity at any given point in time. Think of it as if you are solving a differential equation. You need that initial value before the solution is complete.
Again, the integral of acceleration is velocity + C. So if you know that the final velocity must be zero, then you can correct for that by subtracting the final velocity that was obtained. Subtract that value from the entire curve.
The implicit assumption with cumtrapz is that the INITIAL velocty was zero. If that is not true, this is the source of your problem.

类别

帮助中心File Exchange 中查找有关 Mathematics 的更多信息

产品

版本

R2019b

标签

Community Treasure Hunt

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

Start Hunting!

Translated by