2.1. Software spectrometer (SWspec)

SWspec calculates the spectral power integrated over a selected time of broadband radio signals. It can be utilised to estimate the fast Doppler variations of spacecraft carrier signals as well as single-dish spectral line observations.

2.1.1. Astronomical software spectrometer

SWspec performs the initial detection of the spacecraft carrier signal, subcarriers and ranging tones (if present) at any desired spectral resolution. It is based on the software spectrometer developed for processing observations of the ESA’s Huygens Titan Probe and the SMART-1 Lunar mission Pogrebenko (Reference Pogrebenko2006). SWspec is a versatile and highly configurable spectrometer, written in C++ with a Python wrapper, that handles radio astronomical data recorded by VLBI data acquisition systems (DAS).

SWspec is compatible with most current and legacy VLBI data formats. It utilises a common external library (mark5access)Footnote d, to read the quantised time domain signal recorded by the VLBI DAS and extracts any selected frequency channel. It then performs an accurate window-overlap-add (WOLA) discrete Fourier transform (DFT) per iteration and calculates the time-integrated spectrum. The user can set a wide range of spectral resolution (from sub-Hz to tens of kHz) and integration times (from seconds to minutes).

Both SWspec and SCtracker use window (apodisation) functions to improve the reliability of detections, smooth the spectrum, and reduce the aliasing and cyclicity effects. The choice of apodisation function is a compromise between spectral leakage effects and spectral resolution. The software use a large number of spectral points ( ${\ge}\!10^6$ points) to compensate for the spectral smearing. The available window functions (w[n]) in the software and their corresponding main lobe full width at half maximum (FWHM) are listed in Table 1. Note that the cosine squared and Hann windows are mathematically the same, but they are implemented differently in the software.

Table 1. Available window functions (w[n]) in SWspec and SCtracker, their numerical basis, and main lobe full width at half maximum (FWHM) in spectral bins. Where N is the length of the data vector in samples and n is the sample number (starting with 0).

^*The cosine function becomes a sinus because of our definition of n, which starts with 0.

Each of the window functions provides advantages and disadvantages. For WOLA applications, it is desirable that the window is (a) as close as possible to amplitude 0 at the start and end of an interval and close to 1 at its middle and (b) the sum of the windows shifted by half the length of the data vector be close to flat 1. For example, the cosine window does not obey the rule (b) causing modulation of the output power. We selected the cosine squared (or the Hann) window function to be used by default to calculating the DFT of spacecraft signals after extensive testing.

SWspec is configured by a control file that initialises all the input parameters described in Table 2. A number of these parameters depend on the antenna setup used during the observations: the bandwidth of each critically sampled frequency channel (BandwidthHz), the data format (SourceFormat), and the number of frequency channels (SourceChannels). Other parameters, such as the number of spectral points (FFTpoints), the integration time (FFTIntegrationTime), and the FFT overlap (FFToverlapFactor), are used to adjust the spectral resolution and the accuracy of the carrier detection. These depend mainly on the radial acceleration of the spacecraft. The frequency resolution (dF) is equivalent to the sampling rate (SR) divided by the number of spectral points (dF $\,=\,$ SR $\,\cdot\,$ FFTpoints-1) in Hz, where SR is twice the bandwidth of the frequency channel (SR $\,=\,2\,\cdot\,$ BW). The time resolution (dT) is equal to the integration time of each spectrum in seconds.

Table 2. List of all the parameters available in the control file to initialise and run the software spectrometer.

The number of spectral points and the integration time are the most critical parameters for reliable detection of the spacecraft signal. The selection of the appropriate parameters depends on the receiving antenna system equivalent flux density (SEFD), the spacecraft transmission power, propagation losses, and the relative motion between the ground station and spacecraft. These parameters can differ significantly if the target is an orbiter, a lander or an interplanetary spacecraft. More details on the selection of these parameters are provided in Section 4.1. The spectral output data are written to disk in 32-bit binary format. These spectra are used for modelling the Doppler shift between the ground-station and spacecraft.

2.1.2. Phase polynomial fit

A series of Python (Van Rossum & Drake Jr 1995) libraries have been developed to visualise the spectral data and extract intrinsic measurements of the radio signal(s). These scripts are available in the pysctrack Git repository. CalculateCpp.py allows the user to detect the frequency of the main carrier signal, estimate the rate of the Doppler shift, and build the phase polynomial fit.

CalculateCpp.py searches for the peak of the spacecraft carrier signal and estimates the topocentric frequency detections along the entire frequency band and scan length. When SWspec is used with high spectral resolution the spectrum of the carrier line appears shifted by the Doppler shift. Typical processing configuration consists of $3.2\cdot10^{6}$ FFTpoints for a 16 MHz frequency channel. To compensate for the spectral power being spread out across multiple bins (smearing), the actual location of the carrier tone is estimated via a spectral centroid over five frequency bins around the peak. Thus, the spectral line position is achieved with sub-bin precision.

The time-series of topocentric frequency detections f(t) of the spacecraft carrier signal are used to construct a frequency polynomial regression fit F(t) of m-th using the weighted least mean square method.

(1) \begin{equation} F(t) = C_{\rm fs}(0) + C_{\rm fs}(1) \cdot t + ... + C_{\rm fs}(m-1) \cdot t^{m-1}, \end{equation}

where $C_{\rm fs}$ are the frequency polynomial coefficients. The degree of polynomial (m) ranges between 3 and 7 depending on the length of the scan, the spectral resolution and the rate change of the Doppler shift. We use the signal-to-noise ratio (SNR) as weights per sample.

