source: branches/FACT++_lidctrl_usb/erfa/src/nut00b.c@ 19875

Last change on this file since 19875 was 18711, checked in by tbretz, 10 years ago
Updated to ERFA 1.3.0 (no relevant code change except the leap second at the beginning of 2017)
File size: 17.2 KB
Line 
1#include "erfa.h"
2
3void eraNut00b(double date1, double date2, double *dpsi, double *deps)
4/*
5** - - - - - - - - - -
6** e r a N u t 0 0 b
7** - - - - - - - - - -
8**
9** Nutation, IAU 2000B model.
10**
11** Given:
12** date1,date2 double TT as a 2-part Julian Date (Note 1)
13**
14** Returned:
15** dpsi,deps double nutation, luni-solar + planetary (Note 2)
16**
17** Notes:
18**
19** 1) The TT date date1+date2 is a Julian Date, apportioned in any
20** convenient way between the two arguments. For example,
21** JD(TT)=2450123.7 could be expressed in any of these ways,
22** among others:
23**
24** date1 date2
25**
26** 2450123.7 0.0 (JD method)
27** 2451545.0 -1421.3 (J2000 method)
28** 2400000.5 50123.2 (MJD method)
29** 2450123.5 0.2 (date & time method)
30**
31** The JD method is the most natural and convenient to use in
32** cases where the loss of several decimal digits of resolution
33** is acceptable. The J2000 method is best matched to the way
34** the argument is handled internally and will deliver the
35** optimum resolution. The MJD method and the date & time methods
36** are both good compromises between resolution and convenience.
37**
38** 2) The nutation components in longitude and obliquity are in radians
39** and with respect to the equinox and ecliptic of date. The
40** obliquity at J2000.0 is assumed to be the Lieske et al. (1977)
41** value of 84381.448 arcsec. (The errors that result from using
42** this function with the IAU 2006 value of 84381.406 arcsec can be
43** neglected.)
44**
45** The nutation model consists only of luni-solar terms, but
46** includes also a fixed offset which compensates for certain long-
47** period planetary terms (Note 7).
48**
49** 3) This function is an implementation of the IAU 2000B abridged
50** nutation model formally adopted by the IAU General Assembly in
51** 2000. The function computes the MHB_2000_SHORT luni-solar
52** nutation series (Luzum 2001), but without the associated
53** corrections for the precession rate adjustments and the offset
54** between the GCRS and J2000.0 mean poles.
55**
56** 4) The full IAU 2000A (MHB2000) nutation model contains nearly 1400
57** terms. The IAU 2000B model (McCarthy & Luzum 2003) contains only
58** 77 terms, plus additional simplifications, yet still delivers
59** results of 1 mas accuracy at present epochs. This combination of
60** accuracy and size makes the IAU 2000B abridged nutation model
61** suitable for most practical applications.
62**
63** The function delivers a pole accurate to 1 mas from 1900 to 2100
64** (usually better than 1 mas, very occasionally just outside
65** 1 mas). The full IAU 2000A model, which is implemented in the
66** function eraNut00a (q.v.), delivers considerably greater accuracy
67** at current dates; however, to realize this improved accuracy,
68** corrections for the essentially unpredictable free-core-nutation
69** (FCN) must also be included.
70**
71** 5) The present function provides classical nutation. The
72** MHB_2000_SHORT algorithm, from which it is adapted, deals also
73** with (i) the offsets between the GCRS and mean poles and (ii) the
74** adjustments in longitude and obliquity due to the changed
75** precession rates. These additional functions, namely frame bias
76** and precession adjustments, are supported by the ERFA functions
77** eraBi00 and eraPr00.
78**
79** 6) The MHB_2000_SHORT algorithm also provides "total" nutations,
80** comprising the arithmetic sum of the frame bias, precession
81** adjustments, and nutation (luni-solar + planetary). These total
82** nutations can be used in combination with an existing IAU 1976
83** precession implementation, such as eraPmat76, to deliver GCRS-
84** to-true predictions of mas accuracy at current epochs. However,
85** for symmetry with the eraNut00a function (q.v. for the reasons),
86** the ERFA functions do not generate the "total nutations"
87** directly. Should they be required, they could of course easily
88** be generated by calling eraBi00, eraPr00 and the present function
89** and adding the results.
90**
91** 7) The IAU 2000B model includes "planetary bias" terms that are
92** fixed in size but compensate for long-period nutations. The
93** amplitudes quoted in McCarthy & Luzum (2003), namely
94** Dpsi = -1.5835 mas and Depsilon = +1.6339 mas, are optimized for
95** the "total nutations" method described in Note 6. The Luzum
96** (2001) values used in this ERFA implementation, namely -0.135 mas
97** and +0.388 mas, are optimized for the "rigorous" method, where
98** frame bias, precession and nutation are applied separately and in
99** that order. During the interval 1995-2050, the ERFA
100** implementation delivers a maximum error of 1.001 mas (not
101** including FCN).
102**
103** References:
104**
105** Lieske, J.H., Lederle, T., Fricke, W., Morando, B., "Expressions
106** for the precession quantities based upon the IAU /1976/ system of
107** astronomical constants", Astron.Astrophys. 58, 1-2, 1-16. (1977)
108**
109** Luzum, B., private communication, 2001 (Fortran code
110** MHB_2000_SHORT)
111**
112** McCarthy, D.D. & Luzum, B.J., "An abridged model of the
113** precession-nutation of the celestial pole", Cel.Mech.Dyn.Astron.
114** 85, 37-49 (2003)
115**
116** Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M.,
117** Francou, G., Laskar, J., Astron.Astrophys. 282, 663-683 (1994)
118**
119** Copyright (C) 2013-2016, NumFOCUS Foundation.
120** Derived, with permission, from the SOFA library. See notes at end of file.
121*/
122{
123 double t, el, elp, f, d, om, arg, dp, de, sarg, carg,
124 dpsils, depsls, dpsipl, depspl;
125 int i;
126
127/* Units of 0.1 microarcsecond to radians */
128 static const double U2R = ERFA_DAS2R / 1e7;
129
130/* ---------------------------------------- */
131/* Fixed offsets in lieu of planetary terms */
132/* ---------------------------------------- */
133
134 static const double DPPLAN = -0.135 * ERFA_DMAS2R;
135 static const double DEPLAN = 0.388 * ERFA_DMAS2R;
136
137/* --------------------------------------------------- */
138/* Luni-solar nutation: argument and term coefficients */
139/* --------------------------------------------------- */
140
141/* The units for the sine and cosine coefficients are */
142/* 0.1 microarcsec and the same per Julian century */
143
144 static const struct {
145 int nl,nlp,nf,nd,nom; /* coefficients of l,l',F,D,Om */
146 double ps,pst,pc; /* longitude sin, t*sin, cos coefficients */
147 double ec,ect,es; /* obliquity cos, t*cos, sin coefficients */
148
149 } x[] = {
150
151 /* 1-10 */
152 { 0, 0, 0, 0,1,
153 -172064161.0, -174666.0, 33386.0, 92052331.0, 9086.0, 15377.0},
154 { 0, 0, 2,-2,2,
155 -13170906.0, -1675.0, -13696.0, 5730336.0, -3015.0, -4587.0},
156 { 0, 0, 2, 0,2,-2276413.0,-234.0, 2796.0, 978459.0,-485.0,1374.0},
157 { 0, 0, 0, 0,2,2074554.0, 207.0, -698.0,-897492.0, 470.0,-291.0},
158 { 0, 1, 0, 0,0,1475877.0,-3633.0,11817.0, 73871.0,-184.0,-1924.0},
159 { 0, 1, 2,-2,2,-516821.0, 1226.0, -524.0, 224386.0,-677.0,-174.0},
160 { 1, 0, 0, 0,0, 711159.0, 73.0, -872.0, -6750.0, 0.0, 358.0},
161 { 0, 0, 2, 0,1,-387298.0, -367.0, 380.0, 200728.0, 18.0, 318.0},
162 { 1, 0, 2, 0,2,-301461.0, -36.0, 816.0, 129025.0, -63.0, 367.0},
163 { 0,-1, 2,-2,2, 215829.0, -494.0, 111.0, -95929.0, 299.0, 132.0},
164
165 /* 11-20 */
166 { 0, 0, 2,-2,1, 128227.0, 137.0, 181.0, -68982.0, -9.0, 39.0},
167 {-1, 0, 2, 0,2, 123457.0, 11.0, 19.0, -53311.0, 32.0, -4.0},
168 {-1, 0, 0, 2,0, 156994.0, 10.0, -168.0, -1235.0, 0.0, 82.0},
169 { 1, 0, 0, 0,1, 63110.0, 63.0, 27.0, -33228.0, 0.0, -9.0},
170 {-1, 0, 0, 0,1, -57976.0, -63.0, -189.0, 31429.0, 0.0, -75.0},
171 {-1, 0, 2, 2,2, -59641.0, -11.0, 149.0, 25543.0, -11.0, 66.0},
172 { 1, 0, 2, 0,1, -51613.0, -42.0, 129.0, 26366.0, 0.0, 78.0},
173 {-2, 0, 2, 0,1, 45893.0, 50.0, 31.0, -24236.0, -10.0, 20.0},
174 { 0, 0, 0, 2,0, 63384.0, 11.0, -150.0, -1220.0, 0.0, 29.0},
175 { 0, 0, 2, 2,2, -38571.0, -1.0, 158.0, 16452.0, -11.0, 68.0},
176
177 /* 21-30 */
178 { 0,-2, 2,-2,2, 32481.0, 0.0, 0.0, -13870.0, 0.0, 0.0},
179 {-2, 0, 0, 2,0, -47722.0, 0.0, -18.0, 477.0, 0.0, -25.0},
180 { 2, 0, 2, 0,2, -31046.0, -1.0, 131.0, 13238.0, -11.0, 59.0},
181 { 1, 0, 2,-2,2, 28593.0, 0.0, -1.0, -12338.0, 10.0, -3.0},
182 {-1, 0, 2, 0,1, 20441.0, 21.0, 10.0, -10758.0, 0.0, -3.0},
183 { 2, 0, 0, 0,0, 29243.0, 0.0, -74.0, -609.0, 0.0, 13.0},
184 { 0, 0, 2, 0,0, 25887.0, 0.0, -66.0, -550.0, 0.0, 11.0},
185 { 0, 1, 0, 0,1, -14053.0, -25.0, 79.0, 8551.0, -2.0, -45.0},
186 {-1, 0, 0, 2,1, 15164.0, 10.0, 11.0, -8001.0, 0.0, -1.0},
187 { 0, 2, 2,-2,2, -15794.0, 72.0, -16.0, 6850.0, -42.0, -5.0},
188
189 /* 31-40 */
190 { 0, 0,-2, 2,0, 21783.0, 0.0, 13.0, -167.0, 0.0, 13.0},
191 { 1, 0, 0,-2,1, -12873.0, -10.0, -37.0, 6953.0, 0.0, -14.0},
192 { 0,-1, 0, 0,1, -12654.0, 11.0, 63.0, 6415.0, 0.0, 26.0},
193 {-1, 0, 2, 2,1, -10204.0, 0.0, 25.0, 5222.0, 0.0, 15.0},
194 { 0, 2, 0, 0,0, 16707.0, -85.0, -10.0, 168.0, -1.0, 10.0},
195 { 1, 0, 2, 2,2, -7691.0, 0.0, 44.0, 3268.0, 0.0, 19.0},
196 {-2, 0, 2, 0,0, -11024.0, 0.0, -14.0, 104.0, 0.0, 2.0},
197 { 0, 1, 2, 0,2, 7566.0, -21.0, -11.0, -3250.0, 0.0, -5.0},
198 { 0, 0, 2, 2,1, -6637.0, -11.0, 25.0, 3353.0, 0.0, 14.0},
199 { 0,-1, 2, 0,2, -7141.0, 21.0, 8.0, 3070.0, 0.0, 4.0},
200
201 /* 41-50 */
202 { 0, 0, 0, 2,1, -6302.0, -11.0, 2.0, 3272.0, 0.0, 4.0},
203 { 1, 0, 2,-2,1, 5800.0, 10.0, 2.0, -3045.0, 0.0, -1.0},
204 { 2, 0, 2,-2,2, 6443.0, 0.0, -7.0, -2768.0, 0.0, -4.0},
205 {-2, 0, 0, 2,1, -5774.0, -11.0, -15.0, 3041.0, 0.0, -5.0},
206 { 2, 0, 2, 0,1, -5350.0, 0.0, 21.0, 2695.0, 0.0, 12.0},
207 { 0,-1, 2,-2,1, -4752.0, -11.0, -3.0, 2719.0, 0.0, -3.0},
208 { 0, 0, 0,-2,1, -4940.0, -11.0, -21.0, 2720.0, 0.0, -9.0},
209 {-1,-1, 0, 2,0, 7350.0, 0.0, -8.0, -51.0, 0.0, 4.0},
210 { 2, 0, 0,-2,1, 4065.0, 0.0, 6.0, -2206.0, 0.0, 1.0},
211 { 1, 0, 0, 2,0, 6579.0, 0.0, -24.0, -199.0, 0.0, 2.0},
212
213 /* 51-60 */
214 { 0, 1, 2,-2,1, 3579.0, 0.0, 5.0, -1900.0, 0.0, 1.0},
215 { 1,-1, 0, 0,0, 4725.0, 0.0, -6.0, -41.0, 0.0, 3.0},
216 {-2, 0, 2, 0,2, -3075.0, 0.0, -2.0, 1313.0, 0.0, -1.0},
217 { 3, 0, 2, 0,2, -2904.0, 0.0, 15.0, 1233.0, 0.0, 7.0},
218 { 0,-1, 0, 2,0, 4348.0, 0.0, -10.0, -81.0, 0.0, 2.0},
219 { 1,-1, 2, 0,2, -2878.0, 0.0, 8.0, 1232.0, 0.0, 4.0},
220 { 0, 0, 0, 1,0, -4230.0, 0.0, 5.0, -20.0, 0.0, -2.0},
221 {-1,-1, 2, 2,2, -2819.0, 0.0, 7.0, 1207.0, 0.0, 3.0},
222 {-1, 0, 2, 0,0, -4056.0, 0.0, 5.0, 40.0, 0.0, -2.0},
223 { 0,-1, 2, 2,2, -2647.0, 0.0, 11.0, 1129.0, 0.0, 5.0},
224
225 /* 61-70 */
226 {-2, 0, 0, 0,1, -2294.0, 0.0, -10.0, 1266.0, 0.0, -4.0},
227 { 1, 1, 2, 0,2, 2481.0, 0.0, -7.0, -1062.0, 0.0, -3.0},
228 { 2, 0, 0, 0,1, 2179.0, 0.0, -2.0, -1129.0, 0.0, -2.0},
229 {-1, 1, 0, 1,0, 3276.0, 0.0, 1.0, -9.0, 0.0, 0.0},
230 { 1, 1, 0, 0,0, -3389.0, 0.0, 5.0, 35.0, 0.0, -2.0},
231 { 1, 0, 2, 0,0, 3339.0, 0.0, -13.0, -107.0, 0.0, 1.0},
232 {-1, 0, 2,-2,1, -1987.0, 0.0, -6.0, 1073.0, 0.0, -2.0},
233 { 1, 0, 0, 0,2, -1981.0, 0.0, 0.0, 854.0, 0.0, 0.0},
234 {-1, 0, 0, 1,0, 4026.0, 0.0, -353.0, -553.0, 0.0,-139.0},
235 { 0, 0, 2, 1,2, 1660.0, 0.0, -5.0, -710.0, 0.0, -2.0},
236
237 /* 71-77 */
238 {-1, 0, 2, 4,2, -1521.0, 0.0, 9.0, 647.0, 0.0, 4.0},
239 {-1, 1, 0, 1,1, 1314.0, 0.0, 0.0, -700.0, 0.0, 0.0},
240 { 0,-2, 2,-2,1, -1283.0, 0.0, 0.0, 672.0, 0.0, 0.0},
241 { 1, 0, 2, 2,1, -1331.0, 0.0, 8.0, 663.0, 0.0, 4.0},
242 {-2, 0, 2, 2,2, 1383.0, 0.0, -2.0, -594.0, 0.0, -2.0},
243 {-1, 0, 0, 0,2, 1405.0, 0.0, 4.0, -610.0, 0.0, 2.0},
244 { 1, 1, 2,-2,2, 1290.0, 0.0, 0.0, -556.0, 0.0, 0.0}
245 };
246
247/* Number of terms in the series */
248 const int NLS = (int) (sizeof x / sizeof x[0]);
249
250/*--------------------------------------------------------------------*/
251
252/* Interval between fundamental epoch J2000.0 and given date (JC). */
253 t = ((date1 - ERFA_DJ00) + date2) / ERFA_DJC;
254
255/* --------------------*/
256/* LUNI-SOLAR NUTATION */
257/* --------------------*/
258
259/* Fundamental (Delaunay) arguments from Simon et al. (1994) */
260
261/* Mean anomaly of the Moon. */
262 el = fmod(485868.249036 + (1717915923.2178) * t, ERFA_TURNAS) * ERFA_DAS2R;
263
264/* Mean anomaly of the Sun. */
265 elp = fmod(1287104.79305 + (129596581.0481) * t, ERFA_TURNAS) * ERFA_DAS2R;
266
267/* Mean argument of the latitude of the Moon. */
268 f = fmod(335779.526232 + (1739527262.8478) * t, ERFA_TURNAS) * ERFA_DAS2R;
269
270/* Mean elongation of the Moon from the Sun. */
271 d = fmod(1072260.70369 + (1602961601.2090) * t, ERFA_TURNAS) * ERFA_DAS2R;
272
273/* Mean longitude of the ascending node of the Moon. */
274 om = fmod(450160.398036 + (-6962890.5431) * t, ERFA_TURNAS) * ERFA_DAS2R;
275
276/* Initialize the nutation values. */
277 dp = 0.0;
278 de = 0.0;
279
280/* Summation of luni-solar nutation series (smallest terms first). */
281 for (i = NLS-1; i >= 0; i--) {
282
283 /* Argument and functions. */
284 arg = fmod( (double)x[i].nl * el +
285 (double)x[i].nlp * elp +
286 (double)x[i].nf * f +
287 (double)x[i].nd * d +
288 (double)x[i].nom * om, ERFA_D2PI );
289 sarg = sin(arg);
290 carg = cos(arg);
291
292 /* Term. */
293 dp += (x[i].ps + x[i].pst * t) * sarg + x[i].pc * carg;
294 de += (x[i].ec + x[i].ect * t) * carg + x[i].es * sarg;
295 }
296
297/* Convert from 0.1 microarcsec units to radians. */
298 dpsils = dp * U2R;
299 depsls = de * U2R;
300
301/* ------------------------------*/
302/* IN LIEU OF PLANETARY NUTATION */
303/* ------------------------------*/
304
305/* Fixed offset to correct for missing terms in truncated series. */
306 dpsipl = DPPLAN;
307 depspl = DEPLAN;
308
309/* --------*/
310/* RESULTS */
311/* --------*/
312
313/* Add luni-solar and planetary components. */
314 *dpsi = dpsils + dpsipl;
315 *deps = depsls + depspl;
316
317 return;
318
319}
320/*----------------------------------------------------------------------
321**
322**
323** Copyright (C) 2013-2016, NumFOCUS Foundation.
324** All rights reserved.
325**
326** This library is derived, with permission, from the International
327** Astronomical Union's "Standards of Fundamental Astronomy" library,
328** available from http://www.iausofa.org.
329**
330** The ERFA version is intended to retain identical functionality to
331** the SOFA library, but made distinct through different function and
332** file names, as set out in the SOFA license conditions. The SOFA
333** original has a role as a reference standard for the IAU and IERS,
334** and consequently redistribution is permitted only in its unaltered
335** state. The ERFA version is not subject to this restriction and
336** therefore can be included in distributions which do not support the
337** concept of "read only" software.
338**
339** Although the intent is to replicate the SOFA API (other than
340** replacement of prefix names) and results (with the exception of
341** bugs; any that are discovered will be fixed), SOFA is not
342** responsible for any errors found in this version of the library.
343**
344** If you wish to acknowledge the SOFA heritage, please acknowledge
345** that you are using a library derived from SOFA, rather than SOFA
346** itself.
347**
348**
349** TERMS AND CONDITIONS
350**
351** Redistribution and use in source and binary forms, with or without
352** modification, are permitted provided that the following conditions
353** are met:
354**
355** 1 Redistributions of source code must retain the above copyright
356** notice, this list of conditions and the following disclaimer.
357**
358** 2 Redistributions in binary form must reproduce the above copyright
359** notice, this list of conditions and the following disclaimer in
360** the documentation and/or other materials provided with the
361** distribution.
362**
363** 3 Neither the name of the Standards Of Fundamental Astronomy Board,
364** the International Astronomical Union nor the names of its
365** contributors may be used to endorse or promote products derived
366** from this software without specific prior written permission.
367**
368** THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
369** "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
370** LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
371** FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
372** COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
373** INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
374** BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
375** LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
376** CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
377** LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
378** ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
379** POSSIBILITY OF SUCH DAMAGE.
380**
381*/
Note: See TracBrowser for help on using the repository browser.