Analyze Motion Combining Gyrometer, Accelerometer, and Magnetometer Data


Analyze iPhone Motion And Orientation By Combining Gyrometer,
Acceleromerter, and Magnetometer Sensors


ABSTRACT

Gyroscope or Gyrometer sensors provide rotation data that must be used to extract accurate motion and orientation from the data provided by accelerometer sensors. By incorporating gyrometer (rotation) with accelerometer data we can accurately calculate acceleration and orientation, then subtract the acceleration of gravity, and the remaining acceleration can be integrated once to obtain speeds along the x, y ,z axises, and integrate again to obtain displacements. The spreadsheet provided here may be helpful / useful to some, who may not need to figure out how to combine gyrometer with accelerometer data, but would get a relatively simple view of how it is done. Then use gyrometer, accelerometer, and magnetometer sensor data with this spreadsheet to obtain relatively accurate useful results for fun projects.

Download this article (pdf):  Click Here.

INTRODUCTION

The following will explain a spreadsheet that has been created with the incorporation of gyroscope data to extract relatively accurate motion and orientation for a device that contains gyrometers, accelerometers, and magnetometers, like the iPhone. It is well known that accelerometers by themselves cannot be utilized to accurately characterize motion and orientation of a device containing them. Analyzing data from accelerometers only cannot properly distinguish effects from rotation changes. Meaning, that accelerometers by themselves can give “apparent” orientation changes (rotations) in pitch and roll. And rotations (changes) in pitch and roll can appear as “apparent” accelerations. This can be seen in fig 6 of a previous publication. That being said, this publication is a follow up with a more encompassing analyses for “Calculating Everything From iPhone Sensors”. The previous publication mentioned above was done when iPhone's did not have gyroscopes. Today gyroscopes, accelerometers and magnetometers are common in so called smart phones. Hence this new spreadsheet supersedes the one in the previous publication.

Smart phones with the above sensors can be easily utilized in simple physics experiments, like attach an iPhone (using a cover) with velcro on a radio controlled toy, … all up to the imagination.

To utilize the data from the sensors one must have the ability to extract the data from the sensors. There are lots of “apps” in the “App Store” of Apple's iTune that can provide the data from these sensors. The app must provide a downloadable file with the sensor outputs in a format easily imported to a spreadsheet, such as CVS formatted. You will also need some or develop some proficiency with utilizing a spreadsheet. The spreadsheet that is provided here has been automated so that after entering the sensors data, the calculations are automatically performed. The explanations provided here are performed within the context of utilizing the spreadsheet that is provided. The spreadsheet is that of the free OpenOffice suit. It is called OpenOffice Calc.  (ClickHere) to download the spreadsheet: Sensor Analysis Master EF yN 3.ods, it should automatically download, so look in your download folder. The Excel version is no longer provided since it does not work correctly. To use the spreadsheet you will need to download the open source, free OpenOffice Suit. You may Google, OpenOffice. Or use this link:  http://www.openoffice.org.

If you don’t find the charts that are in the spreadsheets useful, change/delete them and make your own. Make a copy of the spreadsheet and keep the original as the master file, make changes as needed on the copy. Furthermore, even though a cell may look empty, it may contain a function, so check before you delete cells (row and columns). Click on a blank cell to see if it is empty or not. In the worksheets usually starting with row 4, we have used a function that starts as follows

