00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #include <GeographicLib/MagneticCircle.hpp>
00011 #include <fstream>
00012 #include <sstream>
00013 #include <GeographicLib/Geocentric.hpp>
00014
00015 namespace GeographicLib {
00016
00017 using namespace std;
00018
00019 void MagneticCircle::Field(real lon, bool diffp,
00020 real& Bx, real& By, real& Bz,
00021 real& Bxt, real& Byt, real& Bzt) const {
00022 real clam, slam;
00023 CircularEngine::cossin(lon, clam, slam);
00024 real M[Geocentric::dim2_];
00025 Geocentric::Rotation(_sphi, _cphi, slam, clam, M);
00026 real BX0, BY0, BZ0, BX1, BY1, BZ1;
00027 _circ0(clam, slam, BX0, BY0, BZ0);
00028 _circ1(clam, slam, BX1, BY1, BZ1);
00029 if (_interpolate) {
00030 BX1 = (BX1 - BX0) / _dt0;
00031 BY1 = (BY1 - BY0) / _dt0;
00032 BZ1 = (BZ1 - BZ0) / _dt0;
00033 }
00034 BX0 += _t1 * BX1;
00035 BY0 += _t1 * BY1;
00036 BZ0 += _t1 * BZ1;
00037 if (diffp) {
00038 Geocentric::Unrotate(M, BX1, BY1, BZ1, Bxt, Byt, Bzt);
00039 Bxt *= - _a;
00040 Byt *= - _a;
00041 Bzt *= - _a;
00042 }
00043 Geocentric::Unrotate(M, BX0, BY0, BZ0, Bx, By, Bz);
00044 Bx *= - _a;
00045 By *= - _a;
00046 Bz *= - _a;
00047 }
00048
00049 }