Index: trunk/SophyaExt/XAstroPack/Makefile
===================================================================
--- trunk/SophyaExt/XAstroPack/Makefile	(revision 1671)
+++ trunk/SophyaExt/XAstroPack/Makefile	(revision 1672)
@@ -3,7 +3,7 @@
 all: $(LIB)libXAstroPack.a
 clean:
-	rm -f $(OBJ)xastropack.o
+	rm -f $(OBJ)xastropack.o $(OBJ)magfield.o
 	rm -f $(LIB)libXAstroPack.a
-$(LIB)libXAstroPack.a : $(OBJ)xastropack.o
+$(LIB)libXAstroPack.a : $(OBJ)magfield.o $(OBJ)xastropack.o
 	$(ARCXX) $(ARCXXFLAGS) $@ $?
 ifeq ($(CXX),cxx)
@@ -12,5 +12,5 @@
  
 $(OBJ)xastropack.o: xastropack.cc xastropack.h \
- $(INC)machdefs.h \
- $(EXTINCPATH)XAstro/P_.h \
- $(EXTINCPATH)XAstro/astro.h
+  $(INC)machdefs.h \
+  $(EXTINCPATH)XAstro/P_.h \
+  $(EXTINCPATH)XAstro/astro.h
Index: trunk/SophyaExt/XAstroPack/magfield.c
===================================================================
--- trunk/SophyaExt/XAstroPack/magfield.c	(revision 1672)
+++ trunk/SophyaExt/XAstroPack/magfield.c	(revision 1672)
@@ -0,0 +1,781 @@
+/*  module magfield.c */
+
+/* Module to calculate magnetic variation and field given position,
+**               altitude, and date
+** Implements the NIMA (formerly DMA) WMM and IGRF models
+**
+**    http://www.nima.mil/GandG/ngdc-wmm2000.html
+**    For WMM2000 coefficients:
+**    ftp://ftp.ngdc.noaa.gov/Solid_Earth/Mainfld_Mag/DoD_Model/wmm.cof
+**    For IGRF/DGRF coefficients:
+**    http://swdcdb.kugi.kyoto-u.ac.jp/igrf/coef/igrfall.d
+**
+** Copyright (C) 2000  Edward A Williams <Ed_Williams@compuserve.com>
+**
+**  The routine uses a spherical harmonic expansion of the magnetic
+** potential up to twelfth order, together with its time variation, as
+** described in Chapter 4 of "Geomagnetism, Vol 1, Ed. J.A.Jacobs,
+** Academic Press (London 1987)". The program first converts geodetic
+** coordinates (lat/long on elliptic earth and altitude) to spherical
+** geocentric (spherical lat/long and radius) coordinates. Using this,
+** the spherical (B_r, B_theta, B_phi) magnetic field components are
+** computed from the model. These are finally referred to surface (X, Y,
+** Z) coordinates.
+**
+**   Fields are accurate to better than 200nT, variation and dip to
+** better than 0.5 degrees, with the exception of the declination near
+** the magnetic poles (where it is ill-defined) where the error may reach
+** 4 degrees or more.
+**
+**   Variation is undefined at both the geographic and  
+** magnetic poles, even though the field itself is well-behaved. To
+** avoid the routine blowing up, latitude entries corresponding to
+** the geographic poles are slightly offset. At the magnetic poles,
+** the routine returns zero variation.
+**
+** HISTORY
+** Adapted from EAW Excel 3.0 version 3/27/94 EAW
+** Recoded in C++ by Starry Chan
+** WMM95 added and rearranged in ANSI-C EAW 7/9/95
+** Put shell around program and made Borland & GCC compatible EAW 11/22/95
+** IGRF95 added 2/96 EAW
+** WMM2000 IGR2000 added 2/00 EAW
+** Released under GPL  3/26/00 EAW
+** Adaptions and modifications for the SimGear project  3/27/2000 CLO
+** Removed all pow() calls and made static roots[][] arrays to
+** save many sqrt() calls on subsequent invocations
+** 3/28/2000  Norman Vine -- nhv@yahoo.com
+** Put in some bullet-proofing to handle magnetic and geographic poles.
+** 3/28/2000 EAW
+** Added missing comment close, the lack of which caused the altitude 
+** correction to be omitted.
+** 01/31/01 Jim Seymour (jseymour@LinxNet.com)
+**
+** This program is free software; you can redistribute it and/or
+** modify it under the terms of the GNU General Public License as
+** published by the Free Software Foundation; either version 2 of the
+** License, or (at your option) any later version.
+**
+** This program is distributed in the hope that it will be useful, but
+** WITHOUT ANY WARRANTY; without even the implied warranty of
+** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+** General Public License for more details.
+**
+** You should have received a copy of the GNU General Public License
+** along with this program; if not, write to the Free Software
+** Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+**
+*/ 
+
+#include	<stdio.h>
+#include	<stdlib.h>
+#include	<math.h>
+#include        "magfield.h"
+
+#define	max(a,b)	(((a) > (b)) ? (a) : (b))
+
+static const double pi = 3.14159265358979;
+static const double a = 6378.16;	/* major radius (km) IAU66 ellipsoid */
+static const double f = 1.0 / 298.25;	/* inverse flattening IAU66 ellipsoid */
+static const double b = 6378.16 * (1.0 -1.0 / 298.25 );
+	/* minor radius b=a*(1-f) */
+static const double r_0 = 6371.2;	/* "mean radius" for spherical harmonic expansion */
+
+/*IGRF90 constants */
+static const double gnm_igrf90[13][13] = 
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{-29775.4, -1851, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{-2135.8, 3058.2, 1693.2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{1314.6, -2240.2, 1245.6, 806.5, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{938.9, 782.3, 323.9, -422.7, 141.7, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{-211, 352.5, 243.8, -110.8, -165.6, -37, 0, 0, 0, 0, 0, 0, 0}, 
+	{60.7, 63.9, 60.4, -177.5, 2, 16.7, -96.3, 0, 0, 0, 0, 0, 0}, 
+	{76.6, -64.2, 3.7, 27.5, 0.9, 5.7, 9.8, -0.5, 0, 0, 0, 0, 0}, 
+	{22.4, 5.1, -0.9, -10.8, -12.4, 3.8, 3.8, 2.6, -6, 0, 0, 0, 0}, 
+	{4.4, 9.9, 0.8, -12, 9.3, -3.9, -1.4, 7.3, 1.5, -5.5, 0, 0, 0}, 
+	{-3.6, -3.9, 2.4, -5.3, -2.4, 4.4, 3, 1.2, 2.2, 2.9, 0, 0, 0},
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}};
+static const double hnm_igrf90[13][13] = 
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{0, 5410.9, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -2277.7, -380, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -286.5, 293.3, -348.5, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 248.1, -239.5, 87, -299.4, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 47.2, 153.5, -154.4, -69.2, 97.7, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -15.8, 82.7, 68.3, -52.5, 1.8, 26.9, 0, 0, 0, 0, 0, 0},
+	{0, -81.1, -27.3, 0.6, 20.4, 16.4, -22.6, -5, 0, 0, 0, 0, 0}, 
+	{0, 9.7, -19.9, 7.1, -22.1, 11.9, 11, -16, -10.7, 0, 0, 0, 0}, 
+	{0, -20.8, 15.4, 9.5, -5.7, -6.4, 8.6, 9.1, -6.6, 1.9, 0, 0, 0}, 
+	{0, 1.3, 0.4, 3.1, 5.6, -4.2, -0.5, -1.5, 3.8, -0.5, -6.2, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}};
+static const double gtnm_igrf90[13][13] = 
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{18, 10.6, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{-12.9, 2.4, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{3.3, -6.7, 0, -5.9, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0.5, 0.6, -7, 0.5, -5.5, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0.6, -0.1, -1.6, -3.1, 0, 2.3, 0, 0, 0, 0, 0, 0, 0}, 
+	{1.3, -0.2, 1.8, 1.3, -0.2, 0.1, 1.2, 0, 0, 0, 0, 0, 0}, 
+	{0.6, -0.5, -0.3, 0.6, 1.6, 0.2, 0.2, 0.3, 0, 0, 0, 0, 0}, 
+	{0.2, -0.7, -0.2, 0.1, -1.1, 0, 0, -0.5, -0.6, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}};
+static const double htnm_igrf90[13][13] =
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{0, -16.1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -15.8, -13.8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 4.4, 1.6, -10.6, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 2.6, 1.8, 3.1, -1.4, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -0.1, 0.5, 0.4, 1.7, 0.4, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0.2, -1.3, 0, -0.9, 0.5, 1.2, 0, 0, 0, 0, 0, 0}, 
+	{0, 0.6, 0.2, 0.8, -0.5, -0.2, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0.5, -0.2, 0.3, 0.3, 0.4, -0.5, -0.3, 0.6, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}};
+/* IGRF95 */
+static const double gnm_igrf95[13][13] =
+{
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+  {-29682.0, -1789.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+  {-2197.0, 3074.0, 1685.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+  {1329.0, -2268.0, 1249.0, 769.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+  {941.0, 782.0, 291.0, -421.0, 116.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+  {-210.0, 352.0, 237.0, -122.0, -167.0, -26.0, 0.0, 0.0, 0.0, 0.0,0.0,0.0,0.0},
+  {66.0, 64.0, 65.0, -172.0, 2.0, 17.0, -94.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {78.0, -67.0, 1.0, 29.0, 4.0, 8.0, 10.0, -2.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {24.0, 4.0, -1.0, -9.0, -14.0, 4.0, 5.0, 0.0, -7.0, 0.0,0.0, 0.0, 0.0},
+  {4.0, 9.0, 1.0, -12.0, 9.0, -4.0, -2.0, 7.0, 0.0, -6.0, 0.0, 0.0, 0.0},
+  {-3.0, -4.0, 2.0, -5.0, -2.0, 4.0, 3.0, 1.0, 3.0, 3.0, 0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0}
+ };
+
+static const double hnm_igrf95[13][13]=
+{
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 5318.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, -2356.0, -425.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, -263.0, 302.0, -406.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 262.0, -232.0, 98.0, -301.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 44.0, 157.0, -152.0, -64.0, 99.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, -16.0, 77.0, 67.0, -57.0, 4.0, 28.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, -77.0, -25.0, 3.0, 22.0, 16.0, -23.0, -3.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 12.0, -20.0, 7.0, -21.0, 12.0, 10.0, -17.0, -10.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, -19.0, 15.0, 11.0, -7.0, -7.0, 9.0, 7.0, -8.0, 1.0, 0.0, 0.0, 0.0},
+  {0.0, 2.0, 1.0, 3.0, 6.0, -4.0, 0.0, -2.0, 3.0, -1.0, -6.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0}
+};
+
+static const double gtnm_igrf95[13][13]=
+{
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {17.6, 13.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {-13.2, 3.7, -0.8, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {1.5, -6.4, -0.2, -8.1, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.8, 0.9, -6.9, 0.5, -4.6, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.8, 0.1, -1.5, -2.0, -0.1, 2.3, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.5, -0.4, 0.6, 1.9, -0.2, -0.2, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {-0.2, -0.8, -0.6, 0.6, 1.2, 0.1, 0.2, -0.6, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.3, -0.2, 0.1, 0.4, -1.1, 0.3, 0.2, -0.9, -0.3, 0.0,0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0}
+ };
+
+static const double htnm_igrf95[13][13]=
+{
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, -18.3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, -15.0, -8.8, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 4.1, 2.2, -12.1, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 1.8, 1.2, 2.7, -1.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.2, 1.2, 0.3, 1.8, 0.9, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.3, -1.6, -0.2, -0.9, 1.0, 2.2, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.8, 0.2, 0.6, -0.4, 0.0, -0.3, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.4, -0.2, 0.2, 0.7, 0.0, -1.2, -0.7, -0.6, 0.0,0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0},
+  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,0.0, 0.0, 0.0, 0.0}
+};
+
+/*WMM85 constants */
+static const double gnm_wmm85[13][13] = 
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{-29879.8, -1903.3, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{-2070.6, 3040.8, 1696.7, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{1303.9, -2203, 1241.7, 839.4, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{933.8, 781.8, 359, -424.5, 164.5, 0, 0, 0, 0, 0, 0, 0, 0},
+	{-216.4, 353, 254.3, -93.7, -157.5, -45.2, 0, 0, 0, 0, 0, 0, 0}, 
+	{53.2, 63.8, 51.3, -188.4, 3.3, 20.3, -101.7, 0, 0, 0, 0, 0, 0}, 
+	{76.9, -60.7, 0.7, 25.4, -8.1, 6.9, 7, -4.4, 0, 0, 0, 0, 0}, 
+	{18.4, 5.1, 1.2, -12, -9.1, 0.1, 4.7, 6.5, -9.5, 0, 0, 0, 0}, 
+	{5.7, 10.9, 0.9, -12.2, 9.5, -3.3, -1, 6.5, 1.5, -4.8, 0, 0, 0}, 
+	{-3.4, -4.7, 2.5, -5.5, -2.1, 4.6, 3.2, 0.6, 1.9, 2.8, -0.2, 0, 0}, 
+	{2.3, -0.8, -2, 2.1, 0.2, -0.4, -0.4, 1.6, 1.5, -0.7, 2.3, 3.5, 0}, 
+	{-1.8, 0, 0.1, -0.3, 0.5, 0.5, -0.6, -0.4, 0, -0.5, 0, 0.7, -0.2}};
+static const double hnm_wmm85[13][13] = 
+	{
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{0, 5490.5, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -2189.1, -312, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -310.3, 282.6, -299.2, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 227.2, -246.7, 72.5, -299.1, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 43.4, 148.2, -154.8, -71.8, 91.5, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -12.3, 87.9, 67.8, -51.1, -4, 20.8, 0, 0, 0, 0, 0, 0}, 
+	{0, -80.1, -25.9, -0.9, 21.6, 18.5, -20, -7.7, 0, 0, 0, 0, 0}, 
+	{0, 3.8, -20.2, 5, -24.2, 12.2, 7.6, -16.3, -10.9, 0, 0, 0, 0}, 
+	{0, -20.8, 15.8, 9, -5, -6.4, 9.1, 9.9, -5.8, 2.3, 0, 0, 0}, 
+	{0, 1.2, 0.4, 2.5, 5.6, -4.4, -0.5, -1.6, 3.7, -0.5, -6.1, 0, 0}, 
+	{0, 1.3, 2, -1.1, -2.8, 0.7, -0.1, -2.4, -0.4, -1.5, -1.5, 0.7, 0},
+	{0, 0.3, 0.6, 2.5, -1.7, 0.3, 0.2, -0.1, 0.1, 0.1, -1.4, 0.4, 0.7}};
+static const double gtnm_wmm85[13][13] = 
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{21.9, 10.6, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{-11.2, 1.8, 9.3, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{8.3, -2, -0.6, 2.4, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{-1.2, 0.1, -9.7, -1.7, -9.3, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{1.4, -0.5, -1.2, -2.2, 0.9, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{3.1, 0, 1.8, -0.2, -0.4, 2.4, 1.8, 0, 0, 0, 0, 0, 0}, 
+	{-0.1, -0.8, -1.2, 1.1, 0, 0.6, -1.8, -1.2, 0, 0, 0, 0, 0},
+	{0.2, 0, 0.7, 0.1, 0.2, -0.3, -0.1, 0.2, -2.2, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}};
+static const double htnm_wmm85[13][13] = 
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{0, -31.5, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -9.7, -19.9, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 6.1, 1.3, -13, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{0, 1.3, 3.6, 2.5, 0.6, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -0.9, 0.6, 0.3, 2.4, -1.4, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0.7, -2.1, -1.4, -4.3, -0.7, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 1.2, 2, 2.6, 0.9, 0.8, 0.4, 0, 0, 0, 0, 0}, 
+	{0, -0.6, -1.5, 0.1, -1.1, 0.4, -2, 0.9, 1.5, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}};
+
+/* wmm90 constants */
+static const double gnm_wmm90[13][13] =
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{-29780.5, -1851.7, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{-2134.3, 3062.2, 1691.9, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{1312.9, -2244.7, 1246.8, 808.6, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{933.5, 784.9, 323.5, -421.7, 139.2, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{-208.3, 352.2, 246.5, -110.8, -162.3, -37.2, 0, 0, 0, 0, 0, 0, 0}, 
+	{59, 63.7, 60, -181.3, 0.4, 15.4, -96, 0, 0, 0, 0, 0, 0}, 
+	{76.1, -62.1, 1.3, 30.2, 4.7, 7.9, 10.1, 1.9, 0, 0, 0, 0, 0}, 
+	{22.9, 2.3, -1.2, -11.7, -17.5, 2.2, 5.7, 3, -7, 0, 0, 0, 0}, 
+	{3.6, 9.5, -0.9, -10.7, 10.7, -3.2, -1.4, 6.3, 0.8, -5.5, 0, 0, 0},
+	{-3.3, -2.6, 4.5, -5.6, -3.6, 3.9, 3.2, 1.7, 3, 3.7, 0.7, 0, 0}, 
+	{1.3, -1.4, -2.5, 3.2, 0.2, -1.1, 0.3, -0.3, 0.9, -1.1, 2.4, 3, 0}, 
+	{-1.3, 0.1, 0.5, 0.7, 0.4, -0.2, -1.1, 0.9, -0.6, 0.8, 0.2, 0.4, 0.2}};
+static const double hnm_wmm90[13][13] = 
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{0, 5407.2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -2278.3, -384.3, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -284.9, 291.7, -352.4, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 249.4, -232.7, 91.3, -296.5, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 40.8, 148.7, -154.6, -67.6, 97.4, 0, 0, 0, 0, 0, 0, 0},
+	{0, -14.7, 82.2, 70, -56.2, -1.4, 24.6, 0, 0, 0, 0, 0, 0}, 
+	{0, -78.6, -26.7, 0.1, 19.9, 17.9, -21.5, -6.8, 0, 0, 0, 0, 0}, 
+	{0, 9.7, -19.3, 6.6, -20.1, 13.4, 9.8, -19, -9.1, 0, 0, 0, 0}, 
+	{0, -21.9, 14.3, 9.5, -6.7, -6.4, 9.1, 8.9, -8, 2.1, 0, 0, 0}, 
+	{0, 2.6, 1.2, 2.6, 5.7, -4, -0.4, -1.7, 3.8, -0.8, -6.5, 0, 0}, 
+	{0, 0, 1, -1.6, -2.2, 1.1, -0.7, -1.7, -1.5, -1.3, -1.1, 0.6, 0}, 
+	{0, 0.7, 0.7, 1.3, -1.5, 0.3, 0.2, -1.1, 1.2, -0.2, -1.3, 0.6, 0.6}};
+static const double gtnm_wmm90[13][13] = 
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{16, 9.3, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{-11.7, 3.7, 1.8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{2.1, -7.6, 0, -5.8, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{-0.8, 1, -7.4, 0.8, -6.4, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{1.7, 0, 0, -2.7, 0, 3, 0, 0, 0, 0, 0, 0, 0}, 
+	{0.8, 0, 1.5, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0.5, 0, -0.9, 1.5, 2.7, -1, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -1.1, 0, 0, -2.1, 0, 1, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}};
+static const double htnm_wmm90[13][13] = 
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{0, -13.8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -12.8, -14.9, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 3.1, 0.8, -11.3, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 3.3, 3.7, 2.8, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, -2.1, 1.2, 1.2, 0.6, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, -0.6, -0.6, 0, -2.3, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0.6, 0.8, 0, 0, 0, 0.4, 0, 0, 0, 0, 0, 0}, 
+	{0, 0.4, -0.8, 0.5, 0.3, 0.5, 0, -0.7, 0, 0, 0, 0, 0},
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, 
+	{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}};
+/* wmm95 constants */
+static const double gnm_wmm95[13][13] =
+        {
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{-29682.1,-1782.2,0,0,0,0,0,0,0,0,0,0,0},
+        {-2194.7,3078.6,1685.7,0,0,0,0,0,0,0,0,0,0},
+	{1318.8,-2273.6,1246.9,766.3,0,0,0,0,0,0,0,0,0},
+        {940,782.9,290.9,-418.9,113.8,0,0,0,0,0,0,0,0},
+        {-209.5,354,238.2,-122.1,-162.8,-23.3,0,0,0,0,0,0,0},
+        {68.5,65.6,64.1,-169.1,-0.5,16.5,-91,0,0,0,0,0,0},
+        {78,-68.1,0.1,29.6,6,8.7,9.2,-2.4,0,0,0,0,0},
+        {24.7,3.4,-1.5,-9.6,-16.5,2.6,3.6,-4.9,-8.5,0,0,0,0},
+        {2.9,7.5,0.4,-10.3,9.7,-2.3,-2.4,6.8,-0.5,-6.5,0,0,0},
+        {-2.9,-3.3,2.8,-4.3,-3.1,2.4,2.8,0.7,4.1,3.6,0.6,0,0},
+        {1.7,-1.6,-3.6,1.2,-0.6,0.1,-0.7,-0.8,1.3,-0.3,2.2,4.2,0},
+        {-1.8,0.9,-0.1,-0.5,0.8,0.2,0.5,0.4,-0.4,0.3,0.2,0.4,0.6}};
+static const double hnm_wmm95[13][13] =
+	{
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+        {0,5315.6,0,0,0,0,0,0,0,0,0,0,0},
+        {0,-2359.1,-418.6,0,0,0,0,0,0,0,0,0,0},
+        {0,-261.1,301,-416.5,0,0,0,0,0,0,0,0,0},
+        {0,259.4,-230.9,99.8,-306.1,0,0,0,0,0,0,0,0},
+        {0,43.7,157.6,-150.1,-59.2,104.4,0,0,0,0,0,0,0},
+        {0,-15.2,74.3,69.4,-55.3,3,33.3,0,0,0,0,0,0},
+        {0,-76.1,-24.5,1.6,20,16.5,-23.6,-6.8,0,0,0,0,0},
+        {0,14.9,-19.5,6.3,-20.4,12.2,7,-19,-8.8,0,0,0,0},
+        {0,-19.8,14.6,10.9,-7.5,-6.8,9.3,7.7,-8.1,2.6,0,0,0},
+	{0,3.2,1.7,2.9,5.6,-3.4,-0.7,-2.9,2.3,-1.6,-6.6,0,0},
+        {0,0.3,1,-3.6,-1.4,1.9,0.2,-1.3,-2.4,-0.6,-2.2,1.3,0},
+        {0,0.3,1.4,0.8,-3,0.7,0.5,-0.8,0.6,0.1,-1.3,-0.4,0.9}};
+static const double gtnm_wmm95[13][13] =
+        {
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+        {17.6,13.2,0,0,0,0,0,0,0,0,0,0,0},
+        {-13.7,4,-0.3,0,0,0,0,0,0,0,0,0,0},
+        {0.8,-6.6,-0.5,-8.5,0,0,0,0,0,0,0,0,0},
+        {1.2,1.1,-6.8,0.3,-4.5,0,0,0,0,0,0,0,0},
+        {0.9,0.5,-1.4,-1.7,0,2.1,0,0,0,0,0,0,0},
+	{0.4,-0.3,0.3,2.1,0,-0.4,-0.4,0,0,0,0,0,0},
+        {-0.3,-1.1,-0.5,0.5,1.3,0.1,0,-0.9,0,0,0,0,0},
+        {0.1,0,0.4,0.3,-1.3,0.5,0.4,-0.9,0.1,0,0,0,0},
+        {0,0,0,0,0,0,0,0,0,0,0,0,0},
+        {0,0,0,0,0,0,0,0,0,0,0,0,0},
+        {0,0,0,0,0,0,0,0,0,0,0,0,0},
+        {0,0,0,0,0,0,0,0,0,0,0,0,0}};
+static const double htnm_wmm95[13][13] =
+        {
+        {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
+	{0,-18,0,0,0,0,0,0,0,0,0,0,0},
+	{0,-14.6,-7.2,0,0,0,0,0,0,0,0,0,0},
+        {0,4,2.2,-12.6,0,0,0,0,0,0,0,0,0},
+        {0,1.3,1,2.5,-1.2,0,0,0,0,0,0,0,0},
+        {0,0.5,1.5,0.6,1.7,0.6,0,0,0,0,0,0,0},
+        {0,0.7,-1.5,-0.5,-0.7,1.1,2.6,0,0,0,0,0,0},
+        {0,0.3,0,0.7,-0.6,0.1,-0.6,-0.4,0,0,0,0,0},
+        {0,0.4,-0.3,0.1,0.8,-0.1,-1.3,-0.9,-1.1,0,0,0,0},
+        {0,0,0,0,0,0,0,0,0,0,0,0,0},
+        {0,0,0,0,0,0,0,0,0,0,0,0,0},
+        {0,0,0,0,0,0,0,0,0,0,0,0,0},
+        {0,0,0,0,0,0,0,0,0,0,0,0,0}};
+
+static const double gnm_wmm2000[13][13] =
+        {
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-29616.0, -1722.7, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-2266.7, 3070.2, 1677.6, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{1322.4, -2291.5, 1255.9, 724.8, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{932.1, 786.3, 250.6, -401.5, 106.2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-211.9, 351.6, 220.8, -134.5, -168.8, -13.3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{73.8, 68.2, 74.1, -163.5, -3.8, 17.1, -85.1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{77.4, -73.9, 2.2, 35.7, 7.3, 5.2, 8.4, -1.5, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{23.3, 7.3, -8.5, -6.6, -16.9, 8.6, 4.9, -7.8, -7.6, 0.0, 0.0, 0.0, 0.0},
+	{5.7, 8.5, 2.0, -9.8, 7.6, -7.0, -2.0, 9.2, -2.2, -6.6, 0.0, 0.0, 0.0},
+	{-2.2, -5.7, 1.6, -3.7, -0.6, 4.1, 2.2, 2.2, 4.6, 2.3, 0.1, 0.0, 0.0},
+	{3.3, -1.1, -2.4, 2.6, -1.3, -1.7, -0.6, 0.4, 0.7, -0.3, 2.3, 4.2, 0.0},
+	{-1.5, -0.2, -0.3, 0.5, 0.2, 0.9, -1.4, 0.6, -0.6, -1.0, -0.3, 0.3, 0.4},
+	};
+
+static const double hnm_wmm2000[13][13]=
+        {
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 5194.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -2484.8, -467.9, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -224.7, 293.0, -486.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 273.3, -227.9, 120.9, -302.7, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 42.0, 173.8, -135.0, -38.6, 105.2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -17.4, 61.2, 63.2, -62.9, 0.2, 43.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -62.3, -24.5, 8.9, 23.4, 15.0, -27.6, -7.8, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 12.4, -20.8, 8.4, -21.2, 15.5, 9.1, -15.5, -5.4, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -20.4, 13.9, 12.0, -6.2, -8.6, 9.4, 5.0, -8.4, 3.2, 0.0, 0.0, 0.0},
+	{0.0, 0.9, -0.7, 3.9, 4.8, -5.3, -1.0, -2.4, 1.3, -2.3, -6.4, 0.0, 0.0},
+	{0.0, -1.5, 0.7, -1.1, -2.3, 1.3, -0.6, -2.8, -1.6, -0.1, -1.9, 1.4, 0.0},
+	{0.0, -1.0, 0.7, 2.2, -2.5, -0.2, 0.0, -0.2, 0.0, 0.2, -0.9, -0.2, 1.0},
+	};
+
+static const double gtnm_wmm2000[13][13]=
+        {
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{14.7, 11.1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-13.6, -0.7, -1.8, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.3, -4.3, 0.9, -8.4, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-1.6, 0.9, -7.6, 2.2, -3.2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-0.9, -0.2, -2.5, -2.7, -0.9, 1.7, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{1.2, 0.2, 1.7, 1.6, -0.1, -0.3, 0.8, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-0.4, -0.8, -0.2, 1.1, 0.4, 0.0, -0.2, -0.2, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-0.3, 0.6, -0.8, 0.3, -0.2, 0.5, 0.0, -0.6, 0.1, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	};
+
+static const double htnm_wmm2000[13][13]=
+        {
+        {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -20.4, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -21.5, -9.6, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 6.4, -1.3, -13.3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 2.3, 0.7, 3.7, -0.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 2.1, 2.3, 3.1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -0.3, -1.7, -0.9, -1.0, -0.1, 1.9, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 1.4, 0.2, 0.7, 0.4, -0.3, -0.8, -0.1, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -0.5, 0.1, -0.2, 0.0, 0.1, -0.1, 0.3, 0.2, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+        };
+
+static const double gnm_igrf2000[13][13] =
+        {
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-29615.0, -1728.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-2267.0, 3072.0, 1672.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{1341.0, -2290.0, 1253.0, 715.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{935.0, 787.0, 251.0, -405.0, 110.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-217.0, 351.0, 222.0, -131.0, -169.0, -12.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{72.0, 68.0, 74.0, -161.0, -5.0, 17.0, -91.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{79.0, -74.0, 0.0, 33.0, 9.0, 7.0, 8.0, -2.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{25.0, 6.0, -9.0, -8.0, -17.0, 9.0, 7.0, -8.0, -7.0, 0.0, 0.0, 0.0, 0.0},
+	{5.0, 9.0, 3.0, -8.0, 6.0, -9.0, -2.0, 9.0, -4.0, -8.0, 0.0, 0.0, 0.0},
+	{-2.0, -6.0, 2.0, -3.0, 0.0, 4.0, 1.0, 2.0, 4.0, 0.0, -1.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	};
+
+static const double hnm_igrf2000[13][13]=
+        {
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 5186.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -2478.0, -458.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -227.0, 296.0, -492.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 272.0, -232.0, 119.0, -304.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 44.0, 172.0, -134.0, -40.0, 107.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -17.0, 64.0, 65.0, -61.0, 1.0, 44.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -65.0, -24.0, 6.0, 24.0, 15.0, -25.0, -6.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 12.0, -22.0, 8.0, -21.0, 15.0, 9.0, -16.0, -3.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -20.0, 13.0, 12.0, -6.0, -8.0, 9.0, 4.0, -8.0, 5.0, 0.0, 0.0, 0.0},
+	{0.0, 1.0, 0.0, 4.0, 5.0, -6.0, -1.0, -3.0, 0.0, -2.0, -8.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	};
+
+static const double gtnm_igrf2000[13][13]=
+        {
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{14.6, 10.7, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-12.4, 1.1, -1.1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.7, -5.4, 0.9, -7.7, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-1.3, 1.6, -7.3, 2.9, -3.2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -0.7, -2.1, -2.8, -0.8, 2.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{1.0, -0.4, 0.9, 2.0, -0.6, -0.3, 1.2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-0.4, -0.4, -0.3, 1.1, 1.1, -0.2, 0.6, -0.9, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{-0.3, 0.2, -0.3, 0.4, -1.0, 0.3, -0.5, -0.7, -0.4, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	};
+
+static const double htnm_igrf2000[13][13]=
+        {
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -22.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -20.6, -9.6, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 6.0, -0.1, -14.2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 2.1, 1.3, 5.0, 0.3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -0.1, 0.6, 1.7, 1.9, 0.1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, -0.2, -1.4, 0.0, -0.8, 0.0, 0.9, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 1.1, 0.0, 0.3, -0.1, -0.6, -0.7, 0.2, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.1, 0.0, 0.0, 0.3, 0.6, -0.4, 0.3, 0.7, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+	};
+
+
+
+static const int nmax = 12;
+
+double P[13][13];
+double DP[13][13];
+double gnm[13][13];
+double hnm[13][13];
+double sm[13];
+double cm[13];
+
+static double root[13];
+static double roots[13][13][2];
+
+/* Convert date to Julian day    1950-2049 */
+unsigned long int yymmdd_to_julian_days(int yy,int mm,int dd)
+{
+  unsigned long jd;
+ 
+  yy = (yy < 50) ? (2000 + yy) : (1900 + yy);
+  jd = dd - 32075L + 1461L * (yy + 4800L + (mm - 14) / 12 ) / 4;
+  jd = jd + 367L * (mm - 2 - (mm - 14) / 12*12) / 12;
+  jd = jd - 3 * ((yy + 4900L + (mm - 14) / 12) / 100) / 4;
+
+  return(jd);
+} 
+
+/* Convert degrees to radians */
+double deg_to_rad(double deg)
+{
+return deg*pi/180.;
+}
+
+/* Convert radians to degrees */
+double rad_to_deg(double rad)
+{
+return rad*180./pi;
+}
+	     
+
+
+
+/* 
+ * return variation (in radians) given geodetic latitude (radians), longitude
+ * (radians) ,height (km) and (Julian) date
+ * model=1 is IGRF90, 2 is WMM85, 3 is WMM90, 4 is WMM95, 
+ * 5 is IGRF95, 6 is WMM2000, 7 is IGRF2000
+ * N and E lat and long are positive, S and W negative
+*/
+
+double SGMagVar( double lat, double lon, double h, long dat, int model, double* field )
+{
+    /* output field B_r,B_th,B_phi,B_x,B_y,B_z */
+    int n,m;
+    double yearfrac,sr,r,theta,c,s,psi,fn,fn_0,B_r,B_theta,B_phi,X,Y,Z;
+    double sinpsi, cospsi, inv_s;
+
+    static int been_here = 0;
+	
+    double sinlat = sin(lat);
+    double coslat = cos(lat);
+
+    /* convert to geocentric */ 
+    sr = sqrt(a*a*coslat*coslat + b*b*sinlat*sinlat);
+    /* sr is effective radius */
+    theta = atan2(coslat * (h*sr + a*a), sinlat * (h*sr + b*b));
+		  
+    /* theta is geocentric co-latitude */
+
+    r = h*h + 2.0*h * sr +
+	(a*a*a*a - ( a*a*a*a - b*b*b*b ) * sinlat*sinlat ) / 
+	(a*a - (a*a - b*b) * sinlat*sinlat );
+
+    r = sqrt(r);
+
+    /* r is geocentric radial distance */
+    c = cos(theta);
+    s = sin(theta);
+      /* protect against zero divide at geographic poles */
+    inv_s =  1.0 / (s + (s == 0.)*1.0e-8); 
+
+    /*zero out arrays */
+    for ( n = 0; n <= nmax; n++ ) {
+      for ( m = 0; m <= n; m++ ) {
+	P[n][m] = 0;
+	DP[n][m] = 0;
+      }
+    }
+
+    /* diagonal elements */
+    P[0][0] = 1;
+    P[1][1] = s;
+    DP[0][0] = 0;
+    DP[1][1] = c;
+    P[1][0] = c ;
+    DP[1][0] = -s;
+
+    /* these values will not change for subsequent function calls */
+    if( !been_here ) {
+      for ( n = 2; n <= nmax; n++ ) {
+	root[n] = sqrt((2.0*n-1) / (2.0*n));
+      }
+	
+      for ( m = 0; m <= nmax; m++ ) {
+	double mm = m*m;
+	for ( n = max(m + 1, 2); n <= nmax; n++ ) {
+	  roots[m][n][0] = sqrt((n-1)*(n-1) - mm);
+	  roots[m][n][1] = 1.0 / sqrt( n*n - mm);
+	}
+      }
+      been_here = 1;
+    }
+	
+    for ( n=2; n <= nmax; n++ ) {
+      /*  double root = sqrt((2.0*n-1) / (2.0*n)); */
+      P[n][n] = P[n-1][n-1] * s * root[n];
+      DP[n][n] = (DP[n-1][n-1] * s + P[n-1][n-1] * c) * root[n];
+    }
+
+    /* lower triangle */
+    for ( m = 0; m <= nmax; m++ ) {
+      /*  double mm = m*m;  */
+      for ( n = max(m + 1, 2); n <= nmax; n++ ) {
+	/* double root1 = sqrt((n-1)*(n-1) - mm); */
+	/* double root2 = 1.0 / sqrt( n*n - mm);  */
+	P[n][m] = (P[n-1][m] * c * (2.0*n-1) -
+		   P[n-2][m] * roots[m][n][0]) * roots[m][n][1];
+	DP[n][m] = ((DP[n-1][m] * c - P[n-1][m] * s) *
+		    (2.0*n-1) - DP[n-2][m] * roots[m][n][0]) * roots[m][n][1];
+      }
+    }
+
+    /* compute gnm, hnm at dat */
+    switch(model) {
+    case 1:	/* IGRF90 */
+      yearfrac = (dat - yymmdd_to_julian_days(90,1,1)) / 365.25;
+      for (n=1;n<=nmax;n++)
+        for (m = 0;m<=nmax;m++) {
+	  gnm[n][m] = gnm_igrf90[n][m] + yearfrac * gtnm_igrf90[n][m];
+	  hnm[n][m] = hnm_igrf90[n][m] + yearfrac * htnm_igrf90[n][m];
+	}
+      break;
+
+    case 2:	/* WMM85 */
+      yearfrac = (dat - yymmdd_to_julian_days(85,1,1)) / 365.25;
+      for (n=1;n<=nmax;n++)
+        for (m = 0;m<=nmax;m++) {
+	  gnm[n][m] = gnm_wmm85[n][m] + yearfrac * gtnm_wmm85[n][m];
+	  hnm[n][m] = hnm_wmm85[n][m] + yearfrac * htnm_wmm85[n][m];
+	}
+      break;
+
+    case 3:	/* WMM90 */
+      yearfrac = (dat - yymmdd_to_julian_days(90,1,1)) / 365.25;
+      for (n=1;n<=nmax;n++)
+	for (m = 0;m<=nmax;m++) {
+	  gnm[n][m] = gnm_wmm90[n][m] + yearfrac * gtnm_wmm90[n][m];
+	  hnm[n][m] = hnm_wmm90[n][m] + yearfrac * htnm_wmm90[n][m];
+	}
+      break;
+
+    case 4:	/* WMM95 */
+      yearfrac = (dat - yymmdd_to_julian_days(95,1,1)) / 365.25;
+      for (n=1;n<=nmax;n++)
+        for (m = 0;m<=nmax;m++) {
+           gnm[n][m] = gnm_wmm95[n][m] + yearfrac * gtnm_wmm95[n][m];
+	   hnm[n][m] = hnm_wmm95[n][m] + yearfrac * htnm_wmm95[n][m];
+	}
+      break;
+    case 5:	/* IGRF95 */
+      yearfrac = (dat - yymmdd_to_julian_days(95,1,1)) / 365.25;
+      for (n=1;n<=nmax;n++)
+        for (m = 0;m<=nmax;m++) {
+           gnm[n][m] = gnm_igrf95[n][m] + yearfrac * gtnm_igrf95[n][m];
+	   hnm[n][m] = hnm_igrf95[n][m] + yearfrac * htnm_igrf95[n][m];
+	}
+      break;
+    case 6:      /* WMM2000 */
+      yearfrac = (dat - yymmdd_to_julian_days(0,1,1)) / 365.25;
+       for (n=1;n<=nmax;n++)
+         for (m = 0;m<=nmax;m++) {
+           gnm[n][m] = gnm_wmm2000[n][m] + yearfrac * gtnm_wmm2000[n][m];
+	   hnm[n][m] = hnm_wmm2000[n][m] + yearfrac * htnm_wmm2000[n][m];
+	 }
+       break;
+    case 7:      /* IGRF2000 */
+      yearfrac = (dat - yymmdd_to_julian_days(0,1,1)) / 365.25;
+      for (n=1;n<=nmax;n++)
+	for (m = 0;m<=nmax;m++) {
+	  gnm[n][m] = gnm_igrf2000[n][m] + yearfrac * gtnm_igrf2000[n][m];
+	  hnm[n][m] = hnm_igrf2000[n][m] + yearfrac * htnm_igrf2000[n][m];
+	}
+      break;  
+    }
+
+    /* compute sm (sin(m lon) and cm (cos(m lon)) */
+    for (m = 0;m<=nmax;m++) {
+      sm[m] = sin(m * lon);
+      cm[m] = cos(m * lon);
+    }
+
+    /* compute B fields */
+    B_r = 0.0;
+    B_theta = 0.0;
+    B_phi = 0.0;
+    fn_0 = r_0/r;
+    fn = fn_0 * fn_0;
+	
+    for ( n = 1; n <= nmax; n++ ) {
+      double c1_n=0;
+      double c2_n=0;
+      double c3_n=0;
+      for ( m = 0; m <= n; m++ ) {
+	double tmp = (gnm[n][m] * cm[m] + hnm[n][m] * sm[m]); 
+	c1_n += tmp * P[n][m];
+	c2_n += tmp * DP[n][m];
+	c3_n +=  m * (gnm[n][m] * sm[m] - hnm[n][m] * cm[m]) * P[n][m];
+      }
+      /* fn=pow(r_0/r,n+2.0);   */
+      fn *= fn_0;
+      B_r += (n + 1) * c1_n * fn;
+      B_theta -= c2_n * fn;
+      B_phi += c3_n * fn * inv_s;
+    }
+
+
+
+    /* Find geodetic field components: */
+    psi = theta - (pi / 2.0 - lat);
+    sinpsi = sin(psi);
+    cospsi = cos(psi);
+    X = -B_theta * cospsi - B_r * sinpsi;
+    Y = B_phi;
+    Z = B_theta * sinpsi - B_r * cospsi;
+
+    field[0]=B_r;
+    field[1]=B_theta;
+    field[2]=B_phi;
+    field[3]=X;
+    field[4]=Y;
+    field[5]=Z;   /* output fields */
+    /* find variation in radians */
+    /* return zero variation at magnetic pole X=Y=0. */
+    /* E is positive */
+    return (X != 0. || Y != 0.) ? atan2(Y, X) : (double) 0.; 
+}
+
Index: trunk/SophyaExt/XAstroPack/magfield.h
===================================================================
--- trunk/SophyaExt/XAstroPack/magfield.h	(revision 1672)
+++ trunk/SophyaExt/XAstroPack/magfield.h	(revision 1672)
@@ -0,0 +1,66 @@
+/* magfield.h -- compute local magnetic variation given position,
+**               altitude, and date
+**
+** This is an implementation of the NIMA WMM 2000
+**
+**    http://www.nima.mil/GandG/ngdc-wmm2000.html
+**    For WMM2000 coefficients:
+**    ftp://ftp.ngdc.noaa.gov/Solid_Earth/Mainfld_Mag/DoD_Model/wmm.cof
+**    For IGRF/DGRF coefficients:
+**    http://swdcdb.kugi.kyoto-u.ac.jp/igrf/coef/igrfall.d
+**
+** Copyright (C) 2000  Edward A Williams <Ed_Williams@compuserve.com>
+**
+** Adapted from Excel 3.0 version 3/27/94 EAW
+** Recoded in C++ by Starry Chan
+** WMM95 added and rearranged in ANSI-C EAW 7/9/95
+** Put shell around program and made Borland & GCC compatible EAW 11/22/95
+** IGRF95 added 2/96 EAW
+** WMM2000 IGR2000 added 2/00 EAW
+** Released under GPL 3/26/00 EAW
+** Adaptions and modifications for the SimGear project  3/27/2000 CLO
+**
+** This program is free software; you can redistribute it and/or
+** modify it under the terms of the GNU General Public License as
+** published by the Free Software Foundation; either version 2 of the
+** License, or (at your option) any later version.
+**
+** This program is distributed in the hope that it will be useful, but
+** WITHOUT ANY WARRANTY; without even the implied warranty of
+** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+** General Public License for more details.
+**
+** You should have received a copy of the GNU General Public License
+** along with this program; if not, write to the Free Software
+** Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+**
+*/
+#ifndef MAGFIELD_H
+#define MAGFIELD_H
+
+#ifdef __cplusplus
+extern "C" {  /* extern "C" */
+#endif
+
+/* Convert date to Julian day    1950-2049 */
+unsigned long int yymmdd_to_julian_days( int yy, int mm, int dd );
+
+/* Convert degrees to radians */
+double deg_to_rad( double deg );
+
+/* Convert radians to degrees */
+double rad_to_deg( double rad );
+
+/* return variation (in radians) given geodetic latitude (radians), longitude
+(radians) ,height (km), (Julian) date and field model
+model=1 is IGRF90, 2 is WMM85, 3 is WMM90, 4 is WMM95, 5 is IGRF95, 
+6 is WMM2000, 7 is IGRF2000
+N and E lat and long are positive, S and W negative
+*/
+double SGMagVar( double lat, double lon, double h, long dat, int model, double* field );
+
+#ifdef __cplusplus
+}             /* extern "C" */
+#endif
+
+#endif
Index: trunk/SophyaExt/XAstroPack/objlist.list
===================================================================
--- trunk/SophyaExt/XAstroPack/objlist.list	(revision 1671)
+++ trunk/SophyaExt/XAstroPack/objlist.list	(revision 1672)
@@ -1,1 +1,2 @@
+magfield.o
 xastropack.o