=IF(ISBLANK(Data.A4);"";           for Analysis 1, and

=IF(A4>=(COUNT(Data.A$3:Data.A$65556));"";           for Analysis 2.

Notice this is the part of the functions in the cells that set cells to a blank after the last data row, up to row 3000.

Two coordinates systems are utilize for understanding the motion and attitude of a device that contains gyrometers, accelerometers, and magnetometers. We have the “body frame” with it's fixed x, y, z axis defined for an iPhone as follows:



Figure 1

Note the +Z-axis is pointing away from the front and is always perpendicular to the face of the iPhone. The X-axis is along the short length, and the Y-axis along the long length. The +Y-axis is toward the “on” button on the top right of the iPhone. The orientation of the axes follows the Right-Hand Rule. As these axes are fixed to the mobile devices Body, rotating the device will rotate the axes to maintain their configuration as per figure 1.

We define the “Earth Frame”, that of the grounded observer as follows: A unit vector along the x-axis in the horizontal plane and pointing East; then a unit vector along the y-axis also in the horizontal plane and pointing in the direction of magnetic North; then following the Right-Hand Rule the unit vector along the z-axis is perpendicular to the horizontal plane and pointing up. Motion and attitude are calculated from this fixed reference frame or initial x, y, z Earth Frame position. Gyrometers, accelerometers, and magnetometers make measurements with respect to the “Body Frame” on which they are attached, so they move and rotate with the device body. We then convert the measurements made in the “Body Frame” to values that would be observed in the “Earth Frame”.
SOME MATHEMATICS

With understanding of kinematics, the transformation of a vector measured in the Body Frame of reference, to it's representation in the Earth Frame of reference is defined by a 3 x 3 rotation matrix will call R. Where each row is a unit vector representation of the Earth Frame x, y, z axes respectively, as measured in the Body Frame of reference. Hence:

        . 1

And for the transformation of a vector from the Earth Frame of reference back to the Body Frame one must use the inverse rotation matrix of R. Which is created by applying the transpose of R.

     , 1a

where vEF(t) is the velocity vector in the Earth Frame.

Applying the convention defined above for the Earth Frame axises, R at an initial time t0 can be written as

     2

where we use the unit vectors components of the Earth Frame x, y, z axises, E for East, N for North, and U for up, each is measured in the Body Frame. The quantities for equation 2 can be calculated from the accelerometers and magnetometers as vector products as follows:

     3


where the vectors G (gravity) and B (magnetic north) are measured in the Body Frame by the accelerometers and magnetometers under initial static conditions. Initial static conditions is defined at t0 as the device (iPhone) is at rest, and the Earth's magnetic field is the only magnetic field present. No metals nor current carrying wires near the device as they will give an incorrect direction for magnetic north. However, the equations above will still work, but the calculated east and north will not be in the directions of true east and true magnetic north. But they will be the initial or starting conditions from which all subsequent calculations are performed.

The accelerometers gives the values for the acceleration aBF at a time t along the x, y, and z axes of the Body Frame. Then to calculate the acceleration aEF as would be measured or observed in the Earth Frame, we must track the changes of the rotation matrix R(t) starting from the initial static condition R(t0). That is we need the values of R(t+dt) based on knowing the values at R(t), where dt is a very small time interval which equates to very small angle rotations. Without going thru all the math, we write for the small angle approximation

       4

Where

    5

is the magnitude of the rotation vector ω at time t+dt as measured by the gyrometers in the Body Frame, multiplied by dt. I and B are:

    6

The expression of equation 4, in brackets, is a matrix that transforms the rotation matrix from R(t) to values at R(t+dt). By substituting equations 2, 5, & 6 into equation 4, we can track the changing rotation matrix using the gyro data, starting from the initial static condition at time t0 to a time t0+dt after a very small time interval dt. Then iterating for each subsequent time interval dt to update the previous rotation matrix values.

The iteration process usingof equation 4 with equation 1 is relatively easy to implement for calculations using a spreadsheet. We then can subtract the acceleration of gravity, and the remaining acceleration can be integrated (or numerically integrated) once to obtain speeds along each of the Earth Frame coordinate axises, and integrated again to obtain displacements also along each the Earth Frame axises. Multiply the acceleration in units of g's by Standard Gravity for Desired Units: 9.80665 m/s2, or 32.174 ft/s2.


     7


The integrations in equations 7 with iteration can be expressed with the following algorithm for use with a spreadsheet.


8

Using numerically integration the results can also be expressed as

    9


for very small dt. Either equations 8 or 9 are easy to implement with a spreadsheet.


THE SPREADSHEET

The first two rows of each worksheet has descriptions that hopefully are sufficiently clear.

A) Calibration

At least in early iPhone's it was common for the accelerometers to require calibration. This could be observed by measuring the value of the acceleration of gravity at various orientations, namely in the + and – axis directions. By pointing the various axises in the down direction, values of 1 or -1 should be obtained. For some iPhone they could be more than 30% off. The procedure for calibrating the accelerometers can be found elsewhere and therefore will not be discussed here further. Note: the worksheets for performing the calibrations are also contained in this new spreadsheet and can be used as described in the calibration paper referenced above. The spreadsheet is setup to automatically correct the entered accelerometer values from the calibration procedure.

B) Data

