practical_astronomy_rust/
binary.rs

1use crate::binarydata as pa_bd;
2use crate::macros as pa_m;
3use crate::util as pa_u;
4
5/// Calculate orbital data for binary star.
6///
7/// ## Arguments
8/// * `greenwich_date_day` -- Greenwich date (day)
9/// * `greenwich_date_month` -- Greenwich date (month)
10/// * `greenwich_date_year` -- Greenwich date (year)
11/// * `binary_name` -- Abbreviated name of binary
12///
13/// ## Returns
14/// * `position_angle_deg` -- Position angle (degrees)
15/// * `separation_arcsec` -- Separation of binary members (arcseconds)
16pub fn binary_star_orbit(
17    greenwich_date_day: f64,
18    greenwich_date_month: u32,
19    greenwich_date_year: u32,
20    binary_name: String,
21) -> (f64, f64) {
22    let (binary_info, _binary_info_status) = pa_bd::get_binary_info_vector(binary_name);
23
24    let y_years = (greenwich_date_year as f64
25        + (pa_m::cd_jd(
26            greenwich_date_day,
27            greenwich_date_month,
28            greenwich_date_year,
29        ) - pa_m::cd_jd(0.0, 1, greenwich_date_year))
30            / 365.242191)
31        - binary_info.epoch_peri;
32    let m_deg = 360.0 * y_years / binary_info.period;
33    let m_rad = (m_deg - 360.0 * (m_deg / 360.0).floor()).to_radians();
34    let eccentricity = binary_info.ecc;
35    let true_anomaly_rad = pa_m::true_anomaly(m_rad, eccentricity);
36    let r_arcsec = (1.0 - eccentricity * (pa_m::eccentric_anomaly(m_rad, eccentricity)).cos())
37        * binary_info.axis;
38    let ta_peri_rad = true_anomaly_rad + binary_info.long_peri.to_radians();
39
40    let y = (ta_peri_rad).sin() * ((binary_info.incl).to_radians()).cos();
41    let x = (ta_peri_rad).cos();
42    let a_deg = pa_m::degrees(y.atan2(x));
43    let theta_deg1 = a_deg + binary_info.pa_node;
44    let theta_deg2 = theta_deg1 - 360.0 * (theta_deg1 / 360.0).floor();
45    let rho_arcsec =
46        r_arcsec * (ta_peri_rad).cos() / ((theta_deg2 - binary_info.pa_node).to_radians()).cos();
47
48    let position_angle_deg = pa_u::round_f64(theta_deg2, 1);
49    let separation_arcsec = pa_u::round_f64(rho_arcsec, 2);
50
51    return (position_angle_deg, separation_arcsec);
52}