rfa 0.5.9

A port ERFA to Rust.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88

use crate::rfam::*;
use crate::vector_matrix::vector_ops::{pxp::*, pn::*};
use crate::prec_nut::{ltpequ::*, ltpecl::*};
///  ICRS equatorial to ecliptic rotation matrix, long-term.
///
///  # Given:
///    * epj Julian epoch (TT)
///
///  # Returned:
///   * rm  ICRS to ecliptic rotation matrix
///
///  # Notes:
///
///  1) The matrix is in the sense
///
///        E_ep = rm x P_ICRS,
///
///     where P_ICRS is a vector with respect to ICRS right ascension
///     and declination axes and E_ep is the same vector with respect to
///     the (inertial) ecliptic and equinox of epoch epj.
///
///  2) P_ICRS is a free vector, merely a direction, typically of unit
///     magnitude, and not bound to any particular spatial origin, such
///     as the Earth, Sun or SSB.  No assumptions are made about whether
///     it represents starlight and embodies astrometric effects such as
///     parallax or aberration.  The transformation is approximately that
///     between mean J2000.0 right ascension and declination and ecliptic
///     longitude and latitude, with only frame bias (always less than
///     25 mas) to disturb this classical picture.
///
///  3) The Vondrak et al. (2011, 2012) 400 millennia precession model
///     agrees with the IAU 2006 precession at J2000.0 and stays within
///     100 microarcseconds during the 20th and 21st centuries.  It is
///     accurate to a few arcseconds throughout the historical period,
///     worsening to a few tenths of a degree at the end of the
///     +/- 200,000 year time span.
///  References:
///
///    Vondrak, J., Capitaine, N. and Wallace, P., 2011, New precession
///    expressions, valid for long time intervals, Astron.Astrophys. 534,
///    A22
///
///    Vondrak, J., Capitaine, N. and Wallace, P., 2012, New precession
///    expressions, valid for long time intervals (Corrigendum),
///    Astron.Astrophys. 541, C1
///
///  This revision:  2021 May 11
pub fn ltecm(epj: f64, rm:&mut [[f64; 3];3])
{
/* Frame bias (IERS Conventions 2010, Eqs. 5.21 and 5.33) */
    const DX: f64 = -0.016617 * URSA_DAS2R;
    const DE: f64 = -0.0068192 * URSA_DAS2R;
    const DR: f64 = -0.0146 * URSA_DAS2R;

/* Equator pole. */
   let mut p = [0.;3];

   ltpequ(epj, &mut p);

/* Ecliptic pole (bottom row of equatorial to ecliptic matrix). */
   let mut z = [0.;3];

   ltpecl(epj, &mut z);

/* Equinox (top row of matrix). */
   let mut w = [0.;3];

   pxp(&p, &z, &mut w);
   let mut s = 0.;
   let mut x = [0.;3];

   pn(&w, &mut s, &mut x);

/* Middle row of matrix. */
   let mut y = [0.;3];
   pxp(&z, &x, &mut y);

/* Combine with frame bias. */
   rm[0][0] =   x[0]    - x[1]*DR + x[2]*DX;
   rm[0][1] =   x[0]*DR + x[1]    + x[2]*DE;
   rm[0][2] = - x[0]*DX - x[1]*DE + x[2];
   rm[1][0] =   y[0]    - y[1]*DR + y[2]*DX;
   rm[1][1] =   y[0]*DR + y[1]    + y[2]*DE;
   rm[1][2] = - y[0]*DX - y[1]*DE + y[2];
   rm[2][0] =   z[0]    - z[1]*DR + z[2]*DX;
   rm[2][1] =   z[0]*DR + z[1]    + z[2]*DE;
   rm[2][2] = - z[0]*DX - z[1]*DE + z[2];
}