Magnetometer tilt compensation without accelerometer

I need to use the magnetometer in a device to accurately determine how many degrees the user twists (rotates in x-y plane) it. I don't care about absolute orientation, just the angular delta. Of course the user can hold the device in any orientation. Essentially, this is like a doorknob that isn't attached to anything.

I tried atan2(x, y) but the resulting anglular delta varies depending on the device's orientation. Some research suggests that tilt compensation is needed, but all solutions are for correcting compasses based on knowing the gravity vector.

In this application, accelerometer, gravity, gyro data are not available -- just the magnetometer. This means that the pitch & roll used in typical solutions aren't available. However, I'm not looking for true headings, just rotation in x-y.

Can anybody explain the math & solution for this? I realize this is strictly speaking a programming question, but this is stuff that only an EE would know.

The application is typically used such that the x-y plane is rarely, if ever, perpendicular to the earth's field, so sqrt(x2+y2) > 0.

The magnetometer is an AK8975 3-axis MEMS soldered to the device with no gimbal-type leveling arrangement. Hard-iron calibration is already taken care of. Sampling rate is 40 Hz whereas the fastest you can twist your hand back & forth is less than 5 Hz. The magnetometer's orientation is totally unknown since the magnetometer's axis is fixed to the device and the device can be in any orientation. Nothing is known or can be assumed about gravity nor the device's orientation or location.

[Later] Did an experiment where I took the device and lay it flat on a swivel chair facing North, then rotated it thru 360. The chair seat is level with the ground. Here's the chart: x=blue, y=dark green, z=red, sqrt(x^2+y^2+z^2)=aqua The x-axis is elapsed seconds, y-axis is uT.

Two things jump out at me:

1. Shouldn't the aqua-colored magnitude be flat? The trough is 30 vs. 45 for the flat part. Otherwise, is this a tipoff that this sensor is grossly mis-calibrated? Or is this an indication that tilt compensation is needed?
2. The x & Y peaks are slightly offset and have slightly different ranges, but not by that much. This, I'm willing to shrug off as error, as opposed to the 30:45 difference for the magnitude which seems sufficient to mess up any calculations.
3. The Z is flat as expected
• What if the vector of magnetic intensity points in the direction of z-axis? – venny Sep 16 '14 at 19:40
• What kind of magnetometer? A two axis flux gate? Is the orientation of the z axis fixed? (assumed parallel to the earth's g-field) Do you know the magnetic field vector where you are? Which direction does it point? Will it be used inside a building? – George Herold Sep 16 '14 at 19:44
• +1 to GeorgeHerold. Also, can the magnetometer be sampled at a much higher rate than it will move, i.e. many times higher than the Nyquist rate. How slowly might you need to detect? (I'm wondering about drift). What constraint, if any, might you be able to relax or remove if the answer is "it can't be done". – gbulmer Sep 16 '14 at 20:25
• It is impossible because rotation along magnetic field vector is undetectable. – venny Sep 16 '14 at 20:53
• Quantisation error is the last thing I would be afraid of. But zero-gauss offset is $\pm 1000\,\mathrm{LSB}$ and sensitivity may vary by $\pm 5\,\mathrm{\%}$. And if you want to ensure that projection into XY is non-zero, you are restricting yourself to horizontal position. – venny Sep 16 '14 at 21:47

If the orientation of the device is fixed during the measurement then the computation is relatively straightforward as Dave Tweed has pointed out.

The main issue is calibration. A quick way would to establish at an approximate min, max for each x,y,z direction. Hopefully the zeros will lie roughly on the averages of the min and max. That is, this will give an $x_\min$, $x_\max$ and similarly for the other axes. From this estimate $x_\text{zero} = {1 \over 2} (x_\min+x_\max)$, and $r_x = {1 \over 2} (x_\max-x_\min)$, etc.

Then for a reading, compute $x_\text{est} = { x_\text{meas}-x_\text{zero} \over r_x}$, etc. Hopefully this will be in the range $\pm 1$.

I would check that $|z_\text{est}|$ is (1) no larger than some number ($<1$ here) so that you get a meaningful direction, and (2) that it remains fairly constant for the duration of a measurement.

To compute an angle, use $\text{atan2}(y_\text{est} , x_\text{est})$ (check your particular API). To compute the difference with another (hopefully nearby) angle, just be careful with the $\pm\pi$ boundary.

This is obviously a quick overview with many optimisations and sanity checks omitted.

• As an aside: The following discusses the mathematics behind one class of approaches to calibration ulb.ac.be/assoc/bms/Bulletin/sup962/gander.pdf. – copper.hat Sep 17 '14 at 14:45
• As referenced by venny, Freescale's application note seems more robust for limited time frames where you can't guarantee the user has swept through the full range of angles in all dimensions. – MilesK Sep 17 '14 at 18:42
• @MilesK: From what I see in the app. note, they use least squares to fit to their model (a more parameterized ellipse) that models declination as well. You can make a quick calibration using measurements and a SVD that estimates the parameters of the ellipse. I wrote the above so you could (hopefully) make sense out of your readings. – copper.hat Sep 17 '14 at 18:57

If the orientation is absolutely unknown, you are trying to solve a 3D problem with a 2D computation. This cannot work. You need to solve the issue in 3D.

Concerning the calibration, @copper.hat solution is correct, but you have to do the calibration in 3D, which is significantly more difficult to execute < you must either (1: more movements) follow random 3D movements that cover the sphere of possible positions, or (2: complex mechanical setup) start by defining one preferred direction that is aligned with the 3D earth field, then calibrate 2D around this axis, then calibrate the third dimension around one perpendicular axis, or (3: complex computation) select an random orthogonal base and reverse the trigonometric computation to get your calibration variables - this is the most noise sensitive approach>. (1) is what is done by simple consumer objects (like the PS3 controller), (2) is what is done in industrial environment, (3) is what is done when the other two are not available.