source: trunk/FACT++/sofa/src/nut80.c@ 18346

Last change on this file since 18346 was 18346, checked in by tbretz, 9 years ago
File size: 15.9 KB
Line 
1#include "sofa.h"
2
3void iauNut80(double date1, double date2, double *dpsi, double *deps)
4/*
5** - - - - - - - - -
6** i a u N u t 8 0
7** - - - - - - - - -
8**
9** Nutation, IAU 1980 model.
10**
11** This function is part of the International Astronomical Union's
12** SOFA (Standards Of Fundamental Astronomy) software collection.
13**
14** Status: canonical model.
15**
16** Given:
17** date1,date2 double TT as a 2-part Julian Date (Note 1)
18**
19** Returned:
20** dpsi double nutation in longitude (radians)
21** deps double nutation in obliquity (radians)
22**
23** Notes:
24**
25** 1) The TT date date1+date2 is a Julian Date, apportioned in any
26** convenient way between the two arguments. For example,
27** JD(TT)=2450123.7 could be expressed in any of these ways,
28** among others:
29**
30** date1 date2
31**
32** 2450123.7 0.0 (JD method)
33** 2451545.0 -1421.3 (J2000 method)
34** 2400000.5 50123.2 (MJD method)
35** 2450123.5 0.2 (date & time method)
36**
37** The JD method is the most natural and convenient to use in
38** cases where the loss of several decimal digits of resolution
39** is acceptable. The J2000 method is best matched to the way
40** the argument is handled internally and will deliver the
41** optimum resolution. The MJD method and the date & time methods
42** are both good compromises between resolution and convenience.
43**
44** 2) The nutation components are with respect to the ecliptic of
45** date.
46**
47** Called:
48** iauAnpm normalize angle into range +/- pi
49**
50** Reference:
51**
52** Explanatory Supplement to the Astronomical Almanac,
53** P. Kenneth Seidelmann (ed), University Science Books (1992),
54** Section 3.222 (p111).
55**
56** This revision: 2013 June 18
57**
58** SOFA release 2015-02-09
59**
60** Copyright (C) 2015 IAU SOFA Board. See notes at end.
61*/
62{
63 double t, el, elp, f, d, om, dp, de, arg, s, c;
64 int j;
65
66/* Units of 0.1 milliarcsecond to radians */
67 const double U2R = DAS2R / 1e4;
68
69/* ------------------------------------------------ */
70/* Table of multiples of arguments and coefficients */
71/* ------------------------------------------------ */
72
73/* The units for the sine and cosine coefficients are 0.1 mas and */
74/* the same per Julian century */
75
76 static const struct {
77 int nl,nlp,nf,nd,nom; /* coefficients of l,l',F,D,Om */
78 double sp,spt; /* longitude sine, 1 and t coefficients */
79 double ce,cet; /* obliquity cosine, 1 and t coefficients */
80 } x[] = {
81
82 /* 1-10 */
83 { 0, 0, 0, 0, 1, -171996.0, -174.2, 92025.0, 8.9 },
84 { 0, 0, 0, 0, 2, 2062.0, 0.2, -895.0, 0.5 },
85 { -2, 0, 2, 0, 1, 46.0, 0.0, -24.0, 0.0 },
86 { 2, 0, -2, 0, 0, 11.0, 0.0, 0.0, 0.0 },
87 { -2, 0, 2, 0, 2, -3.0, 0.0, 1.0, 0.0 },
88 { 1, -1, 0, -1, 0, -3.0, 0.0, 0.0, 0.0 },
89 { 0, -2, 2, -2, 1, -2.0, 0.0, 1.0, 0.0 },
90 { 2, 0, -2, 0, 1, 1.0, 0.0, 0.0, 0.0 },
91 { 0, 0, 2, -2, 2, -13187.0, -1.6, 5736.0, -3.1 },
92 { 0, 1, 0, 0, 0, 1426.0, -3.4, 54.0, -0.1 },
93
94 /* 11-20 */
95 { 0, 1, 2, -2, 2, -517.0, 1.2, 224.0, -0.6 },
96 { 0, -1, 2, -2, 2, 217.0, -0.5, -95.0, 0.3 },
97 { 0, 0, 2, -2, 1, 129.0, 0.1, -70.0, 0.0 },
98 { 2, 0, 0, -2, 0, 48.0, 0.0, 1.0, 0.0 },
99 { 0, 0, 2, -2, 0, -22.0, 0.0, 0.0, 0.0 },
100 { 0, 2, 0, 0, 0, 17.0, -0.1, 0.0, 0.0 },
101 { 0, 1, 0, 0, 1, -15.0, 0.0, 9.0, 0.0 },
102 { 0, 2, 2, -2, 2, -16.0, 0.1, 7.0, 0.0 },
103 { 0, -1, 0, 0, 1, -12.0, 0.0, 6.0, 0.0 },
104 { -2, 0, 0, 2, 1, -6.0, 0.0, 3.0, 0.0 },
105
106 /* 21-30 */
107 { 0, -1, 2, -2, 1, -5.0, 0.0, 3.0, 0.0 },
108 { 2, 0, 0, -2, 1, 4.0, 0.0, -2.0, 0.0 },
109 { 0, 1, 2, -2, 1, 4.0, 0.0, -2.0, 0.0 },
110 { 1, 0, 0, -1, 0, -4.0, 0.0, 0.0, 0.0 },
111 { 2, 1, 0, -2, 0, 1.0, 0.0, 0.0, 0.0 },
112 { 0, 0, -2, 2, 1, 1.0, 0.0, 0.0, 0.0 },
113 { 0, 1, -2, 2, 0, -1.0, 0.0, 0.0, 0.0 },
114 { 0, 1, 0, 0, 2, 1.0, 0.0, 0.0, 0.0 },
115 { -1, 0, 0, 1, 1, 1.0, 0.0, 0.0, 0.0 },
116 { 0, 1, 2, -2, 0, -1.0, 0.0, 0.0, 0.0 },
117
118 /* 31-40 */
119 { 0, 0, 2, 0, 2, -2274.0, -0.2, 977.0, -0.5 },
120 { 1, 0, 0, 0, 0, 712.0, 0.1, -7.0, 0.0 },
121 { 0, 0, 2, 0, 1, -386.0, -0.4, 200.0, 0.0 },
122 { 1, 0, 2, 0, 2, -301.0, 0.0, 129.0, -0.1 },
123 { 1, 0, 0, -2, 0, -158.0, 0.0, -1.0, 0.0 },
124 { -1, 0, 2, 0, 2, 123.0, 0.0, -53.0, 0.0 },
125 { 0, 0, 0, 2, 0, 63.0, 0.0, -2.0, 0.0 },
126 { 1, 0, 0, 0, 1, 63.0, 0.1, -33.0, 0.0 },
127 { -1, 0, 0, 0, 1, -58.0, -0.1, 32.0, 0.0 },
128 { -1, 0, 2, 2, 2, -59.0, 0.0, 26.0, 0.0 },
129
130 /* 41-50 */
131 { 1, 0, 2, 0, 1, -51.0, 0.0, 27.0, 0.0 },
132 { 0, 0, 2, 2, 2, -38.0, 0.0, 16.0, 0.0 },
133 { 2, 0, 0, 0, 0, 29.0, 0.0, -1.0, 0.0 },
134 { 1, 0, 2, -2, 2, 29.0, 0.0, -12.0, 0.0 },
135 { 2, 0, 2, 0, 2, -31.0, 0.0, 13.0, 0.0 },
136 { 0, 0, 2, 0, 0, 26.0, 0.0, -1.0, 0.0 },
137 { -1, 0, 2, 0, 1, 21.0, 0.0, -10.0, 0.0 },
138 { -1, 0, 0, 2, 1, 16.0, 0.0, -8.0, 0.0 },
139 { 1, 0, 0, -2, 1, -13.0, 0.0, 7.0, 0.0 },
140 { -1, 0, 2, 2, 1, -10.0, 0.0, 5.0, 0.0 },
141
142 /* 51-60 */
143 { 1, 1, 0, -2, 0, -7.0, 0.0, 0.0, 0.0 },
144 { 0, 1, 2, 0, 2, 7.0, 0.0, -3.0, 0.0 },
145 { 0, -1, 2, 0, 2, -7.0, 0.0, 3.0, 0.0 },
146 { 1, 0, 2, 2, 2, -8.0, 0.0, 3.0, 0.0 },
147 { 1, 0, 0, 2, 0, 6.0, 0.0, 0.0, 0.0 },
148 { 2, 0, 2, -2, 2, 6.0, 0.0, -3.0, 0.0 },
149 { 0, 0, 0, 2, 1, -6.0, 0.0, 3.0, 0.0 },
150 { 0, 0, 2, 2, 1, -7.0, 0.0, 3.0, 0.0 },
151 { 1, 0, 2, -2, 1, 6.0, 0.0, -3.0, 0.0 },
152 { 0, 0, 0, -2, 1, -5.0, 0.0, 3.0, 0.0 },
153
154 /* 61-70 */
155 { 1, -1, 0, 0, 0, 5.0, 0.0, 0.0, 0.0 },
156 { 2, 0, 2, 0, 1, -5.0, 0.0, 3.0, 0.0 },
157 { 0, 1, 0, -2, 0, -4.0, 0.0, 0.0, 0.0 },
158 { 1, 0, -2, 0, 0, 4.0, 0.0, 0.0, 0.0 },
159 { 0, 0, 0, 1, 0, -4.0, 0.0, 0.0, 0.0 },
160 { 1, 1, 0, 0, 0, -3.0, 0.0, 0.0, 0.0 },
161 { 1, 0, 2, 0, 0, 3.0, 0.0, 0.0, 0.0 },
162 { 1, -1, 2, 0, 2, -3.0, 0.0, 1.0, 0.0 },
163 { -1, -1, 2, 2, 2, -3.0, 0.0, 1.0, 0.0 },
164 { -2, 0, 0, 0, 1, -2.0, 0.0, 1.0, 0.0 },
165
166 /* 71-80 */
167 { 3, 0, 2, 0, 2, -3.0, 0.0, 1.0, 0.0 },
168 { 0, -1, 2, 2, 2, -3.0, 0.0, 1.0, 0.0 },
169 { 1, 1, 2, 0, 2, 2.0, 0.0, -1.0, 0.0 },
170 { -1, 0, 2, -2, 1, -2.0, 0.0, 1.0, 0.0 },
171 { 2, 0, 0, 0, 1, 2.0, 0.0, -1.0, 0.0 },
172 { 1, 0, 0, 0, 2, -2.0, 0.0, 1.0, 0.0 },
173 { 3, 0, 0, 0, 0, 2.0, 0.0, 0.0, 0.0 },
174 { 0, 0, 2, 1, 2, 2.0, 0.0, -1.0, 0.0 },
175 { -1, 0, 0, 0, 2, 1.0, 0.0, -1.0, 0.0 },
176 { 1, 0, 0, -4, 0, -1.0, 0.0, 0.0, 0.0 },
177
178 /* 81-90 */
179 { -2, 0, 2, 2, 2, 1.0, 0.0, -1.0, 0.0 },
180 { -1, 0, 2, 4, 2, -2.0, 0.0, 1.0, 0.0 },
181 { 2, 0, 0, -4, 0, -1.0, 0.0, 0.0, 0.0 },
182 { 1, 1, 2, -2, 2, 1.0, 0.0, -1.0, 0.0 },
183 { 1, 0, 2, 2, 1, -1.0, 0.0, 1.0, 0.0 },
184 { -2, 0, 2, 4, 2, -1.0, 0.0, 1.0, 0.0 },
185 { -1, 0, 4, 0, 2, 1.0, 0.0, 0.0, 0.0 },
186 { 1, -1, 0, -2, 0, 1.0, 0.0, 0.0, 0.0 },
187 { 2, 0, 2, -2, 1, 1.0, 0.0, -1.0, 0.0 },
188 { 2, 0, 2, 2, 2, -1.0, 0.0, 0.0, 0.0 },
189
190 /* 91-100 */
191 { 1, 0, 0, 2, 1, -1.0, 0.0, 0.0, 0.0 },
192 { 0, 0, 4, -2, 2, 1.0, 0.0, 0.0, 0.0 },
193 { 3, 0, 2, -2, 2, 1.0, 0.0, 0.0, 0.0 },
194 { 1, 0, 2, -2, 0, -1.0, 0.0, 0.0, 0.0 },
195 { 0, 1, 2, 0, 1, 1.0, 0.0, 0.0, 0.0 },
196 { -1, -1, 0, 2, 1, 1.0, 0.0, 0.0, 0.0 },
197 { 0, 0, -2, 0, 1, -1.0, 0.0, 0.0, 0.0 },
198 { 0, 0, 2, -1, 2, -1.0, 0.0, 0.0, 0.0 },
199 { 0, 1, 0, 2, 0, -1.0, 0.0, 0.0, 0.0 },
200 { 1, 0, -2, -2, 0, -1.0, 0.0, 0.0, 0.0 },
201
202 /* 101-106 */
203 { 0, -1, 2, 0, 1, -1.0, 0.0, 0.0, 0.0 },
204 { 1, 1, 0, -2, 1, -1.0, 0.0, 0.0, 0.0 },
205 { 1, 0, -2, 2, 0, -1.0, 0.0, 0.0, 0.0 },
206 { 2, 0, 0, 2, 0, 1.0, 0.0, 0.0, 0.0 },
207 { 0, 0, 2, 4, 2, -1.0, 0.0, 0.0, 0.0 },
208 { 0, 1, 0, 1, 0, 1.0, 0.0, 0.0, 0.0 }
209 };
210
211/* Number of terms in the series */
212 const int NT = (int) (sizeof x / sizeof x[0]);
213
214/*--------------------------------------------------------------------*/
215
216/* Interval between fundamental epoch J2000.0 and given date (JC). */
217 t = ((date1 - DJ00) + date2) / DJC;
218
219/* --------------------- */
220/* Fundamental arguments */
221/* --------------------- */
222
223/* Mean longitude of Moon minus mean longitude of Moon's perigee. */
224 el = iauAnpm(
225 (485866.733 + (715922.633 + (31.310 + 0.064 * t) * t) * t)
226 * DAS2R + fmod(1325.0 * t, 1.0) * D2PI);
227
228/* Mean longitude of Sun minus mean longitude of Sun's perigee. */
229 elp = iauAnpm(
230 (1287099.804 + (1292581.224 + (-0.577 - 0.012 * t) * t) * t)
231 * DAS2R + fmod(99.0 * t, 1.0) * D2PI);
232
233/* Mean longitude of Moon minus mean longitude of Moon's node. */
234 f = iauAnpm(
235 (335778.877 + (295263.137 + (-13.257 + 0.011 * t) * t) * t)
236 * DAS2R + fmod(1342.0 * t, 1.0) * D2PI);
237
238/* Mean elongation of Moon from Sun. */
239 d = iauAnpm(
240 (1072261.307 + (1105601.328 + (-6.891 + 0.019 * t) * t) * t)
241 * DAS2R + fmod(1236.0 * t, 1.0) * D2PI);
242
243/* Longitude of the mean ascending node of the lunar orbit on the */
244/* ecliptic, measured from the mean equinox of date. */
245 om = iauAnpm(
246 (450160.280 + (-482890.539 + (7.455 + 0.008 * t) * t) * t)
247 * DAS2R + fmod(-5.0 * t, 1.0) * D2PI);
248
249/* --------------- */
250/* Nutation series */
251/* --------------- */
252
253/* Initialize nutation components. */
254 dp = 0.0;
255 de = 0.0;
256
257/* Sum the nutation terms, ending with the biggest. */
258 for (j = NT-1; j >= 0; j--) {
259
260 /* Form argument for current term. */
261 arg = (double)x[j].nl * el
262 + (double)x[j].nlp * elp
263 + (double)x[j].nf * f
264 + (double)x[j].nd * d
265 + (double)x[j].nom * om;
266
267 /* Accumulate current nutation term. */
268 s = x[j].sp + x[j].spt * t;
269 c = x[j].ce + x[j].cet * t;
270 if (s != 0.0) dp += s * sin(arg);
271 if (c != 0.0) de += c * cos(arg);
272 }
273
274/* Convert results from 0.1 mas units to radians. */
275 *dpsi = dp * U2R;
276 *deps = de * U2R;
277
278 return;
279
280/*----------------------------------------------------------------------
281**
282** Copyright (C) 2015
283** Standards Of Fundamental Astronomy Board
284** of the International Astronomical Union.
285**
286** =====================
287** SOFA Software License
288** =====================
289**
290** NOTICE TO USER:
291**
292** BY USING THIS SOFTWARE YOU ACCEPT THE FOLLOWING SIX TERMS AND
293** CONDITIONS WHICH APPLY TO ITS USE.
294**
295** 1. The Software is owned by the IAU SOFA Board ("SOFA").
296**
297** 2. Permission is granted to anyone to use the SOFA software for any
298** purpose, including commercial applications, free of charge and
299** without payment of royalties, subject to the conditions and
300** restrictions listed below.
301**
302** 3. You (the user) may copy and distribute SOFA source code to others,
303** and use and adapt its code and algorithms in your own software,
304** on a world-wide, royalty-free basis. That portion of your
305** distribution that does not consist of intact and unchanged copies
306** of SOFA source code files is a "derived work" that must comply
307** with the following requirements:
308**
309** a) Your work shall be marked or carry a statement that it
310** (i) uses routines and computations derived by you from
311** software provided by SOFA under license to you; and
312** (ii) does not itself constitute software provided by and/or
313** endorsed by SOFA.
314**
315** b) The source code of your derived work must contain descriptions
316** of how the derived work is based upon, contains and/or differs
317** from the original SOFA software.
318**
319** c) The names of all routines in your derived work shall not
320** include the prefix "iau" or "sofa" or trivial modifications
321** thereof such as changes of case.
322**
323** d) The origin of the SOFA components of your derived work must
324** not be misrepresented; you must not claim that you wrote the
325** original software, nor file a patent application for SOFA
326** software or algorithms embedded in the SOFA software.
327**
328** e) These requirements must be reproduced intact in any source
329** distribution and shall apply to anyone to whom you have
330** granted a further right to modify the source code of your
331** derived work.
332**
333** Note that, as originally distributed, the SOFA software is
334** intended to be a definitive implementation of the IAU standards,
335** and consequently third-party modifications are discouraged. All
336** variations, no matter how minor, must be explicitly marked as
337** such, as explained above.
338**
339** 4. You shall not cause the SOFA software to be brought into
340** disrepute, either by misuse, or use for inappropriate tasks, or
341** by inappropriate modification.
342**
343** 5. The SOFA software is provided "as is" and SOFA makes no warranty
344** as to its use or performance. SOFA does not and cannot warrant
345** the performance or results which the user may obtain by using the
346** SOFA software. SOFA makes no warranties, express or implied, as
347** to non-infringement of third party rights, merchantability, or
348** fitness for any particular purpose. In no event will SOFA be
349** liable to the user for any consequential, incidental, or special
350** damages, including any lost profits or lost savings, even if a
351** SOFA representative has been advised of such damages, or for any
352** claim by any third party.
353**
354** 6. The provision of any version of the SOFA software under the terms
355** and conditions specified herein does not imply that future
356** versions will also be made available under the same terms and
357** conditions.
358*
359** In any published work or commercial product which uses the SOFA
360** software directly, acknowledgement (see www.iausofa.org) is
361** appreciated.
362**
363** Correspondence concerning SOFA software should be addressed as
364** follows:
365**
366** By email: sofa@ukho.gov.uk
367** By post: IAU SOFA Center
368** HM Nautical Almanac Office
369** UK Hydrographic Office
370** Admiralty Way, Taunton
371** Somerset, TA1 2DN
372** United Kingdom
373**
374**--------------------------------------------------------------------*/
375}
Note: See TracBrowser for help on using the repository browser.