Magnetic Declination by odder
http://en.wikipedia.org/wiki/Magnetic_declination#mediaviewer/File:Magnetic_declination.svg
The short story is that I don’t have enough scientific knowledge to know how to calculate magnetic declination, but enough to know that it can be done, using just a position on the surface of the earth and some information about the earth’s magnetic field. This is such an important topic that I was pretty sure there was open source code out there to do it, and I was right.
I found a pretty excellent implementation in C by Edward A Williams (right here) that is being maintained with the most recent magnetic field data, and ported it to C#. I present the port in the spirit of the GPL License under which the code was originally written. The nicest thing of all is that this code doesn’t depend on anything except system libraries, so you should be able to integrate it into almost anything.
The code was tested in current implementations of TrueNorth, but I’d appreciate any feedback in the comments or you can contact me directly.
I tried to Use your class, but was given a totally different answer to the known correct answer. I may be doing something wrong though. Are you still maintaining this code?inf
Yes I am maintaining and using it.
What are you using as your “known correct” reference? I only ask because declination varies over time, and if you’re using a map that’s a few years old it will have changed.
How far off is it from the “known correct” answer? The code uses a global magnetic model, and will have some error. It also uses the global magnetic field, and it is quite possible to have local deposits of magnetic material that can change the declination. These will be recorded on local maps, but cannot be predicted using a global model like this.
The final things to watch for are that the latitude and longitude are in radians, supplying degrees will result in completely meaningless values, and that the date is correct and converted to Julian days using the supplied method yymmdd_to_julian_days (this was also ported from the c code)
The code is a port of http://williams.best.vwh.net/magvar/magfield.c, and I tested it against the original for completeness.
I also tested it against http://geomag.nrcan.gc.ca/calc/mdcal-eng.php as a second source.
Other than reading your code, this is the only support I can give.
Thanks for the reply. You give me hope! :) If you check the site I am playing with, at listerhome.com, you will see what I am getting. The result of the Magnetic Declination is (temporally) shown UNDER the google map that is displayed. So, if you put in your address, you will see. I believe my result should be around 10.3, but I am getting -2.5. I am 93% sure I am doing something wrong. I am passing the location is decimal value… so, if I was at the Statue of Libery (Using the site I just mentioned), I get:
Location Decimal: 40.6891988, -74.0445167
Those numbers are then used in your code to try get the declination, which I am getting as -3.04092875392727, which I think is wrong.
If I understand correctly, you’re passing in the degrees as you wrote them in the comment.
You need to convert to radians:
double DegToRad(double angle) { return Math.PI * angle / 180.0; }
The output of the function will be in radians, so you’ll need to convert back to degrees
double RadToDeg(double angle) { return angle * (180.0 / Math.PI); }
All trigonometric functions (Sin, Cos, Tan, ATan etc) in C# (and many other languages) take radians instead of degrees as parameters.
Thanks Michael – I am getting a very close match now. Australian Geomagnetic Reference Field Computation says my house should have a value of 10.918 deg, and your code gives me 10.8959717215393. I use 1km height, and the date of today. Maybe they use a slightly different date. But that’s as close as damn it. Thanks for the good work, sir!
Glad it’s working. There may be differences in the geomagnetic models, changing the model parameter might bring it more in line with the one you are using for reference.
Excellent work, this is really great and will save me a ton of time!
One question though, what is the field parameter for in the SGMagVar method?
What should I be putting into that. The rest of it is very straightforward. Thanks!
Alright, I looked at the code and noticed that the array for field is being filled inside the method.
So, I just created an empty double array with 6 spots and dropped that into the method.
Is this the correct way to use the function? and is that array used so that you can see the different components of the equation?
Steve; sorry for not replying earlier.
By passing an appropriately sized array in the field parameter the method will return coefficients that will allow you to calculate magnetic field components at the location specified, according to the magnetic field model you have chosen.
I do not know enough about Geophysics to understand what this can be used for, as a person who writes GIS software I am only interested in the declination, which is a measure of how far the needle will be perturbed when held tangent to the earth’s surface. Perhaps those values are useful for something like tuning a three dimensional compass, or guiding missiles. Who knows.
Michael,
Excellent, thanks very much for the reply and thanks again for the code.
It is working perfectly for me now. Like you, I’m only concerned about the declination at any given point, but good to know that the other stuff is available if needed for a future project.
-Steve
Hello,
case 11: /* IGRF2010 */
yearfrac = (dat – yymmdd_to_julian_days(10, 1, 1)) / 365.25;
nmaxl = 14; // Results in index out of bounds? Andreas
Just a note that the IGRF2010 model is not working. Setting nmaxl = 13 won’t work. Setting nmaxl = 12 feels wrong, since the matrices iterated over are 14×14. No time to solve this tonight, so I’m going with the WMM2010 model. Maybe there is an easy fix?
Andreas, I am using WMM2010 as well and have not tested IGRF2010. The bug could be in the original code, but I will circle back and take a look at it in the next few weeks.
Hello, Can this be used to find the geomagnetic latitude?
I’d never heard of Geomagnetic Latitude before you asked this.
I think it is possible to use this to calculate geomagnetic latitude, but the math is beyond me at the moment.
It would be far easier to use the position of the magnetic north pole, calculate the transform between that and the geographic pole, and use that to calculate your latitude with respect to the new pole.
I’m sorry I don’t have the time and resources to figure this out, sounds like a cool problem.
I assume you want this in order to see the Aurora?
Hi Michale, I would like to use the latest wmm2015-model. How did you generate the gnm_wmm2010, hnm_wmm2010, gtnm_wmm2010, htnm_wmm2010-fields from the WMM.COE-file? Thanks Tom
its ok, figured it out, its all in here:
https://ngdc.noaa.gov/geomag/WMM/data/WMM2015/WMM2015_Report.pdf
Thank you for sharing the c# code
Hi Tom
Traveling for business in South Korea right now but I will look at this when I hey back.
If you’ve found a way to do this, I’d appreciate it if you could contribute it to this repo.