SNR for each sample can be determined as a ratio of detected spectral line power to the rms noise power in the near-by region of spectrum, defined as:

(2) \begin{equation} {\rm SNR} = \frac{P_{\rm sc}}{k \cdot T_{\rm sys} \cdot {\rm dF}}, \end{equation}

where k is the Boltzmann constant, $T_{\rm sys}$ the system temperature of the antenna measured in [K], $P_{\rm sc}$ the spacecraft tone spectral power in [W], and dF the frequency resolution in [Hz] used in SWspec. In the case when it is not possible to compute the rms noise level, weights can be set to 1. See the pysctrack user guide in the GitLab repositoryFootnote e for more information.

The phase polynomial

(3) \begin{equation} P(t) = C_{\rm ps}(0) + C_{\rm ps}(1) \cdot t + ... + C_{\rm ps}(m-1) \cdot t^{m-1}, \end{equation}

is constructed by the phase polynomial coefficients ( $C_{\rm ps}$ ) calculated from the frequency polynomial by integrating its coefficients (for m $\,>\,0$ ) according to:

(4) \begin{equation} C_{\rm ps}(m) = \frac{2 \pi}{m} \cdot C_{\rm fs}(m-1), \end{equation}

where the initial phase (integration constant) is set to 0, the second is the phase rate [rad/s], the third is the acceleration [rad/s $^2$ ], and so on; and $C_{\rm fs}$ are the frequency polynomial coefficients per second being the first component the detected tone frequency [Hz], rate [Hz/s], acceleration [Hz/s $^2$ ]. For computational convenience the phase coefficients are transformed from phase per second to phase per sample ( $C_{\rm pp}$ ) as:

(5) \begin{equation} C_{\rm pp}(m) = C_{\rm ps}(m) \cdot SR^{\,-m}. \end{equation}

where SR is the sampling rate in samples per second.

The procedure CalculateCpp.py performs the tasks described above. It outputs the SNR, the topocentric frequency detections, the stochastic Doppler noise, and stores them in the Fdets file. A detailed description of the output files is presented in Section 2.4. The script also generates the $C_{\rm pp}$ file that contains the phase polynomial coefficients with respect to the initial base band frequency channel in terms of radians per sample and the $C_{\rm fs}$ file with the frequency polynomial coefficients in terms of Hz per second.

The $C_{\rm pp}$ and $C_{\rm fs}$ coefficients per are used as input parameters for SCtracker.

Table 3. Main input files and settings required by SCtracker and the output products after execution. All files are stored in the same directory in which SCtracker runs.

2.2. Multi-tone spacecraft tracker (SCtracker)

SCtracker is the core package for tracking the carrier radio signals transmitted by spacecraft. This program is an evolution of the software developed for tracking of the Huygens Probe Pogrebenko et al. (Reference Pogrebenko, Gurvits, Campbell, Avruch, Lebreton, van’t Klooster and Wilson2004) and is the prototype of the current European VLBI Network (EVN) Super FX correlator (SFXC) Keimpema et al. (Reference Keimpema2015). SCtracker is written in C++ with a light graphical user interface in Python. Table 3 lists the software’s main input and output files. All the configurable parameters in the SCtracker control file are described in the Table A1 at the Appendix.

SCtracker utilises the estimated phase polynomial coefficients ( $C_{\rm ps}$ ) to conduct an initial phase-stop by doing a phase-rotation of the entire recorded bandwidth. The desired carrier frequency tone is then selected to perform signal filtering, extract a narrow band around the tone with the selected bandwidth, and detect the relative phase of the tone. SCtracker is capable of simultaneously tracking the spacecraft carrier signal, its subcarriers and all the ranging tones. This feature is described in more detail in Section 3.1.1.

SCtracker transforms the initial baseband signal S(t) to a complex signal S ^′(t) with the nonlinear part of the phase evolution corrected as expressed in:

(6) \begin{equation} S'(t) = S(t) \cdot e^{-i \phi_{\rm corr(t)}}.\end{equation}

Equation (6) is the analytical form of the double-precision polynomial transformation which is applied to the baseband sample sequence (x[n]), where the applied phase correction ( $\phi_{\rm corr(t)}$ ) corresponds to the polynomial fit coefficients ( $C_{\rm pp}$ ) calculated previously.

(7) \begin{equation} x'[n]=x[n]\exp \left( \pm i \sum\limits_{k=1}^{m-1}C_{\rm pp}(k)\cdot T_{n}^{\,k}\right), \end{equation}

where x ^′[n] are the new samples, x[n] are the original raw samples, $T_{n}$ is the time relative to the beginning of the data segment, and $\pm$ signs take into account that the frequency channel can be upper or lower side band. The calculated time-integrated window-overlapped spectra of the phase stopped baseband signal are saved to disk. These use the same spectral representation as produced with SWspec, with the frequency shift corrected as specified by the polynomial approximation according to Equation (7). The spectral output of SCtracker shows the spacecraft carrier signal as a single spectral line along the entire spectral time series.

A narrow band around each tone is stored for further signal processing at a higher resolution. The user can set the bandwidth of the output tone in the initial control file. Typically, a band of 2 kHz is used. The narrow bands are extracted directly from the stopped baseband signal around each of the frequency tones specified in the tones file. Each tone is associated with a number, 0 being the carrier line, and 1 onwards the subsequent ranging tone(s). They are sorted based on their proximity to the spacecraft carrier line.

