RFI Light-Curve Estimation

Measuring each satellite’s apparent flux over time and frequency straight from the visibilities, by matched-filtering them against the known trajectory phase. No imaging is involved, and the output is the same interchange format an imager would have to produce, so the two are interchangeable as seeds for rfi.est.

The estimator

For a satellite on a known trajectory the RFI contribution to baseline \((p, q)\) is \(A_p A_q^* e^{i(\phi_p - \phi_q)}\) with \(\phi\) the geometric phase, so the unit-modulus template \(T_{pq} = e^{i(\phi_p - \phi_q)}\) de-rotates it and the maximum-likelihood estimate of the source visibility is the inverse-variance-weighted, de-rotated baseline average

\[\hat{S}[f, t] = \frac{\sum_{pq} w_{pq} T_{pq}^{*} V_{pq}}{\sum_{pq} w_{pq}}, \qquad w_{pq} = \frac{1}{\sigma_{pq}^2},\]

with standard error \(1 / \sqrt{\sum w}\) and significance \(z = \mathrm{Re}(\hat{S}) / \text{error}\).

The weights are the noise the MS reports, resolved per baseline and per channel as far as the column resolves it (see tabascal.noise); the template carries no gain. That division of labour is deliberate — see the module docstring below.

Which baselines are coherent

A baseline only helps if the template phase is right on it, and two independent effects put a ceiling on that. A transverse orbit error \(\delta\) at slant range \(r\) shifts the apparent direction by \(\delta / r\); keeping the phase that costs below a radian gives

\[2 \pi \frac{b}{\lambda} \frac{\delta}{r} \le 1 \quad \Longrightarrow \quad b \le \frac{\lambda r}{2 \pi \delta}.\]

A satellite crossing at \(v_\perp\) sweeps the baseline fringe at \((b/\lambda)(v_\perp/r)\), which a model averaging \(N\) fine steps per integration \(\Delta t\) can follow only to the Nyquist rate of its own grid:

\[\frac{b}{\lambda} \frac{v_\perp}{r} \le \frac{N}{2 \Delta t} \quad \Longrightarrow \quad b \le \frac{\lambda r N}{2 \Delta t \, v_\perp}.\]

The smaller of the two binds. On the MWA Cen A field (175 MHz, 567 km range) a 600 m baseline tolerates \(\delta \approx 258\) m while the full 5.3 km array needs \(\approx 29\) m – far tighter than the 0.1-1 km transverse error a Starlink TLE carries. The phase-coherent search over all 9180 baselines therefore failed, the long ones adding with random phase and diluting the statistic, while the same filter over the 1004 baselines under 600 m recovered the satellite at 5.6 sigma.

coherent_baseline_mask() applies both criteria to the lengths from baseline_lengths(), keeping a baseline where \(b \le \min(b_\mathrm{TLE}, b_\mathrm{fringe})\). Both are always in play: n_fine, delta_t and v_perp_m_s are required arguments beside the orbit error, and v_perp_m_s = 0 sends \(b_\mathrm{fringe}\) to infinity, reducing the cut to the TLE criterion alone – the escape hatch for a stationary emitter, or for a caller who wants the orbit ceiling by itself. With soft=True the step becomes Gaussian weights \(e^{-(b/b_\mathrm{coh})^2}\) on the same scale.

Along-track time offset

A TLE’s dominant error is along-track, and along the track an error is very nearly a pure time offset: the satellite is where the elements say it will be \(\tau\) seconds later. One scanned parameter therefore recovers most of the error budget. For each \(\tau\) on a grid the orbit is evaluated at \(t + \tau\), the near-field fringe model is built on \(N\) fine steps inside each integration and averaged over them, and the coherent baselines are summed against the data:

\[z = \sum_{pq} w V M^{*}, \qquad n_1 = \sum_{pq} w |V|^2, \qquad n_2 = \sum_{pq} w |M|^2, \qquad r = \frac{|z|}{\sqrt{n_1 n_2}},\]

with \(r \in [0, 1]\) a per-frame correlation from which the intra-dump fringe smearing divides out. Frames are combined incoherently – the emitter’s own phase is not modelled between integrations – into \(z^2 = \sum_t |z|^2 / (n_1 n_2)\) per channel, and the best cell over \((\tau, \mathrm{channel})\) is the measurement. Only the satellite moves with \(\tau\): the antennas, the sidereal angle and the phase tracking stay at the observation’s own times, because \(\tau\) is an error in the orbit and not in the clock.

The significance comes from a decohered null – the same statistic at the best \(\tau\) with each antenna’s path pushed by an independent \(U(0, 50\,\mathrm{m})\), tens of wavelengths, so every baseline enters with an unrelated phase and the coherent sum collapses to an incoherent one. Two hundred draws give \((z^2_\mathrm{best} - \langle z^2 \rangle_\mathrm{null}) / \sigma_\mathrm{null}\). Nothing about the real distribution of \(z^2\) here is analytic, which is why the null is drawn on the data themselves, carrying their own weights, flagging, baseline set and residual sky.

fit_time_offset() is the whole measurement and tabascal light-curve --fit-offset exposes it. Two caveats travel with the number, both deliberate:

It carries no trials factor. The scan maximises over the whole grid and every channel while the null maximises over channels at the best \(\tau\) alone, so the significance is biased high and grows with the size of the grid searched. The 5 sigma default is a working cut calibrated on the MWA Cen A case, not a false-alarm probability.
The step must resolve the peak. Its half-width scales like \(\lambda r / (2 b_\mathrm{coh} v_\perp)\) – about 0.1 s for a 600 m coherent array at 567 km – so a coarser grid steps over the detection. The 0.25 s default matched the MWA curve, which decays over \(\pm 2\) s because the shortest baselines dominate that sum; a longer coherent array wants a finer step, not a wider grid.

The core (near_field_fringe_model(), matched_filter_sums(), coherence_scores(), tau_scan()) is pure jax.numpy over fixed-shape arrays, walking the grid with lax.map so one compilation covers the whole scan, and is left undecorated so the drivers own the jit and the batched identification search can vmap it over candidates. shift_orbit_record_epoch() closes the loop: an orbit record whose epoch is moved by \(-\tau\) reproduces the measured trajectory through --extra-orbit-dir with no further code.

Searching across candidates

With a TLE snapshot and no prior knowledge of which satellite is in the data, the same scan is run over every candidate at once: enumerate_candidates() screens the records down to the ones that were above the horizon, recording per candidate the frames it was up for and its slant range at maximum elevation; search_candidates() scores them; and select_detections() reads the ranking. Three things make that the same statistic as the single-satellite fit rather than a second one.

One baseline set, each candidate’s own cut. vmap needs static shapes, so the search sums over a single baseline list: the union of the hard coherent sets of the candidates above the geometric horizon. A union, and not the farthest candidate’s set, because the sets are not nested by range – \(b_\mathrm{TLE} \propto r\) alone, while \(b_\mathrm{fringe}\) turns on each candidate’s own transverse speed, so a nearer, slower satellite can steer a baseline that a farther, faster one cannot. Each candidate then applies its own coherence as a per-baseline weight inside the statistic, so another candidate’s excess baselines enter at exactly zero and it is scored over precisely the baselines it could steer. The union is sized from the above-horizon candidates alone: a satellite 13 000 km away, through the Earth, tolerates kilometres of baseline and would otherwise readmit the long ones for everybody – and a search in which no candidate is above the horizon has no honest set to sum over at all, so it returns nothing rather than a ranking.

The cut is sized from one geometry taken at one instant: the mid-window \((r, v_\perp)\) pair over the frames the candidate is in view for, which is the pair fit_time_offset() uses and which the search reports per candidate in fits[i]. The range_m of a candidate, and of a ranking row, is a different number – the closest approach during the pass, reported because it says how near the satellite came. With soft_weights the support is still the hard cut and the Gaussian weights the baselines inside it: the taper is never exactly zero, so a support read off the weights would be every baseline the array has, and would depend on the precision the scan happened to run in.

The horizon inside the statistic. Each candidate’s in-view mask is passed to tau_scan() as frame_mask rather than slicing the arrays, so a satellite that rises or sets mid-observation contributes only its own frames while the batch stays rectangular. Masking with zeros and slicing give the same \(z^2\), so the search and fit_time_offset() report the same detection for the same pass.

One compilation. jax.jit(jax.vmap(tau_scan)) is held at module level and the candidates are fed to it in batches, a ragged last batch padded by repeating its last candidate so every call has one shape. Two arrays per candidate dominate the memory – the fringe model (n_bl, n_freq, n_time, n_fine) complex, one offset at a time, and the paths (n_tau, n_bl, n_time, n_fine) float64 – and max_mem_gb is a budget for their sum, so the batch actually run is the smaller of batch_size and what that budget affords (reported back as batch_size). At MWA scale it is the budget that decides: the union reaches 7704 of the array’s 9180 baselines once candidates come near the horizon, one candidate over 24 channels is then some 2.1 GB, and a batch of eight would ask for 17 GB.

The null is drawn for the top n_null_candidates only: two hundred extra scans per satellite over a whole constellation is the search twice over, spent on candidates nothing will be reported for. That shortlist is taken on raw \(z^2\), which is a sum over in-view frames, so a short pass ranks below a full one at the same per-frame correlation – a caveat on the shortlist rather than a correction to make.

