Measuring Frequency of BNO055 accelerometer using FFT

Hello,

I’ve been trying to measure the frequency of vibration using an accelerometer (BNO055).
First, I wrote the code using the Fast Fourier Transform (fft) library with a function generator as the input.
After making sure that the program gives out the same frequency as the function generator, I tried implementing the accelerometer into the code. However, the results are not correct. I’m outputting the frequencies with their corresponding amplitudes. And at the end it outputs the dominant frequency. Please find code and results below:

#define SAMPLES 128
#define SAMPLING_FREQUENCY 20 //Hz, must be less than 10000 

#include "arduinoFFT.h"
#include <Wire.h>
#include <Adafruit_Sensor.h>
#include <Adafruit_BNO055.h>
#include <utility/imumaths.h>

Adafruit_BNO055 bno = Adafruit_BNO055();
arduinoFFT FFT = arduinoFFT();

unsigned int sampling_period;
unsigned long microseconds;

double vReal[SAMPLES];
double vImag[SAMPLES];

void setup() {
  Serial.begin(115200);
  sampling_period = round(1000000 * (1.0 / SAMPLING_FREQUENCY));
  if (!bno.begin()) {
    Serial.print("Wiring error");
    while (1);
  }
}

void loop() {
  unsigned long t1 = micros();// initiate time to start count

  imu::Vector<3> euler = bno.getVector(Adafruit_BNO055::VECTOR_LINEARACCEL);// calling out acceleration from accelerometer

  //  Serial.print("a = ");
  //  Serial.println(euler.y()); // acceleartion in y direction

  /*SAMPLING*/
  for (int i = 0; i < SAMPLES; i++) {
    microseconds = micros();    //Overflows after around 70 minutes!

    vReal[i] = euler.y();
    vImag[i] = 0;

    while ((micros() - microseconds) < sampling_period); // control constant time-stamp
  }

  /*FFT*/
  FFT.Windowing(vReal, SAMPLES, FFT_WIN_TYP_HAMMING, FFT_FORWARD);
  FFT.Compute(vReal, vImag, SAMPLES, FFT_FORWARD);
  FFT.ComplexToMagnitude(vReal, vImag, SAMPLES);
  double peak = FFT.MajorPeak(vReal, SAMPLES, SAMPLING_FREQUENCY);

  Serial.println("    freq    ampl");
  for (int i = 0; i < (SAMPLES / 2); i++) {
    /*View all these three lines in serial terminal to see which frequencies has which amplitudes*/
    Serial.print(i);
    Serial.print("  ");
    Serial.print((i * 1.0 * SAMPLING_FREQUENCY) / SAMPLES, 1);
    Serial.print(" Hz");
    Serial.print("   ");
    Serial.println(vReal[i], 1);    //View only this line in serial plotter to visualize the bins

  }
  Serial.print("\n");

  /*PRINT RESULTS*/
  Serial.print("f = ");
  Serial.print(peak);     //Print out what frequency is most dominant.
  Serial.println(" Hz");
  Serial.print("\n");

  while ((micros() - t1) < 100000);
}

Results when accelerometer is not moving:

    freq    ampl
