| 1 | #include "slalib.h"
|
|---|
| 2 | #include "slamac.h"
|
|---|
| 3 | void slaPrebn ( double bep0, double bep1, double rmatp[3][3] )
|
|---|
| 4 | /*
|
|---|
| 5 | ** - - - - - - - - -
|
|---|
| 6 | ** s l a P r e b n
|
|---|
| 7 | ** - - - - - - - - -
|
|---|
| 8 | **
|
|---|
| 9 | ** Generate the matrix of precession between two epochs,
|
|---|
| 10 | ** using the old, pre-IAU1976, Bessel-Newcomb model, using
|
|---|
| 11 | ** Kinoshita's formulation (double precision)
|
|---|
| 12 | **
|
|---|
| 13 | ** Given:
|
|---|
| 14 | ** BEP0 double beginning Besselian epoch
|
|---|
| 15 | ** BEP1 double ending Besselian epoch
|
|---|
| 16 | **
|
|---|
| 17 | ** Returned:
|
|---|
| 18 | ** RMATP double[3][3] precession matrix
|
|---|
| 19 | **
|
|---|
| 20 | ** The matrix is in the sense v(bep1) = rmatp * v(bep0)
|
|---|
| 21 | **
|
|---|
| 22 | ** Reference:
|
|---|
| 23 | ** Kinoshita, H. (1975) 'Formulas for precession', SAO Special
|
|---|
| 24 | ** Report No. 364, Smithsonian Institution Astrophysical
|
|---|
| 25 | ** Observatory, Cambridge, Massachusetts.
|
|---|
| 26 | **
|
|---|
| 27 | ** Called: slaDeuler
|
|---|
| 28 | **
|
|---|
| 29 | ** Defined in slamac.h: DAS2R
|
|---|
| 30 | **
|
|---|
| 31 | ** Last revision: 30 October 1993
|
|---|
| 32 | **
|
|---|
| 33 | ** Copyright P.T.Wallace. All rights reserved.
|
|---|
| 34 | */
|
|---|
| 35 | {
|
|---|
| 36 | double bigt, t, tas2r, w, zeta, z, theta;
|
|---|
| 37 |
|
|---|
| 38 | /* Interval between basic epoch B1850.0 and beginning epoch in TC */
|
|---|
| 39 | bigt = ( bep0 - 1850.0 ) / 100.0;
|
|---|
| 40 |
|
|---|
| 41 | /* Interval over which precession required, in tropical centuries */
|
|---|
| 42 | t = ( bep1 - bep0 ) / 100.0;
|
|---|
| 43 |
|
|---|
| 44 | /* Euler angles */
|
|---|
| 45 | tas2r = t * DAS2R;
|
|---|
| 46 | w = 2303.5548 + ( 1.39720 + 0.000059 * bigt ) * bigt;
|
|---|
| 47 | zeta = (w + ( 0.30242 - 0.000269 * bigt + 0.017996 * t ) * t ) * tas2r;
|
|---|
| 48 | z = (w + ( 1.09478 + 0.000387 * bigt + 0.018324 * t ) * t ) * tas2r;
|
|---|
| 49 | theta = ( 2005.1125 + ( - 0.85294 - 0.000365* bigt ) * bigt +
|
|---|
| 50 | ( - 0.42647 - 0.000365 * bigt - 0.041802 * t ) * t ) * tas2r;
|
|---|
| 51 |
|
|---|
| 52 | /* Rotation matrix */
|
|---|
| 53 | slaDeuler ( "ZYZ", -zeta, theta, -z, rmatp );
|
|---|
| 54 | }
|
|---|