Wednesday, September 7, 2011

Improved magnetometer calibration (Part 2)

In Part 1 of this post, I have presented the MagCal software that is used to calibrate a 3-axis magnetometer (or accelerometer).

One problem with MagCal is that one of the results (the matrix A) cannot be used directly. You have to calculate the inverse matrix A-1.  The other problem is that the matrix A elements have a limited number of significant digits, which affects the precision of the inverted result. And finally, the raw measurements numbers cannot be greater than 99.99, so that you may have to scale down all your data to use it.

I present here my alternate implementation called Magneto, which calculates directly the required matrix A-1.  It also presents the inverted A matrix to permit a comparison with MagCal results.

Magneto v1.2 can be downloaded here:
https://sites.google.com/view/sailboatinstruments1





What are the differences between MagCal and Magneto?
 
Both find the equation of an ellipsoid that best fits the raw data points, using however different techniques. MagCal uses an “adaptive least square estimator”, and Magneto uses the “Li's ellipsoid specific fitting algorithm”. Which one is better? I don’t know, but the authors that have developed MagCal suggest in a recent paper that they may now use Li’s method.

Here is the output of both softwares, using the same raw data file provided with MagCal.



As you can see, the results are nearly identical. I plan to add new features and to fully document Magneto on its hosting site when time permits.

20 comments:

  1. Hi there, awesome site. I thought the topics you posted on were very interesting. I tried to add your RSS to my feed reader and it a few. take a look at it, hopefully I can add you and follow.

    ReplyDelete
  2. Hello,
    Is it possible to modify the source code ? Because I would like to implement the algorithm in an embedded system.
    I can't save a lot of data in a file like in a pc, I have to execute the algorithm with only few data.

    ReplyDelete
  3. It is certainly possible, if you can find or develop an alternate math library (instead of Fortran) for matrix calculations. Matrix multiplication and inversion code are readily available, eigenvalues calculations may be harder to come by. I would strongly consider a 32-bit microprocessor for such a task.

    ReplyDelete
  4. For those who have inquired on the possibility of porting the Magneto source code to an embedded system, here is a site with extensive matrix source code in C language:

    http://www.mymathlib.com/matrices/

    Be aware that there is an important difference in the way the matrix elements are stored in memory between Fortran and C. In Fortran, the elements are stored column by column, while in C and C++, they are stored line by line.

    And the Moore-Penrose pseudo-inverse may also be replaced by the normal matrix inverse, as they are normally the same for a square matrix, unless your data set is bad or incomplete.

    ReplyDelete
  5. I have taken a step further in porting the Magneto source code to 100% C-language. It is available at:

    https://sites.google.com/site/sailboatinstruments1/c-language-implementation

    ReplyDelete
  6. Very nice work. Helps me a lot to write a calibration routine for my lowcost MEMS-Magnetometer. Thank you for sharing.

    ReplyDelete
  7. Hi Mr. Merlin

    I am a bit confused about the input data. When you say 'raw data', does it mean the data come from ADC (whose unit is 'counts') ?
    I've looked into your sample data (sample.txt), and i found the value of all data are less than 1 and greater than -1 .
    How did you get this 'raw data' ?
    do we have to calibrate the ADC output first using parameter provided in the datasheet?

    Thank you for your support.

    ReplyDelete
    Replies
    1. Magneto expects to receive raw data in +- format (a value of zero indicating a null field in the current axis), but not necessarily normalized to +-1.0.

      If your sensors have SPI or I2C outputs, they will usually directly produce the required format. For example, the MicroMag3 magnetometer directly produces counts from -3411 to +3411, and the the SCA3000 accelerometer directly produces counts from -1333 to 1333, and Magneto can process direcly these values, without the need to normalize them to +- 1.0. I understand that a normalization may be desirable to avoid machine precision problems, but this has not been the case with these sensors.

      If your sensors produce voltage levels that you have to convert to counts with an ADC, you have indeed to substract a zero field value from the ADC output before using Magneto. You would then normally choose the maximum positive value as input to the 'Norm of Magnetic or Gravitational field'.

      But this norm value is not critical if all you want to calculate later on is a heading (if it is a magnetometer) or a tilt angle (if it is an accelerometer). You can input any reasonable value for the norm, the correction matrix will be different by just a scaling factor, but the calculated heading (or tilt angle) will be the same, as it depends only on the relative value of the field components. The bias values will be unchanged, as they do not depend on the norm.

      Delete
  8. Could I ask you to explain how you find the value "norm of magnetic or gravitational field" please? This has a big effect on the results and I don't know how to find this value for the program. Thank you.

    ReplyDelete
    Replies
    1. I am quoting a part of a previous answer as it is directly relevant to your question:
      "You would then normally choose the maximum positive value as input to the 'Norm of Magnetic or Gravitational field'.

      But this norm value is not critical if all you want to calculate later on is a heading (if it is a magnetometer) or a tilt angle (if it is an accelerometer). You can input any reasonable value for the norm, the correction matrix will be different by just a scaling factor, but the calculated heading (or tilt angle) will be the same, as it depends only on the relative value of the field components. The bias values will be unchanged, as they do not depend on the norm."
      In other words, inspect your output values and use the largest, which will likely represents the field magnitude when the sensor is aligned with it.

      Delete
  9. I am messing around with your C-based calibration tool magneto and am trying to feed it some simple data such as
    1 0 0
    0 1 0
    0 0 1
    -1 0 0
    0 -1 0
    0 0 -1
    .707106781 .707106781 0
    0 .707106781 .707106781
    .707106781 0 .707106781

    I thought these are sampled on a perfect sphere. As expected I get a bias of 0 in all three axis. I would have expected to obtain an A^-1 as the identity. This is not so. What am I doing wrong?

    ReplyDelete
    Replies
    1. Ignore me. You actually do get (almost) the identity. Since the data is on a sphere, the orientation for the elipsoid is not well defined. So it appears to be rotated by some (arbitrary?) orientation but still gives you the expected values h_cal for a choice if input h. A further mistake I made was to ignore the variable hm in the code which encodes for the above mentioned 'norm of magnetic or gravitational field'.

      Delete
  10. i really want to do the same stuff but on the C# can you give me the algorithm of your code it have been 2 weeks i'm working on it and didn't found anything yet.

    Can you help please

    ReplyDelete
  11. I got A-1 and b in the formula, how can I use it in a real application? Is there a related tutorial?

    ReplyDelete
  12. This comment has been removed by the author.

    ReplyDelete
  13. I'm interested in your 100% C-language

    ReplyDelete
  14. I'm interested in your 100% C-language port of the Magneto code, but keep getting errors at the Google URL. Is there any other way to see the file(s)?

    ReplyDelete
    Replies
    1. I am also interested and getting the error back from the link.

      Delete