The spreadsheet contains a worksheet titled “Data”, the data from the sensors is entered in this worksheet in columns A thru J. Note the labels in row A2 thru J2, the correct data must be entered in the appropriate column. In the worksheet (tabs) in the spreadsheet equations have been propagated down for data to row number 3000 where appropriate. If you use the same spreadsheet over and over, be mind full that the if the next time you use less data in columns A thru J, than before, then in the worksheet titled “Data”, you need to delete the content of the excess previous data rows. But delete the content not the formulas.
After the data is entered in the appropriate columns auto calculations are performed in the worksheets: “Data”, “Analysis 1”, and “Analysis 2” down to 3000 rows (or data points).

C) Analysis

The expression
     10

in equation 4 is calculated in the worksheet “Data”. The final result for this expression is in columns “AL” thru “AT”. Remember, equation 10, is for tracking how to change the rotation matrix from R(t) to R(t+dt).

The worksheet “Analysis 1” contains the calculations for correcting the accelerometer values by the calibration procedure. Results are presented in columns E thru G. Furthermore, this worksheet contains the calculations for the representation of the Earth Frame x, y, z axises in terms of unit vectors as are defined in the Body Frame. The corrected accelerometer values and magnetometer values taken as vectors are normalized, and results shown in columns J thru L and N thru P, respectively. When taken under static conditions, the negative of the normalized corrected accelerometer x, y, z values in columns J thru L are defined as the Earth Frame Up z-axis unit vector. Under static conditions, the accelerometers should always read the down Gravity vector x, y, and z values. The Earth Frame x-axis unit vector defined in the Body Frame as the East direction, is presented in columns W thru Y. The Earth Frame y-axis unit vector defined as the North direction, is in columns AF thru AH. These results for “East”, “North”, and “Up” are obtain from applying the vector cross products of equations 3 above.
The data in the Analysis 1 worksheet starts in row 3, but the time here does not start at the 1st time interval. That is because the app being used here for recording the sensors data usually does not read the magnetometer values until after about 1.7 sec. Hence data for times less than 1.7 sec has been ignored and not entered as data in the Data worksheet, cell A3. The “Analysis 1” worksheet also contains the accelerations and speeds along the Body Frame axises. The Body Frame accelerations and speeds have been transformed from the Earth Frame calculations using the inverse (also known as the transpose) rotation matrix. The acceleration results are in columns AJ thru AL; speeds results are in columns AO thru AQ. The displacements along the axises are also presented in columns AO thru AQ, and are calculated by integrating the speeds of columns AO thru AQ. For the accelerations, speeds and displacements charts, the time axis is that defined in column B of worksheet “Analysis 2”, $'Analysis 2'.$B$3:$B$3000.

The worksheet “Analysis 2”, presents the results of transforming the accelerometer values measured in the Body Frame to values that would be observed in the Earth Frame, and the quantities of interest such as pitch, roll, a yaw value, speeds and displacements as would be measured by a stationary observer in the Earth Frame at the initial static condition. Note, the results in the “Analysis 2” worksheet does not contain all the data points (rows) in the “Analysis 1” worksheet. The results presented in the “Analysis 2” worksheets starts from whatever data point (row number) has been designated as the initial starting static condition. This initial starting condition is user set in cell BH1 of the Analysis 1 worksheet, and for this example, the initial starting condition is chosen at 6.23 sec, which correspond to data row 456. Then in worksheet Analysis 2, cells A3 and B3, start with data row 456 of the Analysis 1 worksheet, and time is set start as t0 = 0 sec, respectively. The functions in these cells will not be explained, if you have a reasonable understanding of spreadsheet functions, you will be able to decipher the functions. The functions starting from row number 4 are then propagated down to row 3000.

In worksheet Analysis 2, The updating of the rotation matrix using equation 10 (Analysis 1) starts in cells C3 thru E3, which contains the East unit vector components of the Earth Frame, as calculated in Analysis 1 at the chosen initial starting condition, of t0. Similarly, F3 thru H3, contains the North unit vector, and I3 thru K3, contains the Up unit vector. Then for this columns, starting with row 4, the iterations for updating the rotation matrix, starting from the initial unit vectors defining the Earth Frame axises are performed for tracking the Earth Frame axises after each time interval dt.