select_detections() carries two warnings: a close runner-up within runner_up_ratio of the winner, since satellites in the same train partially match each other’s fringes and a winner that is not clear of the field is a result to look at twice; and a detected candidate whose best \(\tau\) sits on the first or last grid point, whose offset is then a floor rather than a measurement, the remedy being a wider --tau-max. The deliverables are write_config_fragment() – the satellites.norad_ids list, beside the epoch-shifted records it can be replayed from – write_search_results() for the ranking table, and plot_candidate_ranking() for the chart a named satellite is judged against.

Where it is used

rfi.init / rfi.mean: matched-filter seeds the RFI signal model from the visibilities a run has already loaded, through light_curves_from_config().
tabascal light-curve writes the same estimate to an .npz.
tabascal light-curve --fit-offset measures each satellite’s along-track offset first, extracts the curves at it, and records it in the output.
tabascal light-curve -z filters a run’s residual, taken from its results zarr, as a post-fit diagnostic.
tabascal search finds the satellites in an observation from a TLE snapshot alone and emits the satellites.norad_ids a run needs.

Matched-filter light-curve extraction for RFI sources.

Given a set of satellite trajectories and the observed visibilities, this module beam-forms the interferometer toward each satellite and reads off its per-timestep (and per-channel) flux. This is a matched filter in visibility space: for a point source moving along a known trajectory the RFI contribution to baseline (p, q) is:

V_rfi[bl] = A_p conj(A_q) exp(i (phi_p - phi_q))

where phi_a is the geometric phase at each antenna (tabascal.interferometry.get_rfi_phase()). The per-baseline template is the unit-modulus steering vector T_bl = exp(i (phi_p - phi_q)), the data are V_bl = T_bl S + n_bl with per-component noise sigma_bl, and the maximum-likelihood (inverse-variance-weighted) estimate of the source visibility at each (freq, time) is the de-rotated, weighted baseline average:

S_hat[f, t] = sum_bl w_bl conj(T_bl) V_bl / sum_bl w_bl |T_bl|^2,
w_bl = 1 / sigma_bl^2

with variance 1 / sum_bl w_bl |T_bl|^2. Every template here is unit-modulus, so the denominator is just D = sum_bl w_bl and:

error = 1 / sqrt(D),      z = Re(S_hat) / error

z is a calibrated-frame statistic: it reads the real part because a de-rotated real source has no imaginary part to read, which holds only where the antenna gain phases have been taken out. coverage_stats() reports |S_hat| / error beside it as amp_coverage, against a matched Rayleigh threshold. That magnitude is invariant to a phase common to every baseline – an overall offset, or a stable phase on the source itself – which would otherwise turn the signal out of the real part and hide it from z.

It is not a defence against an uncalibrated antenna gain. A gain multiplies each baseline before the average, S_hat = S * sum_bl w g_p conj(g_q) / sum_bl w, so antenna-dependent phases decorrelate the coherent sum itself: the estimate shrinks, and there is nothing left in its magnitude for either statistic to find. Calibrate the phases, or accept that both numbers understate what is there.

The satellite fringe adds coherently after de-rotation while the sky and the noise add incoherently, so S_hat isolates the RFI source visibility. Its magnitude is a per-antenna power estimate; sqrt(|S_hat|) is the per-antenna amplitude used to seed tabascal’s RFI signal model.

The weights carry the calibration; the template does not. w comes from the noise the MS reports, resolved per baseline and per channel wherever the column resolves it that far (SIGMA_SPECTRUM / SIGMA, see tabascal.noise). On an uncalibrated column the source is really g_p conj(g_q) S, so a unit template de-rotates the geometry but not the gains: the per-baseline terms add with a scatter of gain phases and the coherent sum is degraded. That loss is the cost of not calibrating, and it is what the estimate should show. Putting the gain in the template and calibrating the data while carrying the transformed noise (WEIGHT_SPECTRUM = |g|^2 / SIGMA^2, the frame the TAB_* columns are written in) are the same estimator; what is not allowed is calibrating and then weighting uniformly. There is no noise-vs-gain power law anywhere here: the noise is whatever the MS says it is.

Three entry points are provided:

light_curves_from_config() – the in-process path. Reuses the arrays a TabConfig has already loaded (visibilities, noise, antennas, times, orbit records) so tabascal can seed the RFI model without touching the MS again, and returns the curves already ordered to match satellites.norad_ids.
extract_light_curves_from_ms() – the standalone tool. Point it at any MS column and a set of NORAD IDs and it returns / saves the light curves.
extract_light_curves_from_zarr() – the post-fit diagnostic. Matched filters the residual of a run taken straight from its results zarr.

All three share the pure core matched_filter_light_curves().

