rfa/vector_matrix/spherical_cartesian/
pv2s.rs

1
2use crate::utils::*;
3///  Convert position/velocity from Cartesian to spherical coordinates.
4///
5///  Given:
6///   * pv pv-vector
7///
8///  Returned:
9///   * theta    double        longitude angle (radians)
10///   * phi      double        latitude angle (radians)
11///   * r        double        radial distance
12///   * td       double        rate of change of theta
13///   * pd       double        rate of change of phi
14///   * rd       double        rate of change of r
15///
16///  Notes:
17///   * If the position is a pole, theta, td and pd are indeterminate.
18///     In such cases zeroes are returned for all three.
19///
20///  This revision:  2021 May 11
21pub fn pv2s(pv: &[[f64; 3]; 2], theta: &mut f64, phi: &mut f64, r: &mut f64,
22    td: &mut f64, pd: &mut f64 , rd:& mut f64)
23{
24    /* Components of position/velocity vector. */
25    let mut x  = pv[0][0];
26    let mut y  = pv[0][1];
27    let mut z  = pv[0][2];
28    let xd = pv[1][0];
29    let yd = pv[1][1];
30    let zd = pv[1][2];
31    
32    /* Component of r in XY plane squared. */
33    let mut rxy2 = x*x + y*y;
34    
35    /* Modulus squared. */
36    let mut r2 = rxy2 + z*z;
37    
38    /* Modulus. */
39    let rtrue = r2.sqrt();
40    
41    /* If null vector, move the origin along the direction of movement. */
42       let mut rw = rtrue;
43       if rtrue == 0.0 {
44           x = xd;
45           y = yd;
46           z = zd;
47           rxy2 = x*x + y*y;
48           r2 = rxy2 + z*z;
49           rw = r2.sqrt();
50       }
51    
52    /* Position and velocity in spherical coordinates. */
53       let rxy = sqrt(rxy2);
54       let xyp = x*xd + y*yd;
55       if rxy2 != 0.0 {
56           *theta = atan2(y, x);
57           *phi = atan2(z, rxy);
58           *td = (x*yd - y*xd) / rxy2;
59           *pd = (zd*rxy2 - z*xyp) / (r2*rxy);
60       } else {
61           *theta = 0.0;
62           *phi = if z != 0.0{atan2(z, rxy)}else{ 0.0};
63           *td = 0.0;
64           *pd = 0.0;
65       }
66       *r = rtrue;
67       *rd = if rw != 0.0{ (xyp + z*zd) / rw } else { 0.0 };
68    
69    /* Finished. */
70    
71}
rfa/vector_matrix/spherical_cartesian/pv2s.rs

rfa/vector_matrix/spherical_cartesian/
pv2s.rs
