How to integrate acceleration into displacement data correctly?

Question

Alireza 2024-4-2

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2102001-how-to-integrate-acceleration-into-displacement-data-correctly

评论： Paul 2024-4-2

Hi all,

I am trying to extract displacement data from acceleration dataset. Even after detrending and high-pass filtering the result does not look good.

Any hints?

Cod eattached.

Fs = 1000; % sampling freq.

Ts = 1 / Fs; % sampling period or time step

td = 1; % signal duration in seconds

L = Fs * td; % signal length = number of total sampes taken

t = Ts * (0:L-1); % time vector

x = 2.* sin(2*pi*25*t); % cosine function with freq. of: 50 Hz

figure(1);

plot(t, x);

Y = fft(x, L);

P2 = abs(Y/L); % re-scaling

P1 = P2(1 : L/2+1);

P1(2 : end-1) = 2*P1(2 : end-1);

f1 = (Fs/L)*(0:L/2);

figure(2);

plot(f1, P1)

%% Integration

x_detr = detrend(x);

velnf = cumtrapz(t, x);

dispnf = cumtrapz(t, velnf);

% De-trending

x_vel = cumtrapz(t, x_detr);

x_vel_detr = detrend(x_vel);

x_disp = cumtrapz(t, x_vel_detr);

% after de-trending now it is the time for filtering high pass filter

%

figure(3);

plot(t, x, 'k');

hold on;

plot(t, x_disp, 'b');

hold on;

plot(t, dispnf, 'r');

% HP (High Pass) filtering

N = 2; % filter order

cutoff_freq = 1; % cutoff frequency in Hz

[butB, butA] = butter(N, cutoff_freq / (Fs/2), 'high');

accf2 = filter(butB, butA, x_detr); % butter w/ filter

accf3 = filtfilt(butB, butA, x_detr); % butter w/ filtfilt

velf2 = cumtrapz(t, accf2);

velf2 = detrend(velf2);

disp2 = cumtrapz(velf2);

velf3 = cumtrapz(t, accf3);

velf3 = detrend(velf3);

disp3 = cumtrapz(velf3);

figure(4);

plot(t, x, 'k');

hold on;

plot(t, disp2, 'b');

hold on;

plot(t, disp3);

hold off

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Paul 2024-4-2

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2102001-how-to-integrate-acceleration-into-displacement-data-correctly#answer_1435296

在 MATLAB Online 中打开

The integration seems correct based on the implied initial conditions.

Fs = 1000; % sampling freq.
Ts = 1 / Fs; % sampling period or time step 
td = 1; % signal duration in seconds
L = Fs * td; % signal length = number of total sampes taken 
t = Ts * (0:L-1); % time vector 

The acceleration has frequency of 25 Hz, not 50 Hz.

x = 2.* sin(2*pi*25*t); % cosine function with freq. of: 50 Hz

velnf = cumtrapz(t, x);

figure

plot(t,velnf)

dispnf = cumtrapz(t, velnf);

figure

plot(t,dispnf)

We see that cumtrapz starts the integration at 0 for both velnf and dispnf

Let's get the exact solutions for velocity and displacement, don't forget the constants of integration

vel = int(2*sin(2*sym(pi)*25*sym('t')),sym('t')) + sym('C1')

vel =

disp = simplify(int(vel,sym('t'))) + sym('C2')

disp =

In accordance with the cumtrapz, we have to select C1 and C2 such that vel(0) = disp(0) = 0;

C1 = solve(subs(vel,sym('t'),0))

C1 =

C2 = solve(subs(disp,[sym('t') sym('C1')],[0 C1]))

C2 =

0

Now plot

figure

fplot(subs(vel,sym('C1'),C1),[0 1])

figure

fplot(subs(disp,[sym('C1') sym('C2')],[C1 C2]),[0 1])

which are essentially the same plots as using cumtrapz.

4 个评论
显示 2更早的评论隐藏 2更早的评论

Alireza 2024-4-2

ACC.mat

problem is I need to double integrate accelearion data to get displacement data. but looks like there is a drift after integration, besides the fft of displacement is basically empty. not sure if the problem is for the cumtrapz or what ... i have de-trended my dataset and applied high pass filter. imagine instead of x as a sinusoid function i have such data attached:

Paul 2024-4-2

在 MATLAB Online 中打开

As shown, the "drift" in the displacement is because the velocity has a constant offset, called C1 above, and that constant offset in velocity turns into a C1*t in the displacement, which is the drift.

If you think the displacement shouldn't have that drift you can always try to detrend.

Fs = 1000; % sampling freq.

Ts = 1 / Fs; % sampling period or time step

td = 1; % signal duration in seconds

L = Fs * td; % signal length = number of total sampes taken

t = Ts * (0:L-1); % time vector

x = 2.* sin(2*pi*25*t); % cosine function with freq. of: 50 Hz

velnf = cumtrapz(t, x);

dispnf = cumtrapz(t, velnf);

dispnf = detrend(dispnf,1);

figure

plot(t,dispnf)

There might be other options, but, again, I have no idea what the expected result for displacement actually is.

I looked at ACC.mat, but had not idea what all the variables meant or how they are to be processed.

请先登录，再进行评论。

How to integrate acceleration into displacement data correctly?

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

回答（1 个）

4 个评论
显示 2更早的评论隐藏 2更早的评论

另请参阅

类别

标签

Community Treasure Hunt

How to integrate acceleration into displacement data correctly?

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

回答（1 个）

4 个评论 显示 2更早的评论隐藏 2更早的评论

另请参阅

类别

标签

Community Treasure Hunt

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

4 个评论
显示 2更早的评论隐藏 2更早的评论