Which baselines the filter should be summed over is a question of its own, and coherent_baseline_mask() answers it from the orbit accuracy and the fringe rate rather than from a hand-tuned cut; it is a pure array function so the jitted tau-scan and identification search (#190/#191) can apply the same selection inside their own traces.

The along-track offset. A TLE’s dominant error is along-track – kilometres to tens of kilometres of drag mismodelling and unannounced manoeuvres – and an along-track error is very nearly a pure time offset in the trajectory. So one scanned parameter, tau, recovers the bulk of it: evaluate the orbit at t + tau, build the near-field fringe model on a fine grid inside each integration, and coherently correlate it against the data over the baselines coherent_baseline_mask() keeps. fit_time_offset() is the whole measurement – horizon window, coherence cut, scan, best cell and a decohered-antenna null for its significance – and tabascal light-curve --fit-offset exposes it, extracting the curves at the offset it measured and recording that offset in the output. Its core (near_field_fringe_model(), matched_filter_sums(), coherence_scores(), tau_scan()) is pure jax.numpy over fixed-shape arrays, scanning the grid with lax.map and undecorated so the drivers own the jit: one compilation covers the whole scan, and the batched identification search of #191 vmaps the same function over candidates. shift_orbit_record_epoch() is the other end of it – an orbit record moved by -tau, which reproduces the measured trajectory through extra_orbit_dir with no further code.

Which satellite it is. Given a TLE snapshot and nothing else, enumerate_candidates() screens the records down to the ones that were above the horizon, search_candidates() runs that same scan over all of them – jax.vmap over a candidate axis, one jitted program per batch, each candidate’s horizon mask and coherence cut applied inside the statistic so the shapes stay static – and select_detections() reads the ranking against the decohered null. tabascal search is the command, and what it emits is the satellites.norad_ids list a run needs (write_config_fragment()), the ranking table (write_search_results()), and the light curves and shifted orbit records of whatever it named. That is GitHub #191, which produced STARLINK-1765 out of 2551 records on the MWA Cen A dataset.

Only the satellite-trajectory source is implemented. RA/Dec and Alt/Az pointings can be added by constructing rfi_xyz from those and feeding rfi_phase_from_positions().

tabascal.rfi_estimate.DEFAULT_TAU_GRID = array([-4. , -3.75, -3.5 , -3.25, -3. , -2.75, -2.5 , -2.25, -2. , -1.75, -1.5 , -1.25, -1. , -0.75, -0.5 , -0.25, 0. , 0.25, 0.5 , 0.75, 1. , 1.25, 1.5 , 1.75, 2. , 2.25, 2.5 , 2.75, 3. , 3.25, 3.5 , 3.75, 4. ]): Along-track offsets fit_time_offset() scans by default, in seconds: +-4 s in 0.25 s steps, 33 points including 0 and both ends. Wide enough for a day-old Starlink TLE (the MWA Cen A case sat at -2.25 s on a TLE 1.4 h old) and fine enough to resolve the peak an array of a few hundred metres gives it – see fit_time_offset() on sizing the step.

tabascal.rfi_estimate.attach_offset_fits(result: dict, fits: list, threshold_sigma: float) → dict[source]

Add a per-satellite offset fit to a light-curve result.

The output artifact has to record the offset the curves were measured at, or a later run cannot reproduce the trajectory that produced them – so the scan’s answers travel with the curves, stacked along the same source axis, and save_light_curves_npz() writes them into the .npz.

detected is decided here, once, against threshold_sigma, and the threshold is recorded alongside it so a file says what it was judged by.

Parameters:

result (dict) – A light-curve result (see _lc_result()).
fits (list of dict, length n_src) – Fits from fit_time_offset(), in the result’s source order.
threshold_sigma (float) – Significance above the null at which a fit counts as a detection.

Returns:

A copy of result carrying tau_best (n_src,), tau_grid (n_tau,), z2_tau (n_src, n_tau, n_freq), z2_best, best_chan, significance, null_mean, null_std, detected (n_src,) bool, r_best (n_src, n_freq, n_time), offset_threshold_sigma, and the fits themselves under offset_fits for the diagnostics to draw from.

Return type:

dict

Physical separation of each antenna pair, in metres.

Both coherence criteria act on the baseline component perpendicular to the line of sight to the satellite, which changes as it crosses the sky. The physical length bounds that component from above in every direction, so it is the conservative choice and needs no per-timestep geometry. The uv distance is not a substitute: it is the projection toward the phase centre, not toward the satellite, and would admit baselines the satellite sees at full length.

Parameters:

ants_itrf (Array (n_ant, 3)) – Antenna positions in the ITRF (ECEF) frame, in metres.
a1 (Array (n_bl,)) – Antenna indices of each baseline.
a2 (Array (n_bl,)) – Antenna indices of each baseline.

Returns:

Baseline length in metres.

Return type:

Array (n_bl,)

tabascal.rfi_estimate.candidates_from_norad_ids(norad_ids, times_jd: NDArray, extra_orbit_dir=None)[source]

Orbit records for an explicit candidate list, in the order it was asked.

The other way into the search: -n/-np, resolved through the run’s own fetch_orbital_elements(), so extra_orbit_dir, the managed cache and SatChecker take exactly the precedence a run gives them and the search scores the records a run would model.

Returns:: Records, names and NORAD IDs, in the requested order.
Return type:: (list of dict, list of str, list of int)

tabascal.rfi_estimate.candidates_from_orbit_dir(directory: str, times_jd: NDArray, name_filter=None)[source]

Every satellite in a local snapshot directory, one record apiece.

The --tle-dir source: a constellation export dropped in extra_orbit_dir format, read through the same per-ID nearest-epoch policy a run applies to extra_orbit_dir (tabascal.orbit._select_from_extra_dir() with no age ceiling) – so a file carrying several epochs for one satellite contributes the record nearest the observation, rather than whichever row it happened to list first.

That resolver prints a provenance line per satellite, which is the right thing for the handful a run configures and the wrong thing for a whole constellation, so only the lines reporting a rejected record are passed through. The ranking table is the output here, and it should not arrive under 2551 lines of bookkeeping.

Parameters:

directory (str) – Directory of orbit files (TLE or OMM), as --extra-orbit-dir takes.
times_jd (Array (n_time,)) – Observation times as UTC Julian dates; their mean is the epoch records are chosen against.
name_filter (str, optional) – Case-insensitive substring of OBJECT_NAME, e.g. "STARLINK". A snapshot is usually a whole constellation plus whatever else the query dragged in. None keeps all. A record with no name is named by its ID, and so is kept only by a filter that matches the ID.

Returns:

Records, names and NORAD IDs, in ascending ID order.

Return type:

(list of dict, list of str, list of int)

The per-frame correlation and the per-channel score it combines into.

r = |z| / sqrt(n1 n2) is a normalised correlation in [0, 1]: the intra-dump smearing that shrinks |M| on the longer baselines divides out of it, so a frame’s r says how well the data match the trajectory and not how bright the fringe was. z2 = sum_frames |z|^2 / (n1 n2) combines the frames incoherently – the satellite’s own phase is not modelled between integrations – and is therefore bounded by the number of frames in view, which is what makes it comparable against a null.

A cell nothing was measured in has n1 n2 = 0 and contributes exactly zero, with r = 0 there. Guarded with where rather than an epsilon in the denominator: an epsilon shifts every other cell’s value to spare this one.

Parameters:

z (Array (n_freq, n_time)) – The sums from matched_filter_sums().
n1 (Array (n_freq, n_time)) – The sums from matched_filter_sums().
n2 (Array (n_freq, n_time)) – The sums from matched_filter_sums().
frame_mask (Array (n_time,) bool or 0/1, optional) – Frames to combine. None combines all of them. It exists for the batched search of #191, whose shapes must stay static and so cannot slice the in-view window out; masking with zeros and slicing give the same z2, so the two paths report the same detection.

Returns:

r (n_freq, n_time) and z2 (n_freq,).

Return type:

(Array, Array)

Baselines a trajectory is accurate enough to beam-form with.

Two independent effects cap the baseline over which the template phase can be trusted: the orbit error, at b_tle = lam r / (2 pi delta) for a radian of phase (tle_coherence_length()), and the fringe rate against the Nyquist rate of the model’s own fine grid, at b_fringe = lam r n_fine / (2 delta_t v_perp) (fringe_rate_coherence_length()). They are unrelated, so the smaller of the two binds, element by element, and a baseline is kept when b <= min(b_tle, b_fringe). Both are always in play: all three fringe inputs are required, and a stationary emitter – or a caller who wants the orbit ceiling alone – passes v_perp_m_s = 0.0, which sends b_fringe to +inf so that it drops out of the minimum. There is no half-specified call to get wrong.

Beyond the binding length a baseline does not merely stop helping: it enters the coherent sum with an essentially random phase and dilutes the statistic the shorter baselines built. That is what the MWA Cen A case study measured. At 175 MHz and a 567 km slant range a 600 m baseline tolerates delta of about 258 m while the full 5.3 km array needs about 29 m, and Starlink TLEs carry 0.1-1 km of transverse error – so the detection lived entirely in the 1004 baselines under 600 m (5.6 sigma), while the phase-coherent search over all 9180 ranked the true satellite around 34th.

With soft=True the step is replaced by exp(-(b / b_coh)^2), a Gaussian taper on the same scale: unity at zero spacing, 1/e where the hard mask cuts. It down-weights the marginal baselines instead of discarding them, which is the gentler choice when sigma is itself uncertain. It is a weighting, not a wider cut: the Gaussian is never exactly zero, so the drivers keep the hard mask as the baseline support and apply these values inside it (_coherence_weights()). Reading the support off the weights instead would admit every baseline the array has.

One range suffices, taken at mid-observation: b_coh depends on r only linearly and the slant range varies by a few tens of percent over a pass, which moves the cut far less than the order-of-magnitude uncertainty on delta does.

Everything goes through jax.numpy, so this composes inside the jitted, GPU-resident matched-filter core that consumes it. soft is the one static argument – it selects the output dtype and so cannot be traced.

Only the two magnitudes are made sign-safe; freq, range_m, delta_t and n_fine are positive-domain quantities that are taken as given, since a negative frequency or dump time is a caller error rather than a convention to absorb.

Parameters:

bl_len (float or Array) – Baseline lengths in metres (see baseline_lengths()).
freq (float or Array) – See tle_coherence_length(). sigma_transverse_m is read as a magnitude.
range_m (float or Array) – See tle_coherence_length(). sigma_transverse_m is read as a magnitude.
sigma_transverse_m (float or Array) – See tle_coherence_length(). sigma_transverse_m is read as a magnitude.
n_fine (int, float or Array) – See fringe_rate_coherence_length(). Required, not optional: v_perp_m_s = 0.0 is how a caller asks for the TLE ceiling alone. v_perp_m_s is read as a magnitude.
delta_t (int, float or Array) – See fringe_rate_coherence_length(). Required, not optional: v_perp_m_s = 0.0 is how a caller asks for the TLE ceiling alone. v_perp_m_s is read as a magnitude.
v_perp_m_s (int, float or Array) – See fringe_rate_coherence_length(). Required, not optional: v_perp_m_s = 0.0 is how a caller asks for the TLE ceiling alone. v_perp_m_s is read as a magnitude.
soft (bool, default False) – Return Gaussian weights rather than a boolean mask.

Returns:

Boolean mask, or float weights where soft, at the broadcast shape of the inputs.

Return type:

Array

tabascal.rfi_estimate.coverage_stats(result: dict, z_crit: float = 3.0) → dict[source]

Fraction of time-frequency cells consistent with noise, per source.

The z statistic z = Re(S_hat) / error would be ~ N(0, 1) wherever nothing is left after subtraction, so a well-cleaned source has |z| within z_crit almost everywhere. coverage is the fraction of finite (freq, time) cells with |z| <= z_crit; max_z is the peak residual significance.

The z statistic assumes the data are phase calibrated. Re(S_hat) is the whole of a de-rotated real source only when nothing else rotates it, so read it on a calibrated column (CORRECTED_DATA, the TAB_* columns, or a residual against a fitted model).

amp_coverage is the same statistic on |S_hat|/error, against rayleigh_threshold(); its null is analytic – Rayleigh(1) – so it carries no excess column. It is invariant to a rotation common to every baseline: an overall phase offset, or a stable phase on the source, turns S_hat as a whole, which empties Re(S_hat) and spills the source into the imaginary part that null_coverage is measured on – both halves of that comparison then move the wrong way while the magnitude is untouched.

It is not immunity to an uncalibrated antenna gain. A gain multiplies each baseline before the average, S_hat = S * sum_bl w g_p conj(g_q) / sum_bl w, so antenna-dependent phases decorrelate the sum itself: the estimate shrinks and its magnitude with it, and both statistics drift toward “nothing here”. On a raw column neither number is a detection threshold so much as a lower bound. The optimistic-floor caveat below applies to both.

Compare against ``null_coverage``, not against the analytic 2*Phi(z)-1. The floor assumes the de-rotated per-baseline samples are independent. They are not: residual sky is coherent across baselines, so the floor is optimistic and the analytic null over-states the expected coverage. null_coverage is the same statistic on Im(S_hat)/error – after de-rotation a real source sits purely in the real part, so the imaginary part is a matched, source-free null carrying the same noise and the same correlation structure. A source is consistent with noise when its coverage is not significantly below the null; the excess null_coverage - coverage is the part attributable to a real residual.

Parameters:

result (dict) – A driver result, carrying z, light_curves and error.
z_crit (float, default 3.0) – Detection threshold; cells above it are flagged as residual.
whole (A source with no finite cell -- one held out of view for the)
cut (observation by the elevation)
to (or flagged away -- has no coverage)
per_source (report and comes back as nan. It is still listed in)
the (so)
it (table says it was not measured rather than omitting)
overall (but)
were (summarises only the sources that)
no (source nothing was measured for cannot be the worst-fitted one. With)
all (source measured at)
and (every overall coverage metric is nan)
(z_crit (worst_source is None; the thresholds)
are (amp_crit))
populated. (settings rather than measurements and stay)

Returns:

per_source (title, coverage, null_coverage, excess, amp_coverage, max_z, max_amp, n_cells) and overall (pooled coverage, null, amp_coverage, worst source, mean, z_crit, amp_crit). When nothing was measured, worst_source is None and every coverage metric is nan, while z_crit and amp_crit are unchanged.

Return type:

dict

What the statistic scores when the geometry is not real.

The same sums at the same offset, with each antenna’s path pushed by an independent U(0, jitter_m). Tens of metres are tens of wavelengths at metre wavelengths, so every baseline enters with an unrelated phase and the coherent sum collapses to an incoherent one – which is what “no satellite on that trajectory” looks like, measured on these data, with their own weights, flagging, baseline set and residual sky. That is why the null is drawn rather than taken from a chi-squared: nothing about the real distribution of z2 here is analytic.

Drawn with jax.random from PRNGKey(seed), so a significance is reproducible; walked with lax.map rather than vmap because a batched draw would hold n_draws copies of the fringe model at once.

Parameters:

vis (Array) – As tau_scan().
weights (Array) – As tau_scan().
freqs (Array) – As tau_scan().
paths_best (Array (n_bl, n_time, n_fine)) – Path differences at the best offset.
a1 (Array (n_bl,)) – Antenna indices; their maximum sizes the offset vector.
a2 (Array (n_bl,)) – Antenna indices; their maximum sizes the offset vector.
frame_mask (Array (n_time,), optional) – See coherence_scores().
n_draws (int, default 200) – Scrambles to draw.
jitter_m (float, default 50.0) – Upper end of the per-antenna offset, in metres.
seed (int, default 0) – PRNG seed.

Returns:

max-over-channel z2 for each draw.

Return type:

Array (n_draws,)

tabascal.rfi_estimate.enumerate_candidates(records: list, names: list, times_jd: NDArray, ants_itrf: NDArray, min_elevation: float | None = 0.0) → list[source]

Screen a snapshot down to the satellites worth scoring.

Stage one of the identification search of GitHub #191, and the step that turns a constellation into a shortlist: on the MWA Cen A case it took 2551 Starlink records to the 128 that were above the horizon during the 56 s observation. A satellite that never rose is not evidence, and scoring it would cost a scan.

It also decides which frames each candidate is scored over. The mask is computed once, here, and travels with the candidate: the batched scan applies it inside the statistic rather than slicing, because its shapes have to stay static to vmap.

The screen is what keeps the shared baseline set of search_candidates() honest, too. The coherence ceiling grows with the slant range, so a satellite on the far side of the Earth – 13 000 km away, through the ground – would lift it to kilometres and readmit every long baseline. Dropping it here is what stops that; a search whose candidates are all below the horizon has no set to sum over and returns nothing.

Elevations are evaluated once, at tau = 0: the offsets searched for are seconds, which moves a LEO satellite tens of kilometres along its track and its elevation by a fraction of a degree. Nothing in the cut resolves that.

Parameters:

records (list of dict, length n_records) – Orbit records – TLE or OMM, as resolved by tabascal.orbit. The record itself travels with the candidate rather than just its ID, since re-resolving later could pick up a different record for the same satellite.
names (list of str, length n_records) – Catalogue names, aligned with records.
times_jd (Array (n_time,)) – Integration centres as UTC Julian dates.
ants_itrf (Array (n_ant, 3)) – Antenna positions in the ITRF frame, in metres; the mean is the site.
min_elevation (float, optional, default 0.0) – Elevation in degrees at or above which a satellite counts as up, inclusive, as rfi.min_elevation is. A satellite reaching it in at least one integration is kept. None keeps every record over every frame.

Returns:

One per kept record, sorted by descending max_elevation, each with norad_id, name, record, max_elevation, elevation (n_time,) in degrees, frames (n_time,) bool, and range_m – the slant range at the frame of maximum elevation, i.e. the closest approach during the pass.

range_m is for reporting: it says how near the satellite came, which is what a reader of a ranking table wants. It is not what sizes the coherence cut. That is the mid-window (range, speed) pair satellite_range_and_speed() returns for the in-view frames, one geometry taken at one instant, which is what fit_time_offset() uses and what search_candidates() reports in its per-candidate fits. The two are half a pass apart.

Return type:

list of dict

tabascal.rfi_estimate.extract_light_curves_from_ms(ms_path: str, norad_ids: list | None = None, corr: str = 'xx', data_col: str = 'DATA', freq: float | None = None, exclude_autos: bool = True, extra_orbit_dir: str | None = None, min_elevation: float | None = 0.0, max_mem_gb: float = 1.0, offset_fit: dict | None = None) → dict[source]

Extract matched-filter RFI light curves from any column of an MS.

The standalone entry point, used by the tabascal light-curve CLI. Reads the requested data_col (and the MS’s own noise column, through the same tabascal.ms.read_ms() a run uses), propagates the satellites’ orbit records over the MS times, and runs matched_filter_light_curves().

Parameters:

ms_path (str) – Path to the Measurement Set.
norad_ids (list[int]) – NORAD catalogue IDs; their orbit records are resolved through tabascal.orbit, with the same source precedence as a run.
corr (str, default "xx") – Correlation to read (xx/xy/yx/yy).
data_col (str, default "DATA") – MS data column to matched-filter.
freq (float, optional) – If given, use only the single channel nearest this frequency (Hz).
exclude_autos (bool, default True) – Drop autocorrelations from the beam-former.
extra_orbit_dir (str, optional) – Extra local directory of orbit files, searched before the managed cache and SatChecker.
min_elevation (float, optional) – Elevation in degrees below which a satellite is not filtered for. None disables the cut.
max_mem_gb (float, default 1.0) – Memory budget for the matched-filter time-chunk loop.
offset_fit (dict, optional) – Settings for the along-track offset search; see attach_offset_fits(). With it, each satellite’s tau is measured first and the curves are extracted at the offset it found.

Returns:

See _lc_result(). light_curves is (n_src, n_freq, n_time) complex, ordered to match norad_ids.

Return type:

dict

tabascal.rfi_estimate.extract_light_curves_from_zarr(ms_path: str, zarr_path: str, norad_ids: list | None = None, corr: str = 'xx', data_col: str = 'DATA', freq: float | None = None, exclude_autos: bool = True, extra_orbit_dir: str | None = None, min_elevation: float | None = 0.0, max_mem_gb: float = 1.0, offset_fit: dict | None = None) → dict[source]

Matched-filter the residual of a tabascal run, taken from its results zarr.

This is the way to score a run. The MS result columns (TAB_RES_DATA et al.) are overwritten by every tabascal run, so scoring off the MS is only valid if those columns happen to belong to the run meant. The zarr is written once per run and per suffix, so a later run cannot invalidate it.

The residual is formed as data_col - zarr.vis_obs. The zarr’s vis_obs is the model’s own gained prediction, apply_gains(gains, vis_ast + vis_rfi), so this is exactly the residual tabascal.write.write_results_ms() would write, without the MS round trip.

The model is matched to the MS’s channels by frequency, so a freq that narrows the read to one channel still subtracts that channel’s model. See _model_on_ms_channels().

data_col is the reference column the residual is formed against (e.g. DATA), not a residual column. Everything else is as extract_light_curves_from_ms().

tabascal.rfi_estimate.fine_time_offsets(n_fine: int, delta_t: float) → NDArray[source]

Where inside an integration the fringe model is sampled, in seconds.

n_fine equal sub-steps spanning one dump, taken at their midpoints: ((k + 0.5) / n_fine - 0.5) * delta_t. Midpoints rather than edges, because an edge grid samples the boundary between two integrations twice and biases each average by half a sub-step. n_fine = 1 is then exactly the integration centre, which is where the forward model’s own template lives – so the model reduces to rfi_phase_from_records()’s at one step, and the two cannot drift apart.

Parameters:

n_fine (int) – Sub-steps per integration.
delta_t (float) – Integration (dump) time in seconds.

Returns:

Offsets from the integration centre, in seconds.

Return type:

Array (n_fine,) float64

tabascal.rfi_estimate.fit_time_offset(vis: NDArray, record, ants_itrf: NDArray, times_jd: NDArray, phase_centre: dict, freqs: NDArray, a1: NDArray, a2: NDArray, int_time: float, noise=None, flags: NDArray | None = None, taus_s=None, n_fine: int = 40, sigma_transverse_m: float = 300.0, soft_weights: bool = False, min_elevation: float | None = 0.0, exclude_autos: bool = True, n_null: int = 200, null_jitter_m: float = 50.0, seed: int = 0) → dict[source]

Measure one satellite’s along-track time offset from the visibilities.

The single-satellite search of GitHub #190, end to end: window the observation to the frames the satellite is up for, choose the baselines the orbit is accurate enough to beam-form with, score every offset on the grid, and calibrate the best score against a decohered-antenna null.

The statistic, per offset, channel and frame, is the normalised correlation of coherence_scores() over the coherent baselines; frames are combined incoherently into z2 per channel and the best cell is the largest of those over (tau, channel). Its significance is (z2_best - null_mean) / null_std.

Two caveats on that significance, both deliberate.

It carries no trials factor. The scan maximises over the whole offset grid and every channel, while the null is drawn at the best offset and maximises over channels only – so the number is biased high, and grows with the size of the grid it searched. The default threshold of 5 sigma (is_detection()) is calibrated against the MWA Cen A case on the default grid, and is a working cut rather than a false-alarm probability.

And the step has to resolve the peak. Its half-width scales like lam r / (2 b_coh v_perp) – about 0.1 s for a 600 m coherent array at 567 km – so an offset grid coarser than that steps over the detection. The 0.25 s default matched the MWA curve, which decays over about +-2 s because the sum there is dominated by the shortest baselines; a longer coherent array needs a finer step, not a wider grid.

Parameters:

vis (Array (n_bl, n_freq, n_time) complex) – Visibilities to search.
record (dict) – One orbit record.
ants_itrf (Array (n_ant, 3)) – Antenna positions in the ITRF frame, in metres.
times_jd (Array (n_time,)) – Integration centres as UTC Julian dates.
phase_centre (dict) – {"ra": <deg>, "dec": <deg>}.
freqs (Array (n_freq,)) – Channel frequencies in Hz.
a1 (Array (n_bl,)) – Antenna indices of each baseline.
a2 (Array (n_bl,)) – Antenna indices of each baseline.
int_time (float) – Integration (dump) time in seconds.
noise (float or Array, optional) – Per-component noise standard deviation, in any of the shapes matched_filter_light_curves() accepts. None weights every baseline equally.
flags (Array, optional) – True marks samples to exclude. Flagged visibilities are zeroed before the sums, because an MS carries inf and nan in them and 0 * nan is nan.
taus_s (Array (n_tau,), optional) – Offsets to scan, in seconds. Defaults to DEFAULT_TAU_GRID.
n_fine (int, default 40) – Sub-steps per integration in the fringe model.
sigma_transverse_m (float, default 300.0) – Transverse orbit error the coherence cut is sized by, in metres. Around the middle of what a Starlink TLE carries.
soft_weights (bool, default False) – Taper the baselines inside the cut with a Gaussian on the coherence length instead of weighting them all equally. The set summed over is the hard cut either way (_coherence_weights()), so n_bl_used does not depend on this; what changes is how much the marginal baselines are allowed to say.
min_elevation (float, optional, default 0.0) – Elevation in degrees below which the satellite is not searched for, inclusive, as rfi.min_elevation is. None uses every frame. The in-view window is sliced out before the scan, so nothing below the horizon costs anything; the core’s frame_mask is the equivalent for callers whose shapes must stay static.
exclude_autos (bool, default True) – Drop autocorrelations. They carry no path difference and so no fringe.
n_null (int, default 200) – Draws in the decohered null.
null_jitter_m (float, default 50.0) – Per-antenna scramble in the null, in metres.
seed (int, default 0) – PRNG seed for the null.

Returns:

tau_grid (n_tau,), z2_tau (n_tau, n_freq), tau_best, z2_best, best_chan, best_freq, r_best (n_freq, n_time) on the full time axis with nan out of view, frames (n_time,) bool, elevation (n_time,) deg, times_sec (n_time,), null (n_null,), null_mean, null_std, significance, n_bl_used, range_m, v_perp_m_s, b_coh, n_fine and sigma_transverse_m.

The decision is not among them: it needs a threshold, and a dict carrying one would have to guess what the caller means by a detection. See is_detection().

Return type:

dict

tabascal.rfi_estimate.fit_time_offsets(orbit_records: list, norad_ids: list, vis: NDArray, ants_itrf: NDArray, times_jd: NDArray, phase_centre: dict, freqs: NDArray, a1: NDArray, a2: NDArray, int_time: float, threshold_sigma: float = 5.0, **kwargs) → list[source]

Fit and report one along-track offset per satellite, in order.

A loop over fit_time_offset() that prints offset_fit_summary() as each satellite is measured, so a long run says what it found while it is still running. Every other keyword goes straight to the fit.

Longest baseline the model average itself can follow, in metres.

A satellite crossing at transverse speed v_perp sweeps the baseline fringe at (b / lam) (v_perp / r) hertz. A model that averages n_fine sub-steps over an integration of length delta_t samples that fringe on a grid of spacing delta_t / n_fine, so it can follow it only up to the Nyquist rate of its own grid, n_fine / (2 delta_t). Equating the two:

b_fringe = lam r n_fine / (2 delta_t v_perp)

Past it the template decoheres inside the integration against its own discretisation, however good the orbit is; the cure is more fine steps or a shorter dump, not a better TLE.

Parameters:

freq (float or Array) – See tle_coherence_length().
range_m (float or Array) – See tle_coherence_length().
n_fine (int or Array) – Fine sub-steps the model averages over per integration.
delta_t (float or Array) – Integration (dump) time in seconds.
v_perp_m_s (float or Array) – Satellite speed across the line of sight in m/s. Taken as a magnitude, for the same reason as sigma_transverse_m: a pass in the other direction fringes at the same rate. Zero is a stationary emitter, with no fringe to outrun, and returns an infinite length.

Returns:

Coherence length in metres, at the broadcast shape of the inputs.

Return type:

Array

tabascal.rfi_estimate.has_noise_scale(result: dict) → bool[source]

Whether an estimate carries a noise floor, i.e. whether z means anything.

False when the visibilities were filtered with no sigma to weight them by: the light curves are still there, but every error is nan and nothing downstream that divides by one has anything to say.

tabascal.rfi_estimate.is_detection(fit: dict, threshold_sigma: float = 5.0) → bool[source]

Whether a fit clears the null by threshold_sigma.

Kept out of fit_time_offset() on purpose: the fit measures, the caller decides. A pass that was never in view has a nan significance and is not a detection at any threshold.

tabascal.rfi_estimate.light_curves_from_config(tab_config, vis: NDArray | None = None, exclude_autos: bool = True, max_mem_gb: float = 1.0, offset_fit: dict | None = None) → dict[source]

Matched-filter light curves from an already-loaded TabConfig.

The in-process path: reuses the visibilities, noise, antenna positions, times and orbit records the config has loaded, so no second MS read is needed, and the curves come back ordered to match satellites.norad_ids – no title matching, and no light-curve file.

Curves are returned for the real satellites only. Under device sharding the source axis is padded with duplicates of the last satellite, and those rows are re-added as zeros by the seeding code that consumes this.

The elevation cut is the run’s own rfi.min_elevation mask, taken off the config rather than recomputed, so the estimate is masked exactly where the model is.

Parameters:

tab_config (tabascal.config.TabConfig) – A configured object exposing vis_obs, flags, noise, ants_itrf, times_jd, times_mjd, time_scale, freqs, phase_centre, a1, a2, orbit_records and norad_ids.
vis (Array (n_bl, n_freq, n_time), optional) – Visibilities to filter; defaults to tab_config.vis_obs.
exclude_autos (bool, default True)
max_mem_gb (float, default 1.0)
offset_fit (dict, optional) – Settings for the along-track offset search; see attach_offset_fits(). With it, each satellite’s tau is measured first and the curves are extracted at the offset it found.

Returns:

See _lc_result().

Return type:

dict

tabascal.rfi_estimate.matched_filter_light_curves(vis: NDArray, rfi_phase: NDArray, a1: NDArray, a2: NDArray, noise=None, flags: NDArray | None = None, in_view: NDArray | None = None, exclude_autos: bool = True, max_mem_gb: float = 1.0) → Tuple[NDArray, NDArray][source]

Beam-form the data toward each source to estimate its source visibility.

The estimator and its noise floor are the ones this module’s docstring derives. Evaluated in numpy/f64 (a one-shot host-side estimate) with an outer time-chunk loop sized so the (n_bl, n_freq, chunk) per-baseline arrays – the template, the masked visibilities, the weights and their products – stay within max_mem_gb.

max_mem_gb bounds those arrays, not the function’s peak. It does not count the per-antenna exp(-i phi) chunk or the copies fancy indexing makes of a partly-masked block, and it says nothing about the arrays held for the whole call: vis and rfi_phase as given, and the (n_src, n_freq, n_time) accumulator, which grows with the number of sources rather than with the chunk. Lowering it shrinks the loop’s working set and nothing else.

Parameters:

vis (Array (n_bl, n_freq, n_time) complex) – Observed visibilities (any MS data column).
rfi_phase (Array (n_src, n_ant, n_freq, n_time)) – Per-antenna geometric phase (see rfi_phase_from_positions()).
a1 (Array (n_bl,)) – Antenna indices of each baseline.
a2 (Array (n_bl,)) – Antenna indices of each baseline.
noise (float or Array, optional) – Per-component noise standard deviation, as TabConfig resolves it: a scalar, (n_bl,), (n_bl, n_freq), or a three-dimensional array whose axes match the visibilities or are 1. None weights every baseline equally, which under-weights the quiet baselines and over-weights the loud ones, and returns a nan error: without a sigma the weights are a shape rather than a variance, and 1 / sqrt(N) would be asserting a noise of 1 Jy that nobody wrote down. The callers warn when they fall back to it.
flags (Array (n_bl, n_freq, n_time) bool, optional) – True marks samples to exclude from the average.
in_view (Array (n_src, n_time) bool, optional) – False marks (source, timestep) pairs the source is not up for. Those times are skipped entirely – the template is never evaluated there – and come back as an exact zero, the same “no signal known” convention the forward model’s elevation mask uses.
exclude_autos (bool, default True) – Drop autocorrelation baselines (a1 == a2). They carry no fringe to de-rotate, only each antenna’s own power.
max_mem_gb (float, default 1.0) – Approximate cap on the working-array size of the time-chunk loop.

Returns:

Array (n_src, n_freq, n_time) complex – Matched-filter source-visibility estimate. nan where every baseline in a (freq, time) cell is flagged – nothing was measured there, and a zero would read as a measured zero – and exactly 0 where in_view says the source was not up.
Array (n_src, n_freq, n_time) real – The standard error of the (real) flux estimate, 1 / sqrt(sum_bl w): the beam-former’s noise floor, the visibility-space equivalent of a dirty-image aperture standard deviation. nan wherever the estimate is not a measurement, and everywhere when noise is None.

The three weighted inner products the coherence statistic is built from.

z = sum_bl w V conj(M), n1 = sum_bl w |V|^2, n2 = sum_bl w |M|^2, each summed over the baseline axis alone so the result is per channel and per frame: the satellite is coherent within an integration and the frames are combined afterwards, incoherently, by coherence_scores().

A plain weighted sum, with no flag handling of its own: a zero weight removes a baseline exactly, and whether a flagged sample may be nan is the driver’s problem (fit_time_offset() zeroes those before the sums, since 0 * nan is nan and would poison the whole cell).

Parameters:

vis (Array (n_bl, n_freq, n_time) complex) – Visibilities.
model (Array (n_bl, n_freq, n_time) complex) – Fringe model on the same grid.
weights (Array) – Anything broadcastable onto (n_bl, n_freq, n_time) – a per-baseline (n_bl, 1, 1) is the common case.

Returns:

(z, n1, n2), each (n_freq, n_time); z complex, the others real.

Return type:

(Array, Array, Array)

tabascal.rfi_estimate.near_field_baseline_paths(record, ants_itrf: NDArray, times_jd: NDArray, phase_centre: dict, a1: NDArray, a2: NDArray, n_fine: int, delta_t: float, taus_s=0.0) → NDArray[source]

Per-baseline near-field path difference on the fine grid, in metres.

The quantity the whole search is built on. Per antenna the path is the spherical one, |x_sat - x_a| – never the plane-wave projection – plus the phase-tracking term w_a the visibilities are already rotated by, exactly as tabascal.interferometry.get_rfi_phase_numpy() assembles it. The baseline quantity is the difference path_p - path_q, and the fringe model is exp(-2 pi i (path_p - path_q) / lam) averaged over the fine axis (near_field_fringe_model()).

taus_s moves the satellite only: the orbit is propagated to t + offset + tau while the antennas, the sidereal angle and w stay at t + offset. That is what an along-track TLE error is. Shifting the whole geometry instead turns the Earth under the phase tracking as well, which moves the path difference by a tenth of a wavelength at metre wavelengths – a different model, and the wrong one.

Everything here is numpy/f64 on the host, and stays f64 whatever --x64 says: an absolute path is hundreds of kilometres, which f32 resolves to tens of metres, and the model needs a small fraction of a wavelength. The difference handed to the jitted core is at most the array’s diameter, which f32 does hold to well under a wavelength, so the cast happens there and not before.

The antennas’ fine-grid positions are computed once and shared by every tau; only the orbit is re-propagated, in a single vectorised call over the flattened (n_tau, n_time, n_fine) grid.

Sizing: the result is n_tau * n_bl * n_time * n_fine float64. For the MWA case (33 offsets, ~1000 coherent baselines, 27 frames, 40 sub-steps) that is ~290 MB, which is why the coherence cut is applied to the baseline list before the paths are built rather than after.

Parameters:

record (dict) – One orbit record – TLE or OMM, as resolved by tabascal.orbit.
ants_itrf (Array (n_ant, 3)) – Antenna positions in the ITRF (ECEF) frame, in metres.
times_jd (Array (n_time,)) – Integration centres, as UTC Julian dates.
phase_centre (dict) – {"ra": <deg>, "dec": <deg>} phase centre of the visibilities.
a1 (Array (n_bl,)) – Antenna indices of each baseline.
a2 (Array (n_bl,)) – Antenna indices of each baseline.
n_fine (int) – Sub-steps per integration (see fine_time_offsets()).
delta_t (float) – Integration time in seconds.
taus_s (float or Array (n_tau,), default 0.0) – Along-track offsets to evaluate. A scalar still returns a grid of one, so the core is written once, for a grid.

Returns:

Path difference path_p - path_q in metres.

Return type:

Array (n_tau, n_bl, n_time, n_fine) float64

The per-integration fringe model: the template, averaged over its dump.

mean_fine exp(-2 pi i path / lam). Each sub-step is a unit-modulus steering vector, so the average can only shrink: a baseline whose fringe turns a fraction of a cycle inside one integration keeps its modulus, and one that sweeps tens of cycles averages away to nothing. That loss is real – it is in the data too – and modelling it is what lets the normalised statistic of coherence_scores() stay a correlation.

Pure jax.numpy and undecorated: the drivers own the jit, and the batched search of #191 needs it traceable inside its own. Its working precision is the session’s, which is why the paths are built in f64 on the host and only their (short) baseline differences arrive here.

The intermediate is (..., n_bl, n_freq, n_time, n_fine) complex, which is the largest array in the scan; tau_scan() holds one offset’s worth of it at a time.

Parameters:

paths (Array (..., n_bl, n_time, n_fine)) – Baseline path differences in metres (see near_field_baseline_paths()).
freqs (Array (n_freq,)) – Channel frequencies in Hz.

Returns:

Fringe model, of modulus at most 1.

Return type:

Array (…, n_bl, n_freq, n_time) complex

tabascal.rfi_estimate.offset_fit_summary(norad_id, fit: dict, threshold_sigma: float = 5.0) → str[source]

One line per satellite: what was measured, and whether it counts.

DETECTED in upper case and not detected in lower, so a log can be grepped for the one without matching the other.

tabascal.rfi_estimate.plot_candidate_ranking(search: dict, selection: dict, save_path: str, title: str | None = None) → str[source]

The whole field, ranked: the chart a named satellite is judged against.

One bar per candidate in rank order, detections in red, the candidate median as a dotted line and the winner annotated with its channel, its score and the runner-up’s. A single NORAD ID is not evidence; a winner standing clear of a median drawn from every satellite that was up is. On the MWA Cen A case that was 0.0995 against a runner-up of 0.0523 and a median of 0.0446.

Parameters:

search (dict) – A result from search_candidates().
selection (dict) – The matching select_detections() result; its rows are highlighted.
save_path (str) – Output PNG path.
title (str, optional) – Prefix for the figure title.

Returns:

save_path.

Return type:

str

tabascal.rfi_estimate.plot_offset_diagnostics(fit: dict, save_path: str, title: str | None = None) → str[source]

The three panels an along-track detection is judged by eye on.

Top, the frame-by-channel correlation |r| at the best offset, with the satellite’s elevation over it: a real detection is a band that lights up while the satellite is up and goes out when it sets, on one channel, not a scatter of warm cells. Bottom left, the per-channel z2 at that offset against the decohered null’s mean and spread, which is the comparison the significance is. Bottom right, the scan itself for the best channel, with the other channels’ range shaded behind it: a single peak that decays back into that band is what makes the best cell a measurement rather than the largest of many sidelobes.

Parameters:

fit (dict) – A result from fit_time_offset().
save_path (str) – Output PNG path.
title (str, optional) – Prefix for the figure title, usually the NORAD ID.

Returns:

save_path.

Return type:

str

tabascal.rfi_estimate.plot_z_spectrograms(result: dict, save_path: str, z_crit: float = 3.0, vmax: float | None = None) → str[source]

Per-source spectrogram of the z statistic (residual / floor).

One panel per source, time on x, frequency (MHz) on y, colour = signed z = Re(S_hat) / error on a diverging scale (blue = over-subtracted, red = under-subtracted residual). Single-channel data degrades to a z-vs-time line plot with the +/- z_crit band shaded.

Returns:: save_path.
Return type:: str

tabascal.rfi_estimate.rayleigh_threshold(z_crit: float) → float[source]

The |S_hat|/error cut enclosing the same probability as |z| <= z_crit.

Under the complex Gaussian null the real and imaginary parts of S_hat are independent N(0, error^2), so |S_hat|/error is Rayleigh(1) with P(R <= c) = 1 - exp(-c^2/2). Matching that to the two-sided normal probability leaves c = sqrt(-2 ln(erfc(z / sqrt 2))) – 3.44 for the usual 3 sigma – so the two coverages are read on the same scale rather than against thresholds that mean different things.

The tail is taken from erfc rather than as 1 - erf: that subtraction cancels to exactly zero once erf rounds to 1, somewhere past 6 sigma, and the threshold came back infinite – which marks every cell as consistent with noise and reports a coverage of 100% for any data at all. erfc is the two-sided tail, computed without the cancellation.

tabascal.rfi_estimate.rfi_phase_from_positions(rfi_xyz: NDArray, ants_itrf: NDArray, times_jd: NDArray, phase_centre: dict, freqs: NDArray) → NDArray[source]

Per-antenna geometric phase for RFI sources at known ECI positions.

Numpy/f64 host-side computation mirroring tabascal.components.trajectory.FixedOrbit._compute_rfi_phase(), so the estimate is matched to the phase the forward model itself builds.

Parameters:

rfi_xyz (Array (n_src, n_time, 3)) – Source positions over time in the ECI (GCRF) frame, in metres.
ants_itrf (Array (n_ant, 3)) – Antenna positions in the ITRF (ECEF) frame, in metres.
times_jd (Array (n_time,)) – Observation times in Julian date.
phase_centre (dict) – {"ra": <deg>, "dec": <deg>} phase centre of the visibilities.
freqs (Array (n_freq,)) – Channel frequencies in Hz.

Returns:

Geometric phase at each antenna for each source.

Return type:

Array (n_src, n_ant, n_freq, n_time)

tabascal.rfi_estimate.rfi_phase_from_records(orbit_records: list, ants_itrf: NDArray, times_jd: NDArray, phase_centre: dict, freqs: NDArray, time_offsets_s=None) → NDArray[source]

Per-antenna geometric phase for satellites given their orbit records.

Propagates each record – TLE or OMM, as resolved by tabascal.orbit – over times_jd and defers to rfi_phase_from_positions().

Parameters:

orbit_records (sequence of dict, length n_src) – Orbit records, in the order the curves are wanted in.
ants_itrf – See rfi_phase_from_positions().
times_jd – See rfi_phase_from_positions().
phase_centre – See rfi_phase_from_positions().
freqs – See rfi_phase_from_positions().
time_offsets_s (float or sequence of float, optional) – Along-track offset per source, in seconds, as measured by fit_time_offset(): source i’s orbit is evaluated at times_jd + tau_i. A scalar applies to every source. It moves the satellite only – the antennas, the sidereal angle and the phase tracking stay at times_jd, because tau is an error in the orbit, not in the observation’s clock. Re-propagating the whole geometry at t + tau would rotate the Earth under the fringe tracking as well and give a different phase. None (the default) is the behaviour every existing caller has, bit for bit.

Return type:

Array (n_src, n_ant, n_freq, n_time)

tabascal.rfi_estimate.satellite_range_and_speed(record, ants_itrf: NDArray, times_jd: NDArray)[source]

Slant range and line-of-sight-crossing speed, at mid-observation.

The two numbers coherent_baseline_mask() needs, measured from the array centre by finite difference over a second of the pass. Both are taken in the Earth-fixed frame, which is the one the array sits still in: what fringes a baseline is the rate at which the direction to the satellite sweeps across the antennas, and a geostationary emitter – which hangs motionless over the array – moves at three kilometres a second in the inertial frame.

That emitter does not fringe nothing, though. On phase-tracked visibilities the w term keeps turning at the sidereal rate whatever the satellite does, worth a fringe of order Omega_e b / lam – about 0.03 Hz on 600 m at 1.7 m, negligible beside a LEO’s but not zero. This function returns the satellite’s own transverse speed; the term the phase tracking adds is bounded by Omega_e * range_m and fit_time_offset() adds it where the fringe-rate ceiling is sized.

The Earth-fixed direction is recovered by turning the inertial separation back by the sidereal angle. That leaves the precession-nutation rotation in it, which is constant to a part in 1e12 over the second differenced here and so cannot affect a rate; against the full frame transform the speed agrees to a part in 1e4.

It is not the orbital speed. Near the horizon most of a LEO satellite’s motion is along the line of sight and does not fringe, and its range is several times the overhead one – both of which lengthen the coherence ceiling rather than shorten it.

Parameters:

record (dict) – One orbit record.
ants_itrf (Array (n_ant, 3)) – Antenna positions in the ITRF frame, in metres; the mean is the site.
times_jd (Array (n_time,)) – Observation times as UTC Julian dates. The middle one is used: the range varies by tens of percent over a pass, far less than the uncertainty on the orbit error the cut is set by.

Returns:

Slant range in metres, and speed across the line of sight in m/s.

Return type:

(float, float)

tabascal.rfi_estimate.save_light_curves_npz(path: str, result: dict) → None[source]

Save a light-curve result as the rfi.est interchange format.

The four names tabascal.components.rfi_signal.read_light_curves() requires – light_curves (n_src, n_time, n_freq), norad_ids, times (UTC MJD) and freqs (Hz) – so the output of a tabascal light-curve run can be pointed at with rfi.est unchanged.

times is UTC and not the MS’s TIME column as declared. The format states one scale so a curve stays interpretable away from the MS it was measured on, and the reader samples it on the same one: on a TAI-declared MS the declared numbers are 37 s from the instants they name, which would seed a later run with a satellite that brightens at the wrong times. That scale is stamped into the file as time_scale, so a reader never has to assume it – and so an untagged file, which pre-dates the stamp and may have been written on a declared scale, can be told apart and warned about.

light_curves is the magnitude |S_hat|, an apparent flux in Jy: the reader casts to float64, which would silently discard the imaginary part of a complex array. The native complex estimate is kept alongside it under light_curves_complex, together with the noise floor (error), the z statistic and the in-view mask. Readers of the format ignore the extras.

Where an along-track offset was fitted (attach_offset_fits()) its answers are written too – tau_best above all, without which a later run cannot reproduce the trajectory these curves were measured on – with r_best swapped into the file’s (n_src, n_time, n_freq) orientation like the curves themselves.

Parameters:

path (str) – Output .npz path.
result (dict) – A dict from one of the driver functions.

tabascal.rfi_estimate.search_candidates(vis: NDArray, candidates: list, ants_itrf: NDArray, times_jd: NDArray, phase_centre: dict, freqs: NDArray, a1: NDArray, a2: NDArray, int_time: float, noise=None, flags: NDArray | None = None, taus_s=None, n_fine: int = 40, sigma_transverse_m: float = 300.0, soft_weights: bool = False, exclude_autos: bool = True, batch_size: int = 8, max_mem_gb: float = 4.0, n_null: int = 200, null_jitter_m: float = 50.0, n_null_candidates: int = 5, seed: int = 0, progress=None) → dict[source]

Score every candidate satellite against the visibilities, and rank them.

Stage two of GitHub #191: the tau scan of fit_time_offset(), run over a whole snapshot. It is the same statistic – tau_scan() under jax.vmap, one jitted program per batch – because a search that re-derived the filter would be free to disagree with the single-satellite fit about what a detection is. The score per candidate is the largest z2 over the offset grid and the channels: the emission is narrowband, so the maximum over channels is the right statistic and the winning channel is itself a deliverable, telling the user which channels to fit.

One baseline set, each candidate’s own cut. vmap needs static shapes, so the search sums over one baseline list: the union of the hard coherent sets of the candidates above the geometric horizon. A union rather than the farthest candidate’s set, because the sets are not nested by range – the fringe-rate ceiling depends on each candidate’s own transverse speed, so a nearer, slower satellite can steer a baseline a farther, faster one cannot. Each candidate then applies its own coherence as a per-baseline weight inside the statistic, so another candidate’s excess baselines enter at exactly zero and it is scored over precisely the baselines it could steer – the same numbers a search of that candidate alone would give. Only the above-horizon candidates size the union: a satellite 13 000 km away, through the Earth, tolerates kilometres of baseline and would otherwise readmit the long ones for everybody. With soft_weights the support is still the hard union and the taper weights the baselines inside it (_coherence_weights()).

The horizon lives inside the statistic. Each candidate’s in-view mask from enumerate_candidates() is passed to tau_scan() as frame_mask rather than slicing the arrays, so a satellite that rises or sets mid-observation contributes only its own frames while the batch stays rectangular. Masking and slicing give the same z2, so the search and the single-satellite fit report the same detection for the same pass.

Sizing. Two arrays dominate, per candidate of a batch: the fringe model, (n_bl, n_freq, n_time, n_fine) complex, held one offset at a time, and the paths for the whole grid, (n_tau, n_bl, n_time, n_fine) float64 on the host. max_mem_gb is a budget for their sum, and the batch actually run is the smaller of batch_size and what that budget affords – reported back as batch_size, since it is what the sweep did rather than what it was asked for. It is what sets the batch at real scale: on the MWA case the coherent union reaches 7704 of the array’s 9180 baselines once candidates come near the horizon (b_coh_max 2880 m), one candidate over 24 channels is then some 2.1 GB, and a batch of eight would ask for 17 GB. Not counted are the per-candidate weights, at whatever shape the noise and the flags resolve to – per baseline for a scalar sigma, the full (n_bl, n_freq, n_time) once anything is flagged – so this is a sizing heuristic and not a cap. A ragged last batch is padded by repeating its last candidate so every batch has one shape and the kernel compiles once; the padding rows are dropped.

The null is drawn for the top ``n_null_candidates`` only. Two hundred extra scans per satellite over a whole constellation is the search twice over, spent on candidates nothing will be reported for. The rest are scored and ranked but carry no significance. The shortlist is taken on raw z2, as the issue specifies – and z2 is a sum over in-view frames, so a short pass ranks low against a full one at the same per-frame correlation. That is a caveat on the shortlist, not a correction to make: comparing partial passes against full ones is what the score is for.

Parameters:

vis (Array (n_bl, n_freq, n_time) complex) – Visibilities to search.
candidates (list of dict) – Screened candidates from enumerate_candidates().
ants_itrf (Array (n_ant, 3)) – Antenna positions in the ITRF frame, in metres.
times_jd (Array (n_time,)) – Integration centres as UTC Julian dates.
phase_centre (dict) – {"ra": <deg>, "dec": <deg>}.
freqs (Array (n_freq,)) – Channel frequencies in Hz.
a1 (Array (n_bl,)) – Antenna indices of each baseline.
a2 (Array (n_bl,)) – Antenna indices of each baseline.
int_time (float) – Integration (dump) time in seconds.
noise – As fit_time_offset(), and meaning the same thing.
flags – As fit_time_offset(), and meaning the same thing.
taus_s – As fit_time_offset(), and meaning the same thing.
n_fine – As fit_time_offset(), and meaning the same thing.
sigma_transverse_m – As fit_time_offset(), and meaning the same thing.
soft_weights – As fit_time_offset(), and meaning the same thing.
exclude_autos – As fit_time_offset(), and meaning the same thing.
n_null – The decohered null, as fit_time_offset() draws it.
null_jitter_m – The decohered null, as fit_time_offset() draws it.
seed – The decohered null, as fit_time_offset() draws it.
batch_size (int, default 8) – Most candidates to score per jitted call. An efficiency knob and nothing else: the answers do not depend on it.
max_mem_gb (float or None, default 4.0) – Memory budget in gigabytes, which lowers the batch when batch_size candidates would not fit. None is no budget, for a caller who has measured their own. See Sizing above for what it does and does not cover.
n_null_candidates (int, default 5) – How many of the ranked candidates get a decohered null, and so a significance.
progress (callable, optional) – Called (done, total) after each batch. A search over a constellation runs for minutes, and a command that says nothing for all of them cannot be told from a hung one.

Returns:

table, a list of per-candidate dicts in ranked order carrying rank, norad_id, name, max_elevation, range_m, z2_best, tau_best, best_chan, best_freq, r_max, significance, null_mean, null_std and n_frames; the same columns as arrays under norad_ids, z2_best, tau_best, best_chan and significance; the scan curves z2_tau (n_cand, n_tau, n_freq) on tau_grid; frames (n_cand, n_time) bool; n_bl_used and b_coh_max for the shared baseline set; batch_size, the batch the sweep actually ran at; median_z2; fits, one dict per candidate shaped like fit_time_offset()’s, so plot_offset_diagnostics() and attach_offset_fits() take them unchanged; and candidates, the screened candidates themselves in the same ranked order, which is where the orbit record of a named satellite comes from.

A row’s range_m is the candidate’s closest approach, for reporting; the geometry the cut was sized from is the mid-window pair in fits[i]["range_m"] and fits[i]["v_perp_m_s"].

The decision is not among them; see select_detections().

Return type:

dict

tabascal.rfi_estimate.select_detections(search: dict, threshold_sigma: float = 5.0, runner_up_ratio: float = 1.5) → dict[source]

Decide which candidates are detections, and what to warn about.

Ranking is not deciding. A detection is a candidate whose significance was measured and clears threshold_sigma; one outside the null’s shortlist carries nan, and nan >= 5 is false while nan < 5 is false too – only one of those readings is safe to rely on. Multiple detections are simply all returned, in rank order: the fit accepts several satellites already.

Two warnings, both from the issue and both about a ranking that cannot be read at face value:

a close runner-up – the second candidate within runner_up_ratio of the first. Satellites in the same train partially match each other’s fringes, so a winner that is not clear of the field is a result to look at twice rather than a satellite to name. It is a statement about the ranking and does not wait for a detection: two candidates level at the noise floor say the search could not separate them, which is worth knowing even when neither is named.
a scan edge – a detected candidate whose best offset is the first or last point of the grid. The peak may be off the grid, so the offset is at least that large and the number reported is a floor; widening --tau-max is the fix. Only for a detection: an undetected candidate’s best offset is wherever noise happened to peak, and that it did so at the end of the grid says nothing about anything.

Parameters:

search (dict) – A result from search_candidates().
threshold_sigma (float, default 5.0) – Significance above the decohered null at which a candidate counts. It carries no trials factor; see fit_time_offset().
runner_up_ratio (float, default 1.5) – Warn when z2[1] >= z2[0] / runner_up_ratio.

Returns:

detected, the qualifying rows of search["table"] in rank order, and warnings, a list of sentences.

Return type:

dict

tabascal.rfi_estimate.select_sources(result: dict, keep) → dict[source]

A light-curve result restricted to some of its sources.

Threshold-gated saving is a selection along the source axis, and it has to be made in one place: a per-source array left behind would mislabel every curve after the first one dropped.

Parameters:

result (dict) – A light-curve result, with or without an offset fit attached.
keep (Array (n_src,) bool, or Array of int) – A mask or an index array over the sources, read the same way either way.

Returns:

A copy carrying only the selected sources. The coordinates – the frequencies, the times, the offset grid – are not per source and come through untouched.

Return type:

dict

tabascal.rfi_estimate.shift_orbit_record_epoch(record, tau_s: float) → dict[source]

A copy of record whose epoch is moved by -tau_s seconds.

The zero-code-change way to use a fitted offset. tau is measured as the time the satellite is evaluated at, so a positive one means the elements are late – the satellite is where they say it will be tau seconds later – and the record has to become that trajectory: propagating the shifted elements at t reproduces the original at t + tau. Hence the minus.

Written into a directory by write_shifted_orbits(), the result is picked up by --extra-orbit-dir and nothing else in tabascal has to know about the search at all.

A TLE’s line-1 epoch field (columns 19-32, YYDDD.DDDDDDDD) is rewritten and the modulo-10 checksum recomputed, since a rewritten field invalidates it and the parser rejects a bad one – as it should, that being how a single-character corruption is caught. Nothing else on either line moves. The field quantises to 1e-8 days, 0.86 ms, which is about 7 m along a LEO track; an OMM has no fixed-width field and keeps the epoch to the microsecond, so it is the format to shift where there is a choice. Either way the record’s own EPOCH column is moved with it, since a record whose column disagreed with its own elements would be read one way by the age policy and another by the propagator.

Parameters:

record (dict) – A TLE or OMM orbit record.
tau_s (float) – The measured along-track offset in seconds. Zero returns the record’s own epoch unchanged, not a re-encoding of it.

Returns:

A copy; the record handed in is not edited under the caller.

Return type:

dict

Score every along-track offset on the grid: the scan itself.

One program for the whole grid. The offset axis is walked with jax.lax.map, not a Python loop over jitted kernels: a loop re-enters the compiler per step and gives up the point of the design, which is one compilation and then a device-resident sweep. lax.map rather than a vmap over the grid for the same reason the paths are looped on the host – the per-offset (n_bl, n_freq, n_time, n_fine) model is the biggest array in the calculation and only one of them need exist at a time.

The function is pure and of fixed-shape arrays, so it vmaps over a leading candidate axis of paths; that is the contract the multi-satellite search of #191 is built on.

Parameters:

vis (Array (n_bl, n_freq, n_time) complex) – Visibilities, already zeroed wherever their weight is.
weights (Array) – Broadcastable onto vis; see matched_filter_sums().
paths (Array (n_tau, n_bl, n_time, n_fine)) – Baseline path differences per offset (near_field_baseline_paths()).
freqs (Array (n_freq,)) – Channel frequencies in Hz.
frame_mask (Array (n_time,), optional) – See coherence_scores().
ant_offsets (Array (n_ant,), optional) – Per-antenna path offsets in metres, added as paths + off[a1] - off[a2]. This is the decohered null’s hook (decohered_null()); a1 and a2 are then required.
a1 (Array (n_bl,), optional) – Antenna indices, needed only with ant_offsets.
a2 (Array (n_bl,), optional) – Antenna indices, needed only with ant_offsets.

Returns:

z2 (n_tau, n_freq) and r (n_tau, n_freq, n_time). r keeps every frame, masked or not: it says what each frame did, and the mask decides only which are combined.

Return type:

dict

Longest baseline an orbit known to sigma metres can still steer.

A transverse position error delta at slant range r moves the satellite’s apparent direction by delta / r, which costs 2 pi (b / lam) (delta / r) of template phase on a baseline of length b. Holding that to one radian gives:

b_tle = lam r / (2 pi delta)

Only the transverse error enters. The along-track error, which is the larger part of a TLE’s, is very nearly a pure time offset and is absorbed by the tau search instead of by this cut.

The expression is symmetric in b and delta, so it inverts itself: the same call evaluated at sigma_transverse_m = b is the largest orbit error the baseline b tolerates.

Parameters:

freq (float or Array) – Observing frequency in Hz.
range_m (float or Array) – Slant range to the satellite in metres.
sigma_transverse_m (float or Array) – Transverse (across the line of sight) TLE position error in metres. Taken as a magnitude: a caller holding a signed component – an offset measured along some axis – means its size, and a negative ceiling would be met by no baseline at all. Zero is a perfect orbit and returns an infinite coherence length.

Returns:

Coherence length in metres, at the broadcast shape of the inputs.

Return type:

Array

tabascal.rfi_estimate.write_config_fragment(path: str, selection: dict, shifted_orbit_dir=None) → str[source]

Write the detections as a tabascal satellites section.

The deliverable the issue asks for: the norad_ids list a run needs, produced from the data rather than known in advance, ready to merge into a config. Written even when nothing was detected – an empty list with a line saying why is an artifact to point at, where a missing file leaves the reader to guess the run failed.

A bare list of IDs is not auditable, so above it, per satellite and as YAML comments, is what it was detected on: the offset its curves must be extracted at, its score, its significance, and the channel – which is what tells the user which channels to fit. Any warnings from select_detections() are written under them.

With shifted_orbit_dir the section also points at the epoch-shifted records (write_shifted_orbits()) with no age ceiling, which is the whole point of writing them: a later run reproduces the trajectories the search measured, whatever SatChecker serves by then.

Returns:: path.
Return type:: str

tabascal.rfi_estimate.write_search_results(path: str, search: dict, selection: dict, threshold_sigma: float) → str[source]

Save the ranking table, detections or not.

“Nothing above the threshold” is a result about every satellite that was up, and the evidence for it is the table – so this is written either way, and detected records which rows the selection named against the threshold it was named by.

Returns:: path.
Return type:: str

tabascal.rfi_estimate.write_shifted_orbits(directory: str, norad_ids, records: list, taus_s, filename: str = 'shifted_orbits.json') → str[source]

Write epoch-shifted orbit records where a later run can pick them up.

One file in extra_orbit_dir format (tabascal.orbit.save_orbits_for_reuse()), so tabascal run --extra-orbit-dir <directory> reproduces the trajectories the search measured – with the default unlimited age ceiling, and independently of what SatChecker serves by then.

Parameters:

directory (str) – Directory to write into; created if it does not exist.
norad_ids (sequence of int) – Catalogue IDs, aligned with records.
records (sequence of dict) – The orbit records to shift.
taus_s (float or sequence of float) – Offsets in seconds, one per record or one for all of them.
filename (str, default "shifted_orbits.json") – Name of the file inside directory.

Returns:

The path written.

Return type:

str