The current implementation of the software allows an unlimited number of extracted narrow bands, with an arbitrary distribution in the input band. These bands (with the tones in the centre) are filtered out into continuous complex time-domain signals using a $2^{\rm nd}\,$ order WOLA DFT-based algorithm of the Hilbert transform approximation. The extracted signals in temporal domain are written to output files in binary format for further post-processing. All complex data are stored with floating-point precision to avoid losing resolution between the different iterations. The initial timestamp of the processed scan, based on the VLBI data header, and the cut-off frequencies of the narrow band(s) are written to auxiliary output files.

2.3. Digital phase lock loop (dPLL)

dPLL provides a finer detection of the topocentric frequency measurements, an ultra-narrow band of the carrier signal, and the residual phase with respect to the station clock. It repeats with high precision the steps described in Sections 2.1 and 2.2 on the filtered narrow band spacecraft signal. The dPLL.py script that performs these operations, can be found in the pysctrack software package.

The first step in dPLL consists of calculating the new time-integrated overlapped spectra with longer integration time and finer spectral resolution. It estimates a regression line to the frequency detections and calculates a new set of phase polynomial coefficients, $C_{\rm pp}$ (r2i), using the WLMS method. A band-pass filter centred at the tone is applied, followed by a Hilbert transform approximation to discard a negative part of the spectrum. The user can select the output filter of dPLL via command line. By default, the output bandwidth (BW $_{\rm o}$ ) is set to $20\,$ Hz.

The second step in dPLL consists of shifting the spacecraft tone to a lower frequency, near the baseband cut-off. We set F $_{\rm startc}$ as the first element of C $_{\rm fs}$ . The bin position in the initial narrow band is defined as F $_{\rm max}$ . We set the frequency target (F $_{\rm target}$ ) in the middle of the output bandwidth; thus, the carrier signal can be phase rotated with the frequency shift of F $_{\rm rot}$ = (F $_{\rm max}$ - F $_{\rm target}$ ). In order to perform the phase-stop in an ultra narrow band, we apply the same concept as in Equation (6):

(8) \begin{equation} s_{\rm sc}(t) = s(t) \cdot e^{-i 2 \pi F_{\rm rot} \cdot t}.\end{equation}

The phase of the spacecraft signal is retrieved from the signal in the time domain after accounting for the phase unwrap. The phase follows an almost straight trend with a slope proportional to $2 \pi F_{\rm target}$ . Finally, the residual phase of the signal with respect to the station clock is extracted by applying a new phase polynomial fit using again a WLMS method. The order of the third polynomial fit is mostly consistent with the two previous fits.

The reconstructed phase of carrier signal is one of the outputs of dPLL.py. The phase is a function of the residual phase $\delta\Phi(t)$ , the three phase-rotations performed inside dPLL ( $F_{\rm pll}$ = $F_{\rm target}$ + $F_{\rm rot}$ + $F_{\rm startc}$ ) and the SCtracker phase-stop:

(9) \begin{equation}\begin{split} \Phi(t) = - \delta\Phi(t) + 2 \pi F_{\rm pll} \cdot t + \Phi_{\rm sctrack}(t),\end{split}\end{equation}

where $\Phi_{sctrack}(t)$ is the polynomial model obtained from the raw Doppler measurements after running SWspec and SCtracker and can be derived as

(10) \begin{equation} \Phi_{\rm sctrack}(t) = \sum_{\kappa=1}^{N_{\rm pp}} \left[C_{\rm pp,ip} \cdot \left(\frac{t}{T_{\rm span}} \right)^{\kappa} \right],\end{equation}

where $N_{\rm pp}$ is the highest polynomial order used in the SCtracker and dPLL, $C_{\rm pp}$ is the sum of the respective phase polynomial coefficients ( $C_{\rm pp1}$ and $C_{\rm pp2}$ ), and $T_{\rm span}$ is the length of the full scan.

This second step can be repeated multiple times in order to achieve a higher precision in the tone detection by narrowing the final band-pass filter. For studies of interplanetary scintillation, as described in Molera Calvés et al. (Reference Molera Calvés2014), data were filtered down to a bandwidth of $20\,$ Hz. However, when radio signals propagated through the solar corona and were affected by strong scintillation, data were analysed in a larger bandwidth (i.e. $50\,$ Hz). For phase referencing and orbit determination, the resolution of the carrier signal is improved by filtering down to a bandwidth of 5 or $1\,$ Hz.

2.4. Data products

The main results obtained with SDtracker are the topocentric frequency detections of the spacecraft carrier signal and the reconstructed and the residual phase with respect to the station reference clock. An example of the output data files generated in a session conducted on the MEX spacecraft with the Hobart-12 (Hb) radio telescope (Australia) on $2020.02.29$ is shown in Figure 2.

Figure 2. Data products obtained by running all three packages in SDtracker: SWspec, SCtracker and dPLL. The labels are based on a session observed on $2020.02.29$ at Hobart.

calculateCpp.py generates the first set of Doppler detection and the residual frequency. The file is labelled with the code r0i since it is the first iteration. The file contains four header lines with metadata including sky frequency, station, epoch, and 6 columns of values in ASCII format: Modified Julian Date (MJD), sample timestamp [s], SNR, raw spectral power of the carrier peak (normalised to 1), frequency of the carrier [Hz], and the stochastic noise of the polynomial fit [Hz]. The topocentric frequency detections for deep space probes have accuracy better than $0.2\,$ Hz at this stage, using a frequency resolution of $5\,$ Hz and $5\,$ s integration time with an input bandwidth of $8\,$ MHz. calculateCpp.py can also be used to verify the quality of the phase-stop after SCtracker. The generated topocentric frequency detections and polynomial coefficients are then labelled with the code r1i. These are rarely used since it is easier to work with narrow band rather than wide band signals.