The accelerometer values measured in the Body Frame are then transformed to the values that would be observed in the Earth Frame using the rotation matrix described above after each time interval. This are calculated in columns M thru O. And in column P we subtract 1 from the z component of the acceleration. In the Earth Frame, the Gravity vector has only one nonzero component (0, 0, -1), that is down along the z-axis, hence we subtract the z-axis component from -1. This gives motion in the + z-axis as 'Up”.

Column R contains the magnitude of the acceleration calculated as

A' = SQRT [ A'EFx^2 + A'EFy^2 + ( A'EFz  - 1)^2 ]

To proceed, we need to be able to distinguish when the device (iPhone) is not moving. Note: that just because the acceleration can be 0, the speed may not necessarily be 0. The speed maybe a constant value other than 0. Furthermore, an instantaneous speed of 0 does not mean the acceleration is 0. So we need to know when both accelerations and speeds are 0, thus no motion. This should should become clear from performing the integrations, small values, non-zero, are added up in the integration process and distort results. Hence we need to set the magnitude of the acceleration and the speeds along each axis to “0” when the device is not moving. This can be performed by noting the magnitude of the rotation vector shown in column “AY” of the “Data” worksheet. This is plotted, note that for no-motion the magnitude of the rotation vector is less than ~ 0.06. And as with any sensor there is always some noise. Hence we use 0.06 as a coefficient to kind of filter the magnitude of the acceleration and the speeds along the Earth Frame axises. This coefficient can be set in cell AY2 in the Analysis 2 worksheet. You may need to adjust this value. That is, we set the magnitude of the acceleration and the speeds along the Earth Frame axises to “0”, for times where the magnitude of the rotation vectors is less than 0.06. Results for the acceleration magnitude are in column S. For the speeds along the x, y, and z axises results are in columns AE thru AG. The part of the function in columns AE thru AG,

IF(INDIRECT("Data.AY"&A3)<S$2;0;

performs the filtering. The author has found this approach useful, however, you maybe clever, and find a different approach. Note, you may have to adjust this coefficient. If speeds get very low ~3.5 Ft / Sec, you may have to lower this value of this coefficient. Or if you have a noisy initial condition, it may need to be increase. If the value of this coefficient is made to large, it will distort the calculated speeds and hence the displacements. Be mindful, start with a value of “0”, then increase to “0.01”, and study the resulting speeds. Increase till the appropriate results are obtained. The author has found that a value of “0.06” seems to work well for the few experiments that were performed. If you set this filter coefficient to “0”, you might need to use other techniques to remove background form the integration process. Often there will be a linear with time background, due to the initial static condition (no motion) that has a Pitch and or Roll value not equal 0 deg. That is why this author uses the filter coefficient.

D) Results

The results presented in this “Master” Spreadsheet is for an experiment consisting of (while collecting the sensors data) raising the iPhone from the surface of a mattress up to a height, then dropping the iPhone on the mattress. Falling on the mattress it bounced 3 or 4 times and came to rest with the Body Frame z-axis pointing down. Note that the Body Frame +z-axis and the Earth Frame +z-axis are co-linear, and point in the same direction. Furthermore, for this experiment we are only interested in analyzing the up and down motion, the z-axis motion.

The iPhone accelerometers have limitations, they are set to saturate at + and – 2 g's. Meaning they cannot accurately measure acceleration values greater than +2 g's, or less than -2 g's. Dropping the iPhone onto a mattress will exceed the +/- 2 g's limitations on impact (but not before). Though the accelerometers will continue reading accurate accelerations after a fraction of a sec from the impacts, the data lost in that fraction of time is sufficient to make analysis of speeds and displacements useless after the 1st impact. Since acceleration data is lost in the fraction of time during impact, the speeds cannot be tracked during the impact (integrated from the acceleration data); therefore, manually we may “adjust” the speeds and displacements to “0” after the 1st bounce. Since we are only interested in the z-axis motion in this experiment, columns AT and AU are set to z-axis speed and displacement respectively. Then the “adjustment” is performed and shown in columns AT and AU (in the Analysis 2 worksheet) for the z-axis motion. The z-axis speeds and displacements are taken from columns AL and AQ respectively. The results are shown in the graphs (Figures 2 thru 4). The user may need to make adjustment to the appropriate speed results, or just ignore or delete the content in this columns. Some experiments may not need adjustments, hence the user can focus in the appropriate x, y, z axis motion or the magnitudes in columns AE thru AM (Analysis 2 worksheet) for the Earth Frame, and columns AO thru AW (Analysis 1 worksheet) for Body Frame motion.