0  0.0 Hz   6.2
1  0.2 Hz   2.7
2  0.3 Hz   0.0
3  0.5 Hz   0.0
4  0.6 Hz   0.0
5  0.8 Hz   0.0
6  0.9 Hz   0.0
7  1.1 Hz   0.0
8  1.2 Hz   0.0
9  1.4 Hz   0.0
10  1.6 Hz   0.0
11  1.7 Hz   0.0
12  1.9 Hz   0.0
13  2.0 Hz   0.0
14  2.2 Hz   0.0
15  2.3 Hz   0.0
16  2.5 Hz   0.0
17  2.7 Hz   0.0
18  2.8 Hz   0.0
19  3.0 Hz   0.0
20  3.1 Hz   0.0
21  3.3 Hz   0.0
22  3.4 Hz   0.0
23  3.6 Hz   0.0
24  3.7 Hz   0.0
25  3.9 Hz   0.0
26  4.1 Hz   0.0
27  4.2 Hz   0.0
28  4.4 Hz   0.0
29  4.5 Hz   0.0
30  4.7 Hz   0.0
31  4.8 Hz   0.0
32  5.0 Hz   0.0
33  5.2 Hz   0.0
34  5.3 Hz   0.0
35  5.5 Hz   0.0
36  5.6 Hz   0.0
37  5.8 Hz   0.0
38  5.9 Hz   0.0
39  6.1 Hz   0.0
40  6.3 Hz   0.0
41  6.4 Hz   0.0
42  6.6 Hz   0.0
43  6.7 Hz   0.0
44  6.9 Hz   0.0
45  7.0 Hz   0.0
46  7.2 Hz   0.0
47  7.3 Hz   0.0
48  7.5 Hz   0.0
49  7.7 Hz   0.0
50  7.8 Hz   0.0
51  8.0 Hz   0.0
52  8.1 Hz   0.0
53  8.3 Hz   0.0
54  8.4 Hz   0.0
55  8.6 Hz   0.0
56  8.8 Hz   0.0
57  8.9 Hz   0.0
58  9.1 Hz   0.0
59  9.2 Hz   0.0
60  9.4 Hz   0.0
61  9.5 Hz   0.0
62  9.7 Hz   0.0
63  9.8 Hz   0.0

f = 9.89 Hz

Results when accelerometer is moving:

    freq    ampl
0  0.0 Hz   1586.0
1  0.2 Hz   682.7
2  0.3 Hz   3.6
3  0.5 Hz   1.4
4  0.6 Hz   0.7
5  0.8 Hz   0.4
6  0.9 Hz   0.3
7  1.1 Hz   0.2
8  1.2 Hz   0.2
9  1.4 Hz   0.1
10  1.6 Hz   0.1
11  1.7 Hz   0.1
12  1.9 Hz   0.1
13  2.0 Hz   0.1
14  2.2 Hz   0.1
15  2.3 Hz   0.0
16  2.5 Hz   0.0
17  2.7 Hz   0.0
18  2.8 Hz   0.0
19  3.0 Hz   0.0
20  3.1 Hz   0.0
21  3.3 Hz   0.0
22  3.4 Hz   0.0
23  3.6 Hz   0.0
24  3.7 Hz   0.0
25  3.9 Hz   0.0
26  4.1 Hz   0.0
27  4.2 Hz   0.0
28  4.4 Hz   0.0
29  4.5 Hz   0.0
30  4.7 Hz   0.0
31  4.8 Hz   0.0
32  5.0 Hz   0.0
33  5.2 Hz   0.0
34  5.3 Hz   0.0
35  5.5 Hz   0.0
36  5.6 Hz   0.0
37  5.8 Hz   0.0
38  5.9 Hz   0.0
39  6.1 Hz   0.0
40  6.3 Hz   0.0
41  6.4 Hz   0.0
42  6.6 Hz   0.0
43  6.7 Hz   0.0
44  6.9 Hz   0.0
45  7.0 Hz   0.0
46  7.2 Hz   0.0
47  7.3 Hz   0.0
48  7.5 Hz   0.0
49  7.7 Hz   0.0
50  7.8 Hz   0.0
51  8.0 Hz   0.0
52  8.1 Hz   0.0
53  8.3 Hz   0.0
54  8.4 Hz   0.0
55  8.6 Hz   0.0
56  8.8 Hz   0.0
57  8.9 Hz   0.0
58  9.1 Hz   0.0
59  9.2 Hz   0.0
60  9.4 Hz   0.0
61  9.5 Hz   0.0
62  9.7 Hz   0.0
63  9.8 Hz   0.0

f = 9.89 Hz

