Converting acceleration data to displacement data using FFT results in Zero-drift

I have a (n x 1) vector containing acceleration data from a sensor and I am trying to obtain the displacement by:

1. performing an FFT
2. Convert acceleration to displacement in the frequency domain by dividing by -omega^2
3. 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)).';

[Mathieu NOE]

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

[Johann Vormadal]

Hi,

In my code I remove my DC offset before applying the fft.

The requirement is to have zero drift. the fft/ifft seems to be the best way to acheive this (in theory).

Answer 1

[Mathieu NOE]

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)');

[Johann Vormadal]

Thank you for the response

I was not aware of those functions, but it seems there is still a bit of zero drift in your displacement result, which I cannot have.

Since my displacement result doesn't even reach zero it seems to me like I am applying the fft/ifft wrong. It shouldn't be possible for the ifft to produce a constant zero offset when the 0hz bucket is emtpy.

[Mathieu NOE]

a slight improvement , as you can "play" with the order of detrending ...(constant, or nth order polynomial)

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 = detrend(acc,0);    % baseline correction

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
velocity = cumtrapz(dt,acc);
velocity = detrend(velocity,1);    % baseline correction

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
disp     = cumtrapz(dt,velocity);
disp = detrend(disp,3);    % 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)');

[Johann Vormadal]

I should have mentioned in my original post that I will need to run this script for many samples, so I need a solution that doesn't require any manual fiddling.

I also notice that the phasing in the 'disp' graph is slightly misaligned towards the end with the phasing on the 'accel' graph.

I wonder if cropping the raw data to a whole number of cycles would improve the quality of the results.

[Mathieu NOE]

have to try a few other options , like :

Finding the velocity from displacement - MATLAB Answers - MATLAB Central (mathworks.com)

Converting acceleration to displacements records - File Exchange - MATLAB Central (mathworks.com)

there are probably more ...

[Johann Vormadal]

Thank you for those links!

Those also use the fft method, but their code is slightly different and much more elegant. I will try and emulate them, and see if it solves my problem.