dPLL.py generates the second iteration of the topocentric frequency detections r2i. Here, the input bandwidth is $2\,$ kHz, and the integration time varies between 1 and $10\,$ s based on the motion of the spacecraft and SNR of the carrier. For example, we use $10\,$ s for interplanetary plasma scintillation studies Molera Calvés et al. (Reference Molera Calvés2014), but $1\,$ s in the radio occultation experiments and aerobraking experiments Bocanegra-Bahamón et al. (Reference Bocanegra-Bahamón2019). These residuals frequency has generally an accuracy better than $5\,$ mHz. Finally, the topocentric frequency detections r3i are generated after extracting the residual phase of the main carrier tone and reconstructing the frequency from the phase of the signal. The difference between the r2i and r3i revisions is small ( $\leq\! 1\,$ mHz), but noticeable in experiments where ultra-high accuracy is required Litvinov et al. (Reference Litvinov2016); Litvinov et al. (Reference Litvinov2018).

The residual phases file includes two columns with the timestamp and residual phase in radians for the whole scan. This file also contains four header lines summarising the epoch, station name, MJD, initial time, averaged SNR, and residual frequency after the second polynomial fit. The total number of samples (TS) is twice the output bandwidth times the scan length (TS = $2\cdot$ BW $_{\rm o}\cdot T_{\rm span}$ ).

3.1. Software features

The SDtracker has been used over a thousand spacecraft tracking sessions to optimise its features and develop a robust software package.

3.1.2. Radio frequency interference (RFI)

RFI is a severe problem in radio astronomy and spacecraft tracking applications are also affected from it. It can be only a minor concern while working with spacecraft equipped with high-gain antennas, as the received signal is strong. However, it may become problematic in space missions at large distances from Earth or when equipped with low-gain antennas. Identifying the signal that belongs to the spacecraft is achieved with the Python scripts in pysctrack. When the spacecraft tone is strong, the scripts automatically identifies the signal and estimates the intrinsic Doppler shift. In the case of a weak carrier signal, or if severe interference is present (especially at the lower operational band around $2\,$ GHz, the so called S-band), additional visual inspection is required by the user. We distinguish two types of RFI:

1. • Spurious random interferences are present during short duration of a scan. Frequency peaks are randomly encountered throughout the spectra, but only for short periods of time. This type of interference is generally not problematic as the software automatically discards it. Its origin is usually associated with the instrumentation noise and RFI reflections.

2. • Steady frequency tones present constantly for the entire duration of the scan. The origin of these tones is associated with human-made signals, such as near-by tracking stations, communication towers, and radar systems. They could be confused with a spacecraft signal but can be often recognised since their frequency does not shift in time. A similar to an RFI impact on the signal detection can be produced by the phase calibration tones used at VLBI telescopes (see Section 3.1.4).

We have developed a utility that autonomously distinguishes the RFI peaks from any carrier signal transmitted by satellites and spacecraft. The huntSC.py script can also be found in the pysctrack package. The program inspects the entire frequency bandwidth observed for potential RFI peaks. Based on SNR and the Doppler shift, it identifies all spacecraft tones present in the full bandwidth. The spacecraft carrier frequency is a priori known and varies on the order of hundred kHz depending on the spacecraft’s motion. While the search for the tone is not usually needed, it is useful in multiple scenarios: for environments where multiple spacecraft communicate simultaneously with ground-stations (see Figure 4 for example); for the first time observing a particular spacecraft; for spacecraft with weak transmission power, or for locations with severe RFI.

Figure 3. (i) Processing time for the full pipeline of a single 19-min scan. The configuration used here is $8\,$ MHz bandwidth and $3.2\,$ M FFT and $5\,$ s integration time. (ii) Processing time for the full pipeline of six scans processed in parallel using the same observation settings.

3.2. Software performance

SWspec and SCtracker are compact applications intended for deployment on the generic setup of VLBI stations. They have been tested in most common Linux operative systems and macOS. The software takes advantage of the Intel Integrated Performance Primitives (IPP) librariesFootnote f to optimise the code and maximise the performance of computationally intensive tasks. Benchmarks on Linux platforms showed an improvement of $40\,\%$ in CPU performance using the IPP libraries compared to the standard Linux Fastest Fourier Transform in the West package. The software can be used on the same data acquisition system where the data have been recorded, avoiding large data transfers to the processing centre.

The full data processing pipeline for a typical 19-min tracking scan takes around 50–55 min including the observing time. The breakdown of the different tasks is shown in Figure 3 (left panel). In this example, we used SWspec with four-cores and SCtracker in single-core mode. In 2010, the processing of a single 19-min scan took about 2 h of computational time Molera Calves (Reference Molera2012). The latency has been improved by more than $200\,$ % between 2010 and 2020. The right panel shows the total computational time when processing six 19-min scans. Since the data processing is performed in parallel, the largest component to the overall time is now the antenna observing time. Six scans can be computed in parallel within 160 min.