You can see from the results that there are mainly two issues:

  1. It keeps giving 9.89 Hz as the dominant frequency even though it has no amplitude.

  2. When the accelerometer moves, you should expect values that are greater than 1 Hz to be present. However, it seems that the motion only increases the amplitude of 0 Hz and 0.2 Hz. This occurs even when the accelerometer is moving quite fast that the frequency should definitely be above 3 Hz.

I don’t know what to change to make the readings correct. Keep in mind that this code worked with a function generator connected to analogue pin 0.

I would appreciate any help!

Thank you!

I'd dump the measured values vReal[], for both signal sources. Eventually the Serial Plotter can be used to verify the curve form, but the amplitudes must be verified as well.

A 9.89Hz frequency at 20Hz sampling rate is bogus. Do you apply an anti-aliasing filter?

Thank you for replying!

First, what do you mean by dump vReal[], for both signals?

When I Serial plot the acceleration (euler.y), the plot and readings are correct. The amplitude of acceleration plot is usually between 20 and -20 (m/s^2). Keep in mind that the euler.y uses SI units. However the amplitudes fft results are different. I'm not sure if they should be the same.

Lastly, to answer your question, no i'm not using an anti-aliasing filter. I thought that since the sampling frequency is definitely higher than twice the frequency I'm measuring, I wouldn't need one. I even tried increasing the Sampling frequency to a bunch of different values. The peak kept outputting the highest frequency. For example, if the Sampling frequency was 300 Hz, the peak reading would be 149 Hz. If it was 800 Hz, the peak reading would be 399 Hz. And so on and so forth. I've been stuck on this for quite some time. I'm not that experienced in both Arduino nor Fourier.

Looking forward to your reply Thank you in advance

I am not familiar with that FFT library, nor with your experiment. How comes that you expect a sinuosidal shaped acceleration from the IMU?

After some research I suspect that the strange MajorPeak frequency is generated by the Hamming filter, possibly together with a signal offset. See a plot or dump of the values after the window function has been applied.

What if you call DCRemoval before Windowing?

I also could not find any documentation of the library, not even comments on the methods - links anybody?