Results from a second experiment is shown in figure 5 below.



Figure 2. Accelerometer magnitude from x, y , z axis Body Frame readings.


Figure 3. Speed along the z-axis in the Earth Frame. The speed is adjusted to “0” after 1st bounce.



Figure 4. z-axis displacement in the Earth Frame is performed by numerical integration of the z-axis speed in the Earth Frame. The displacement is adjusted “0” after the 1st bounce.


Figure 5. Different experiment: Toss up an iPhone with Body Frame z-axis pointing up; then falls on a mattress. Bounces 3 or 4 times, and mattress is below initial height. Distance (yellow) plotted on right axis. Speeds and displacements are adjusted to “0” after 1st bounce. Note: Careful not exceed accelerometer saturation (+/- 2 g's) during toss up.


CONCLUSION

Combining gyrometer, accelerometer and magnetometer data, makes it possible to conduct simple physic's experiments on motion and orientation using an iPhone to collect the data, or a device with the above sensors, and a spreadsheet for performing the analysis. A spreadsheet is provided here containing all the calculations for performing the analysis. The analysis consists of calculating accelerations, speeds, distances (displacements), pitch and roll, and the orientation of the Body Frame x, y, z axises from magnetic north. The analysis is performed by effecting a transformation (with a rotation matrix) of the accelerometer readings from the Body Frame to an Earth Frame of reference. The inverse rotation matrix is then use to transform results from the Earth Frame back to the Body Frame. The rotation matrix is constructed with the gyrometer data to obtain accurate rotations and accelerations values, that are free from “apparent” artifacts that result from calculating rotations from accelerometer data only. Have fun using the spreadsheet !!

9 comments:

  1. Got absolutely no idea what you are doing but I am some impressed with all them formulas and such. Took me three times to finally pass Algebra I which is probably why I made a living digging ditches.
    Your old high school buddy Bob Gomes

    ReplyDelete
  2. Hello,
    I am trying to do some vibrational analysis on an object. This spreadsheet would be perfect, however once i enable editing in excel, the whole spreadsheet goes into whack. Any help on this? Thanks in advance.

    ReplyDelete
  3. Juan,
    Sorry to read that the Excel version seems to be problematic. I don’t have Excel, therefore, regretfully I can’t help / offer a fix. Consider using the “free” Apache OpenOffice open source suite. It is very good and is continually improved. Don’t know whether you are using Windows or OS X, I built the spreadsheet with the OS X (Mac) version of Apache OpenOffice, likely I used an older version than what is currently available from the Apache Foundation. They have version's for Windows, OS X, and many other operating systems.

    If I think up anything useful, will let you know. !!

    http://www.openoffice.org/download/

    ReplyDelete
  4. Juan, Happy New Year !!
    Have sent you two files, were they able to help you ?

    ReplyDelete
  5. Happy New Year Rouel! Thank you for reaching out. I apologize, I didn't realize you responded on here! Unfortunately I didn't receive those files. If you could please sent them to juancarreon89@yahoo.com that'd be great!

    Thanks!

    ReplyDelete
  6. Hi Rouel,
    Thank you very much for sharing the knowledge and this great post.
    Is it possible to analyze the errors of the IMU for example bias and scaling for accelerometer and gyroscope??
    Please shed some light on this...that would be really helpful.
    Thanks in advance!!

    ReplyDelete
    Replies
    1. Hi Ravindran Ji,
      Regretfully, not able to shed any light on the errors of the IMU. Have not study on it. As you know, you probably have to Google it, “errors of IMU”.

      Good luck finding something useful in the internet.

      Delete
  7. Hi,
    I would join one of the previous comments as well: this text presents such a practical approach which way I can approach this topic with higher confidence.
    I also tried to reach the spreadsheet to take a closer look, but without success. Would you be so kind to send this also to me: grektor1980@gmail.com
    Thank you

    ReplyDelete
  8. My apologies, that you were not able to download the spreadsheet. Will email to you.

    ReplyDelete