Table 4. SWspec performance benchmarks with different input settings and using parallel vs single processing modes. The table includes integration time (dT), number of FFT’s points (FFT’s) [million points], bandwidth (BW) [MHz], FFT overlap (Ovl), number of cores in parallel ( $\#$ c), FFT calculation time (FFT T), total time (Total T).

SCtracker was developed to run on a single-core processors. The parallelisation is done via a Python script included in pysctrack, which processes multiple scans simultaneously. SWspec, on the other hand, is capable of processing a single scan on a multi-core system; this can be selected from the configuration file. A benchmark of the SWspec performance with different input and multi-core settings is presented in Table 4. I/O operations become the main bottleneck when reading the file from a single hard disk and using more than six cores in parallel. In that case, it is recommended to store the files in a RAID system or a vbs-typeFootnote g file disk system.

VLBI generates large volumes of raw data, with typical data rates of several gigabits per second per station resulting in a total volume of accumulated raw data per experiment of hundreds of terabytes. In order to maximise the sensitivity of VLBI observations, observing systems are widening the frequency bandwidth. Current receivers are capable of collecting radio signals continuously from 2 to $14\,$ GHz. This observed bandwidth is then segmented into frequency bands, or frequency channels, which are then recorded by DAS Molera Calves (Reference Molera2012).

Several DAS types are utilised in modern VLBI facilities. The most widely used are the Flexbuff systems developed at the Metsähovi Radio Observatory and the Mark5/6 (see review in Lindquist & Szomuru Reference Lindquist and Szomuru2014). These units record multiple digitised frequency channels in files ranging in size from tens to hundreds of gigabytes. The data format includes the observing mode, antenna code, and the time of the observations, among other metadata. The spacecraft tracking software is compatible with most of the major versions of the data format developed for VLBI in the recent years, including Mark5A/B, RadioAstron Data Format (RDF) and the latest VLBI Data Interchange Format (VDIF). SWspec and SCtracker can also read non-formatted data without timestamps produced by non-standard radio astronomical instrumentation.

4.1. Data processing pipeline

Radio telescopes equipped with VLBI instrumentation provide for acquisition of wide–band data, making possible to capture multiple spacecraft radio signals simultaneously. The transmission band for deep space navigation at X-band, between $8\,400$ and $8\,450\,$ MHz, fits in one or two standard VLBI frequency channels. In the case of Martian missions, where no less than 10 craft are operating in 2021, all spacecraft carrier signals can be detected and processed simultaneously by a single VLBI radio telescope. Figure 4 shows the data recorded in one experiment in February 2020 in which five craft (one lander and four orbiters of Mars) were transmitting simultaneously to ground tracking stations. We can distinguish NASA Mars Insight ( $8\,404.5\,$ MHz), NASA Mars Odyssey ( $8\,407.2\,$ MHz), ESA ExoMars Trace Gas Orbiter ( $8\,409.7\,$ MHz), ESA Mars Express ( $8\,419.9\,$ MHz), and NASA Mars Reconnaissance Orbiter ( $8\,439.7\,$ MHz). The power received, the frequency, and the stability of the phase are different, and they depend on the characteristics of each spacecraft. The topocentric frequency detections in a session can vary by an order of few mHz to several Hz per second due to the Doppler shift.

Figure 4 illustrates the X-band electromagnetic spectrum captured at a VLBI radio telescope. We can run the full software pipeline to any detected planetary spacecraft to extract the intrinsic parameters of the carrier signals. Multiple stages of the processing pipeline are presented in the four panels of Figure 5. Panel (i) shows the results of SWspec, that is the accumulated 227 integrated spectra ( $5\,$ s) for a scan of 19 min of the Mars Express carrier signal acquired with the 30-m Ceduna radio telescope on $2020.02.23$ . The plot shows the Doppler shift range of $\sim\! 4.3\,$ kHz in 19-min ( $3.7\,$ Hz/s or $134\,$ mm/s) and the peak maximum of the normalised spectral (SNR $\sim40\,$ dB).

Figure 4. Radio spectrum at X-band showing signals from several spacecraft operating on the surface of Mars (Mars Insight) and in various aerocentric orbits (Mars Express, Mars Odyssey, Mars Reconnaissance Orbiter, and ExoMars Trace Gas Orbiter).

Figure 5. (i) Spectral power of Mars Express observed with the 30-m Ceduna radio telescope on $2020.02.23$ . The plot shows the Doppler shift on a 19-min scan using dF of 5 Hz and $5\,$ s integration time. (ii) Spectral power of Mars Insight observed with the 65-m Tianma radio telescope on the same epoch. The SNR of Mars Insight was $20.1\,$ dB and a Doppler shift of $400\,$ Hz in 19-min (iii) MEX carrier tone detected in narrow band ( $20\,$ Hz) after $1140\,$ s of integration with a frequency resolution of $0.4\,$ Hz. The SNR of the signal is $51\,$ dB, and the tracking of the signal is within a precision of $0.8\,$ mHz. (iv) Same carrier tone detected in ultra narrow band ( $1\,$ Hz), now the SNR improves to almost 60 dB.

For comparison, panel (ii) shows the signal of Mars Insight after SWspec for the same scan. Instead of using data from Ceduna, the plot presents data from Tianma (65-m) because the lander signal was weak. NASA’s Mars Insight is a lander deployed on the surface of Mars and equipped with a low-gain antenna. Due to the lower transmission power combined with a season of dust storms on Mars, the detection was 20– $30\,$ dB weaker than other standard missions orbiting Mars at the moment. The Doppler shift was also small, around $400\,$ Hz over 19 min ( $0.35\,$ Hz/s). The average SNR in the scan was only $20.1\,$ dB, while the SNR for Ceduna using 10-s integration time was $11\,$ dB (i.e. barely detectable).

The detection of the initial tone is improved by using all three packages in SDtracker to completely phase-stop the spacecraft carrier in a narrow band. In panel (i), the spectral tone is spread across multiple frequency channels due to the Doppler shift introduced by the motions of the ground-station and spacecraft. In panel (iii), the spacecraft carrier signal is corrected for the Doppler shift with a mHz resolution. The phase correction includes the compensation applied first in SCtracker and then in dPLL.

Finally, the panel (iv) shows the result of applying all the possible corrections with the third set of phase polynomials to the carrier tone in a ultra-narrow band of a single Hz. In this case, we can see that all spectral power is distributed in a single spectral bin, for the duration of the entire scan. Demonstrating that the phase-compensation has worked as expected to stop the signal. The SNR now was almost 60 dB.

After SWspec, the first estimate of the apparent topocentric frequency detections (r0i) of each spacecraft signal are calculated. The Doppler shift depends on the spacecraft tracked, especially on their orbits, and on the ground-station on Earth. Examples of apparent topocentric frequency detections are shown in Figure 6 and 7. In the first plot, a dynamic spectrum of the MEX signal shows the total Doppler shift ( $874\,$ Hz) during a 1-min scan as the spacecraft moved away from the line of sight of the receiver. The Doppler shift rate was $15\,$ Hz/s. The single scan was observed on $2020.06.21$ with Ceduna radio telescope. This dynamic range was generated from the apparent topocentric frequency detections calculated after dPLL (r2i)

In Figure 7, the topocentric frequency detections of multiple telescopes of the Mars Express MultiView phase referencing experiment are shown. The experiment was observed on $2020.06.21$ with 10 VLBI radio telescopes. In this experiment, a fast nodding cycle of 1 min on-spacecraft and 1 min on each of the three reference sources was used. The experiment lasted 3 h with the following antennae participating: Hobart-26, Hobart-12, Katherine, Yarragadee, and Ceduna (Australia), Warkworth (New Zealand), Tianma and Urumqi (China), Hartebeesthoek (South Africa), and Zelenchukskaya (Russia). The topocentric frequency detections are highly correlated, since the main contributor to the Doppler shift is the spacecraft motion. Differences are driven by the different motions of each station which depend on their terrestrial coordinates.

The main data products of the SDtracker also include the spacecraft carrier phase and the residual phase with respect to local reference clock. The variations of the residual phase usually range from $100\,$ mrad to $20\,$ rad Molera Calvés et al. (Reference Molera Calvés2014). Figure 8 presents a comparison of the residual phase retrieved from a spacecraft carrier (upper panel) and from a pcal tone (bottom panel) internally generated for the calibration of the system noise. Both tones were present in the same scan recorded at Hartebeesthoek on $2015.03.30$ . MEX spacecraft was at $2.35\,$ AU with a solar elongation angle of $19.5\,$ degrees. The phase of the spacecraft signal fluctuates due to the propagation media between observer and target. The major contribution is due to the interplanetary plasma scintillation. For this particular experiment, the standard deviation ( $\sigma_{\rm sc}$ ) was $0.521\,$ rad. Observations at high solar elongations, at low horizon, and/or during bad weather conditions, the phase fluctuations can be dominated by tropospheric and ionospheric turbulence. The pcal tone provided an estimate of the system phase noise $\sigma_{\rm sn}$ , which was equal to $0.014\,$ rad.

Breaks in the residual phase can occur due to fast changes in the state vectors of the spacecraft, low SNR because of a low transmission power or a insensitive ground antenna, a data gap caused by a communications failure, and recording problems, — among others. Some of these breaks can be fixed by tuning the time lags, spectral resolution, FFT zero-padding, or overlapped factors in SCtracker and dPLL. Other breaks might require the recorded data to be split into multiple sub-scans and for the data processing to be performed on each of them individually. For example, this is the case for the gravitational redshift experiment carried out with the RadioAstron mission, in which the operations mode between spacecraft and ground station can switch from one-way to two-way modes multiple times in a single scan Litvinov et al. (Reference Litvinov2016); Litvinov et al. (Reference Litvinov2018).

SDtracker can process data obtained in an open-loop regime, which can be implemented as one-, two-, or three-way mode. In one-way mode, the spacecraft transmits the signal to the ground station using an onboard ultra stable oscillator (USO) as the frequency reference. In the two-way mode, the ground station sends a signal to the spacecraft that is locked to a local frequency standard (typically a hydrogen maser) and the spacecraft retransmits the signal to the same ground station. In the three-way mode, the receiving and transmitting ground stations are different. The main difference between the one- and the two-/three-way regimes is that, in the latter, the phase stability is higher in most cases. Due to the better stability of the signal, some radio science experiments are preferably conducted in the two-/three-way mode. However, Pogrebenko et al. (Reference Pogrebenko, Gurvits, Campbell, Avruch, Lebreton, van’t Klooster and Wilson2004) demonstrated successful tracking of the Huygens descent in Titan’s atmosphere supported by the Huygens onboard USO. Another example of a successful use of a one-way mode is the radio occultation experiment Bocanegra-Bahamón et al. (Reference Bocanegra-Bahamón2019). In the latter case, the one-way mode makes possible retrieval of both ingress and egress occultation profiles, as opposed to only ingress in the two-/three-way mode.

4.2. Tracking of deep space missions and Earth satellites

The SDtracker software has been used for different types of spacecraft operating at both interplanetary distances as well as Earth satellites at sufficiently high orbits.

Figure 6. Dynamic spectrum of the MEX Doppler shift observed by the Ceduna radio telescope on $2020.06.21$ , with a topocentric frequency shift rate of $15\,$ Hz/s.

Figure 7. Apparent topocentric frequency detections observed with multiple radio telescopes on a MultiView VLBI phase referencing session. The scans were 1-min long and alternated with reference sources.

4.2.1. Deep space planetary missions

The ESA’s Venus Express (VEX) was the first spacecraft used for extensive tests of the tracking software presented here. VEX was launched in 2005 and operated until December 2014. Most of the communication operations were carried out at X-band ( $8\,419\,$ MHz), including some at S-band ( $2\,296\,$ MHz). We conducted a total of 319 VEX sessions using VLBI radio telescopes.

Another extensively observed spacecraft has been the ESA’s Mars Express (MEX), operational since 2003. MEX sessions with VLBI radio telescopes are still being conducted at the time of this writing in 2021. Furthermore, a number of other spacecraft have been tracked by our group in the framework of various science projects since 2009. In addition to their specific science objectives, these observations enabled us to verify and update the software. In particular, we modified the software to maintain its compatibility with the progressing spacecraft communications and VLBI instrumentation. Table 5 lists all the observations conducted between 2009 and 2021 for which the software described in this paper has been used. The main carrier frequency is an approximation as this can vary by several tens and hundreds of kHz.

Table 5. List of spacecraft observed through $2021.03.31$ used in developing SDtracker at X-band, including the spacecraft name, epochs of the sessions (1), approximate frequency of the main carrier (2), station code (3), mean SNR of the carrier (4), stochastic frequency noise (5), and Doppler shift range over a 19-min scan (6). These values depend on the epoch, ground-station, and settings used in SWspec; hence, we used a standard session.

Cd–Ceduna 30-m (Australia); Ho–Hobart 26-m (Australia); Hh–Hartebeesthoek 26-m (South Africa); Mh–Metsähovi 14-m (Finland); T6–Tianma 65-m (China); Ww–Warkworth 30-m (New Zealand); Wz–Wettzell 25-m (Germany); Yg–Yarragadee 12-m (Australia). CNSA–China National Space Administration, JAXA–Japan Aerospace Exploration Agency, MRO–Mars Reconnaissance Orbiter

^*The spacecraft was probably operating in open-loop mode using the on-board oscillator during those observations.

Figure 8. Phase extraction from the Mars Express signal (top) and telescope phase calibration tone (lower panel). The data were collected on $2015.03.30$ with the Hartebeesthoek 26-m telescope in South Africa for interplanetary scintillation studies.

Venus Express spacecraft was used in three different science cases by Duev et al. (Reference Duev, Molera Calvés, Pogrebenko, Gurvits, Cimó and Bocanegra Bahamon2012), Molera Calvés et al. (Reference Molera Calvés2014) and Bocanegra-Bahamón et al. (Reference Bocanegra-Bahamón2019). VEX had a 24-h near-polar elliptical orbit with a Doppler shift ranging between 400 and $2\,500\,$ Hz per 19-min scans ( $0.35$ to $2.3\,$ Hz/s). The spacecraft was equipped with a high-gain antenna at X-band for the radio communications. The standard SNR received with a $\sim\! 20$ -m dish radio telescope was $10\,000$ ( $40\,$ dB). Duev et al. (Reference Duev, Molera Calvés, Pogrebenko, Gurvits, Cimó and Bocanegra Bahamon2012) demonstrated precise estimates of the orbital state vectors of VEX with EVN antennas. We used analytical measurements of the topocentric frequency detections with SDtracker to verify the model. Data were acquired with multiple $16\,$ MHz frequency channels and scans of 2 min. For the SWspec data processing, we used $3.2$ million FFT points and $5\,$ s integration time. These settings yielded to a frequency accuracy in the order of $50\,$ mHz and an SNR of $32\,$ dB. For the dPLL, we used a 6-degree frequency polynomial fit in the finer Doppler compensation and a 9-degree phase polynomial to extract the residual phase of the carrier line. The final output of the carrier tone was filtered down to $1\,$ Hz.

Molera Calvés et al. (Reference Molera Calvés2014) presented the interplanetary scintillation (IPS) studies campaign based on the investigation of VEX spacecraft carrier phase fluctuations caused by the solar wind. We used consistently a dF of $\,5\,$ Hz and a dT of $\,5\,$ s when running SWspec and SCtracker. We kept the settings constant for better comparison of the results between the stations. For the dPLL, we used 6-order polynomial fits to finer compensate the Doppler and the trend of the phase. The phase scintillation caused by the solar wind was evaluated in a $20\,$ Hz output bandwidth.

Bocanegra-Bahamón et al. (Reference Bocanegra-Bahamón2019) investigated the applicability of this setup for conducting planetary atmospheric studies by means of radio occultation experiments. For the data analysis of the occultation scans, the frequency polynomial coefficients are computed up to a few seconds before the ingress starts and a few seconds after the egress has ended. In this way, the fit corresponds to the Doppler shift due to the relative motion of the transmitter with respect to the receiver, and the remaining frequency residuals will correspond mainly to propagation effects of the signal through the target planet’s atmosphere and Earth’s ionosphere and troposphere. The integration time used in SWspec and SCtracker was $1\,$ s and a 9-order polynomial. The processing with the dPLL was done with an integration time of $0.1\,$ s and an output bandwidth of $50\,$ Hz.

Mars Express spacecraft observations further developed the PRIDE techniques as presented in Duev et al. (Reference Duev2016), Molera Calvés et al. (Reference Molera Calvés2017) and Bocanegra-Bahamón et al. (Reference Bocanegra-Bahamón, Molera Calvés, Gurvits, Duev, Pogrebenko, Cimò, Dirkx and Rosenblatt2018). MEX is in a polar orbit with a period of $7.5$ h, which yielded to the Doppler shift measurements per scan larger that for VEX, as seen in Table 5. The antenna type in MEX the same as in VEX, thus, the SNR received at the stations were similar to previous experiments. Duev et al. (Reference Duev2016) and Bocanegra-Bahamón et al. (Reference Bocanegra-Bahamón, Molera Calvés, Gurvits, Duev, Pogrebenko, Cimò, Dirkx and Rosenblatt2018) estimated the state vectors of MEX during a flyby on Phobos in 2015. We used a combination of long (19-min) and short scans (2-min) alternating with a reference source for OD. To cope with the short scans and the fast dynamics of MEX, we used a dF of $40\,$ Hz and dT of $1\,$ s in SWspec and SCtracker. The degree of the phase polynomial used in dPLL varied depending on the scan’s orbital phase, ranging between 6, 9, and up to 15 close to the actual flyby.

The IPS campaign continued in 2014 using the MEX carrier signal. The work is currently under preparation and covers 7 years of data. Typical Doppler shift ranges from 1 to $20\,$ kHz per 19-min scans (1 to $17\,$ Hz/s). Due to the higher radial acceleration of the observer-spacecraft system, we needed to use lower integration times. Typical dT are 5 or $2\,$ s, while dF is 5 or $10\,$ Hz. To prioritise the consistency on the results between MEX and VEX data, we used the same order of polynomials (6) and output bandwidth (20 Hz).

Finally,the detection of a coronal mass ejection (CME) crossing the line of sight between the MEX spacecraft and Earth was presented in Molera Calvés et al. (Reference Molera Calvés2017). Due to the large frequency fluctuations caused by the CME, we had to use $1\,$ s integration time with a frequency resolution of $5\,$ Hz. The final output of the tone was filtered down to $20\,$ Hz as usual, but we used a $8-$ degree phase polynomial to de-trend the phase in dPLL.

4.3. Single-dish spectroscopy for astronomical applications

Pogrebenko et al. (Reference Pogrebenko2009) conducted spectrometric observations of the Saturnian system with the VLBI radio telescopes of Medicina and Metsähovi at the frequency of the water maser emission near $22\,$ GHz. The data processing was conducted with an earlier version of SWspec. The current version of the program can be also utilised as a single-dish spectrometer directly from any VLBI recorder.

SWspec is used to continue the study of astronomical masers by the University of Tasmania using their network of radio telescopes. Long-term monitoring of the 6.7-GHz class II methanol maser spectral profile is necessary for the early detection of maser flares in high-mass star formation regions. These flare events allow for the detection of rare maser transitions, which can then be used to investigate the physical properties in these star-formation regions Breen et al. (Reference Breen, Sobolev, Kaczmarek, Ellingsen, McCarthy and Voronkov2019); Brogan et al. (Reference Brogan2019). In the example shown in Figure 9, one can see the spectral profile of the methanol maser G $323.740$ as observed with the Ceduna radio telescope on $2020.09.21$ . We used SWspec with dF of $\,1\,$ kHz and dT of $\,5\,$ s. The peak of the maser is at topocentric frequency of $6668.518\,$ MHz. The radial velocity of the maser has been corrected for the local standard of rest.

Figure 9. G $323.740$ methanol maser observed on $2020.09.21$ with the radio telescope at Ceduna (Australia) as part of the University of Tasmania’s effort of long-term study of methanol masers. Spectra generated using SWspec with $1\,$ kHz frequency resolution and $5\,$ s integration time.

4.4. Bistatic radar

VLBI-equipped radio telescopes can also be set up as a bistatic receiver in radar applications. In these operations, a large antenna transmits a tone to the target of interest, while the receiver antenna(s) observes the bounced echo. The technique has been successfully applied to determine the spin states of asteroids Busch et al. (Reference Busch, Kulkarni, Brisken, Ostro, Benner, Giorgini and Nolan2010) and surveying the space debris with VLBI technology Yajima et al. (Reference Yajima, Tsuchikawa, Murakami, Katsumoto and Takano2007); Montebugnoli et al. (Reference Montebugnoli, Pupillo, Salerno, Pluchino and di Martino2010). In 2021, we have demonstrated the capability of the Hobart and Katherine 12-metre antennas for tracking the echo reflected by a large space vehicle. In the experiment on $2021.03.21$ , the three DSN dishes of Tidbinilla (Australia) transmitted three powerful tones to the 2015-56B rocket. These independent tones were detected at both antennas. The experiment is conducted in collaboration with NASA, and it is the preparatory tests for exploring the features of VLBI radio telescopes towards Space Domain Awareness (SDA). This technique can also be applied to the planned bistatic radar observations of the ESA’s JUICE mission to study Jupiter’s atmosphere.

Figure 10 shows the echoes of the signal transmitted by DSN34 as received at Ke and Hb radio telescopes. The Doppler shift in 30 min was $10\,$ Hz, and $160\,$ Hz, respectively. The Doppler shift included the movement of stations and the correction for the transmitted frequency (the so-called ramp table) to receive a static tone at Katherine.

Figure 10. A single tone, transmitted by the DSN34 antenna, was reflected by the $2015-56$ B space rocket, and its echoes were received at Hb (blue) and Ke (purple). This bistatic radar SDA experiment was conducted on $2021.03.21$ .