multi_target_radar.py

# multi_target_radar.py - Test multiple target OD.
# Copyright (C) 2021-2023 University of Texas
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.

import json
from orbdetpy import configure, add_station, Constant, MeasurementType
from orbdetpy.conversion import get_J2000_epoch_offset
from orbdetpy.estimation import multi_target_OD
from orbdetpy.propagation import propagate_orbits
from orbdetpy.plotting.estimation import plot

t0, t1 = get_J2000_epoch_offset(("2019-07-10T23:30:00", "2019-07-11T00:00:00"))

config = [configure(prop_start=t0, prop_end=t1, prop_step=60.0, prop_initial_state=(-20000151.1484, 21900.3443, 0.2039, 0, -4464.9556, 0),
                    sim_measurements=True),
          configure(prop_start=t0, prop_end=t1, prop_step=60.0, prop_initial_state=(-20000000, 0, 0, -100, -4464.3029, 0),
                    sim_measurements=True),
          configure(prop_start=t0, prop_end=t1, prop_step=60.0, prop_initial_state=(-20000000, -20000, -20000, -100, -4464.3029, 0),
                    sim_measurements=True)]

# Simulate RA and dec
for cfg in config:
    add_station(cfg, "Arecibo", 0.3202, -1.1651, 497.0)
    add_station(cfg, "Kaena", 0.08683, -2.5872, 460.0)
    add_station(cfg, "Boston", 0.7495, -1.2504, 460.0)
    cfg.measurements[MeasurementType.RIGHT_ASCENSION].error[:] = [Constant.ARC_SECOND_TO_RAD]
    cfg.measurements[MeasurementType.DECLINATION].error[:] = [Constant.ARC_SECOND_TO_RAD]

obs_list = propagate_orbits(config)
merge_optical = [sorted((y for x in obs_list for y in list(x.array)), key=lambda x: f"{x.time}{x.station}")]

# Simulate range and range rate
for cfg in config:
    del cfg.measurements[MeasurementType.RIGHT_ASCENSION]
    del cfg.measurements[MeasurementType.DECLINATION]
    cfg.measurements[MeasurementType.RANGE].error[:] = [10.0]
    cfg.measurements[MeasurementType.RANGE_RATE].error[:] = [0.1]

radar_list = propagate_orbits(config)
merge_radar = [sorted((y for x in radar_list for y in list(x.array)), key=lambda x: f"{x.time}{x.station}")]

# Combine radar and optical data
for o, r in zip(merge_optical[0], merge_radar[0]):
    r.angle_rates[:] = o.values

fit_data = multi_target_OD(config[:2], merge_radar)

print(f"Unassociated observations: {fit_data.unassociated_obs}")
for obs, fit in zip(fit_data.associated_obs, fit_data.est_output):
    assoc = [merge_radar[0][idx] for idx in obs.array]
    if (len(assoc) > 0):
        plot(config[0], assoc, fit.array, interactive=True, estim_param=False)