Assuming this is the library you are using, (https://github.com/kosme/arduinoFFT/blob/master/src/arduinoFFT.cpp) the interpolation formula in "MajorPeak" doesn't appear to be protected against the highest peak being at the endpoints. That is, the formula indexes the array point before and after the peak cell in its calculation which will give indeterminate results.

The reason the first cell has the highest value is presumably because the signal input into the FFT has a non-zero mean (i.e. a DC bias) that dominates everything else, but one would have to see a sample of that data to confirm.

I found a document that claims that the first array entry represents the average power of the signal. Dunno whether this is equivalent to the signal bias. Also the second half of the array does not contain amplitudes, in so far peak values at both ends of the frequency range have questionable values, in both frequency and amplitude.

I thought that since the sampling frequency is definitely higher than twice the frequency I'm measuring, I wouldn't need one.

No that is very wrong.

It is what frequencies are in the signal you are sampling. If those frequencies are higher than half the sample rate then you will get aliasing. It is not what you said, those two statements are in no way equivalent.

So it not the frequency YOU ARE MEASURING but the frequencies CONTAINED in the signal you are sampling.

Thank you all for replying!!

The experiment I’m working on is basically using the accelerometer to measure the frequency of vibration of metal or wooden members when oscillated. I’m trying to measure this frequency to eventually be able to find the stiffness of the material. To get a hint of what I’m trying to do, please visit this website, and scroll down to the accelerometer section.

The libraries and codes I’m using were mainly from the links shown at the bottom of the post.

So far, I have added DCRemoval before Windowing. In addition, I added a low pass filter. It seems that the values I’m getting now are more promising, but still far from the expected values. The peak frequencies are no longer constant. However, they don’t output values that make sense. Please refer to code below:

#define SAMPLES 128 // number of samples
#define SAMPLING_FREQUENCY 100 //Hz, must be less than 10000 

#include "arduinoFFT.h"
#include <Wire.h>
#include <Adafruit_Sensor.h>
#include <Adafruit_BNO055.h>
#include <utility/imumaths.h>

Adafruit_BNO055 bno = Adafruit_BNO055();
arduinoFFT FFT = arduinoFFT();

double sensorValue = 0;    //initialization of sensor variable, equivalent to EMA Y
double EMA_a = 0.1;      //initialization of EMA alpha
double EMA_S = 0;          //initialization of EMA S

unsigned int sampling_period;
unsigned long microseconds;

double vReal[SAMPLES];
double vImag[SAMPLES];

void setup() {
  Serial.begin(115200);
  sampling_period = round(1000000 * (1.0 / SAMPLING_FREQUENCY));
  if (!bno.begin()) {
    Serial.print("Wiring error");
    while (1);
  }
}


void loop() {
 
  /*SAMPLING*/
  for (int i = 0; i < SAMPLES; i++) {
    microseconds = micros();    //Overflows after around 70 minutes!

    imu::Vector<3> euler = bno.getVector(Adafruit_BNO055::VECTOR_LINEARACCEL);// calling out acceleration (m/s/s) from accelerometer

    double ay = euler.y(); // euler.y() is the acceleration in y dirction. 

    //    Serial.print("ay = ");
    //    Serial.println(ay);

    sensorValue = ay;
    EMA_S = (EMA_a * sensorValue) + ((1 - EMA_a) * EMA_S); //run the EMA (LOWPASS FILTER)

    //    Serial.print("EMA = ");
    //    Serial.println(EMA_S);
    //    Serial.print("\n");

    vReal[i] = EMA_S;
    vImag[i] = 0;

    while ((micros() - microseconds) < (sampling_period)); // control constant time-stamp
  }

  /*FFT*/
  FFT.DCRemoval(vReal, SAMPLES);
  FFT.Windowing(vReal, SAMPLES, FFT_WIN_TYP_HAMMING, FFT_FORWARD);
  FFT.Compute(vReal, vImag, SAMPLES, FFT_FORWARD);
  FFT.ComplexToMagnitude(vReal, vImag, SAMPLES);
  double peak = FFT.MajorPeak(vReal, SAMPLES, SAMPLING_FREQUENCY);

  for (int i = 0; i < (SAMPLES / 2); i++) {
    /*View all these three lines in serial terminal to see which frequencies has which amplitudes*/
    Serial.print((i * 1.0 * SAMPLING_FREQUENCY) / SAMPLES, 1);
    Serial.print(" Hz");
    Serial.print("   ");
    Serial.println(vReal[i]);    //View only this line in serial plotter to visualize the bins
  }
  Serial.print("\n");

  /*PRINT RESULTS*/

  Serial.print(peak);     //Print out what frequency is most dominant.
  Serial.print(" Hz");
  Serial.println("\n");

}

I’ve attached the results. You can tell when the accelerometer is moving by how the values change. I’ve deleted some of the Serial.prints I had before since they take a lot of memory. In addition, I’ve attached how the filtered acceleration (EMA_S) plot vs time looks like when I move it up and and down.

I would say that there are still two main problems:

  1. When the accelerometer is not moving, it keeps giving out random values. (refer to attached results).
  2. When the accelerometer moves (by hand up and down), the amplitudes of ALL frequencies change, and the peak value is still outputting inconsistent results.

I truly appreciate your replies.
Thank you in advance!

LINKS:

link of FFT library


link of LP filter

link to accelerometer

link to FFT code

Results.txt (31.2 KB)

Your results don’t look like useful spectres. Excitation by hand is forcing random results.

I’d try either a single stroke and record the decay curve, or provoke a resonant frequency by a generator.

So far, I have added DCRemoval before Windowing.

Waste of time and may introduce noise. The FFT calculates the DC offset and presents that as the zero frequency component of the transform.

In addition, I added a low pass filter.

A simple first order filter won’t be effective, but tell us the knee frequency.

You do NOT want the so-called “linear acceleration”, as the process of subtracting the gravity vector introduces substantial errors. Use the raw accelerometer data instead, and for tests, only for an axis parallel to the applied acceleration.

 imu::Vector<3> euler = bno.getVector(Adafruit_BNO055::VECTOR_LINEARACCEL);

In fact you would be better off using a standard, much simpler three axis accelerometer like the ADXL345, rather than the much more complex BNO055.

Frankly, it appears that you do not have a clear idea what you are doing, and should not expect much from the effort as it now stands.

You need to do some very carefully controlled experiments, where you subject the accelerometer to known accelerations at known frequencies, and verify that you are getting the expected results, before moving on.

jremington: (about DCRemoval) Waste of time and may introduce noise. The FFT calculates the DC offset and presents that as the zero frequency component of the transform.

A look at the first component values makes me disagree. This component continues being the highest "amplitude", even after DCRemoval. The Adafruit tutorials claim it's the "average power" of the signal - what's that?

What does the windowing function with a DC offset present?

How reliable is that FFT library at all?

jremington:
Waste of time and may introduce noise. The FFT calculates the DC offset and presents that as the zero frequency component of the transform.

A concern is that applying the window function spreads the power at DC across a few adjacent bins so the skirts of a large DC offset may swamp whatever signal of interest appears in these bins.

A simple first order filter won’t be effective, but tell us the knee frequency.

I’ve never used the BNO055, but the datasheet indicates that it has internal programmable lowpass filtering (down to 12 Hz). If the library supports using this feature that’s probably the way to go.

You do NOT want the so-called “linear acceleration”, as the process of subtracting the gravity vector introduces substantial errors. Use the raw accelerometer data instead, and for tests, only for an axis parallel to the applied acceleration.

 imu::Vector<3> euler = bno.getVector(Adafruit_BNO055::VECTOR_LINEARACCEL);

In fact you would be better off using a standard, much simpler three axis accelerometer like the ADXL345, rather than the much more complex BNO055.

Indeed, a real concern with accelerometers is differentiating the signal of interest from acceleration due to gravity. The latter will have a time varying component if the accelerometer orientation changes during data collection which is inevitable if one is shaking the sensor by hand. This is at least part of the problem with the results presented in the first post.

Frankly, it appears that you do not have a clear idea what you are doing, and should not expect much from the effort as it now stands.

You need to do some very carefully controlled experiments, where you subject the accelerometer to known accelerations at known frequencies, and verify that you are getting the expected results, before moving on.

In addition the way I would approach this problem in the early stage is to use the Arduino simply to capture the sensor data and send it to a PC. I’d use a suitable interactive signal processing/data visualization tool on the PC (e.g. Matlab (Octave) or Python/Numpy/Matplotlib) to work out the algorithm prior to implementing it in the embedded Arduino environment.

Indeed, a real concern with accelerometers is differentiating the signal of interest from acceleration due to gravity.

Agreed, and it is impossible, for an application like this one, with low quality, inaccurate orientation sensors like the BNO055.

How reliable is that FFT library at all?

Standard code that has been around for at least 50 years. Either it works or it doesn't, and that is trivial to test.

OP: see this explanation for the background that explains the problem: [http://www.chrobotics.com/library/accel-position-velocity[/ur ](http://www.chrobotics.com/library/accel-position-velocity)

DrDiettrich: The Adafruit tutorials claim it's the "average power" of the signal - what's that?

The Adafruit tutorial is in error. The first (zeroth?) bin of the FFT is proportional to the average of the signal. In most cases, including this one, the signal represents an amplitude so it is the average amplitude not the average power.

I've submitted a correction request to Adafruit. I consider their tutorials to be among the best available and a great service to the community, so will be interested to see if they update the article.

One more thought: if the true major frequency is the lowest one in the spectrum, it indicates an effectivly still lower resonant frequency. The setup should be changed until a true major peak stands out, faster than the lowest non-zero frequency.