Converting acceleration data to displacement data using FFT results in Zero-drift
35 次查看(过去 30 天)
显示 更早的评论
I have a (n x 1) vector containing acceleration data from a sensor and I am trying to obtain the displacement by:
- performing an FFT
- Convert acceleration to displacement in the frequency domain by dividing by -omega^2
- performing an IFFT to get back to the time domain.
The idea is to make a conversion which doesn't result in any Zero-drift, but all my attempts so far result in the oscillation being superimposed on a very low frequency wave (See figure attached), and I can't tell where this low frequency oscillation is coming from. The result should be a near sinusoidal wave centred around zero.
Below is a (slightly modified) snippet from my code. The sensor data (Ch5) is attached as a .mat file.
Tinterval = 0.000625;
Length = length(Ch5);
SamplingFrequency = 1/Tinterval;
Frequency = SamplingFrequency*(0:Length-1)/Length;
Omega = 2*pi*Frequency;
ConversionFactor = -(Omega.^2).';
ConversionFactor(1,:) = 1; %To avoid dividing by zero
% Base Data
RawData = Ch5;
Data = RawData - mean(RawData);
FourierTransform = fft(Data);
Y_accel_base = ifft(FourierTransform).';
dY_base = abs(ifft(FourierTransform./ConversionFactor)).';
2 个评论
Mathieu NOE
2024-5-27
hello
doing it by fft/ifft is a requirement or just an idea ?
also notice that the raw data is not zero mean, so you're going to have drift on first integration then a parabolic trend on the 2nd integration
回答(1 个)
Mathieu NOE
2024-5-27
FYI, a very simple time domain approach - assuming we don't care about the DC values , only the dynamic portion of the data is of interest
results may be even refined by adding some high pass filtering if you really need it
load('Ch5 data.mat')
dt = 0.000625;
fs = 1/dt; % sampling rate
acc = Ch5;
samples = numel(acc);
t = (0:samples-1)*dt;
% remove DC value
acc = acc - mean(acc);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
velocity = cumtrapz(dt,acc);
velocity = detrend(velocity); % baseline correction
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
disp = cumtrapz(dt,velocity);
disp = detrend(disp); % baseline correction
%
figure(1)
subplot(311),plot(t,acc);
title('accel'); ylabel('m/s²');xlabel('time (s)');
subplot(312),plot(t,velocity);
title('velocity'); ylabel('m/s');xlabel('time (s)');
subplot(313),plot(t,disp);
title('disp'); ylabel('m');xlabel('time (s)');
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!