Orbit Determination Before Detect: Orbital Parameter Matched Filtering for Uncued Detection

The motivation for this paper is to further develop techniques for the surveillance of space with the Murchison Widefield Array (MWA) using passive radar. The paper is particularly concerned with developing techniques for uncued detection over a wide field of regard. The MWA is a low frequency (70 - 300 MHz), wide field-of-view, radio telescope located in Western Australia [1]. The MWA has demonstrated the incoherent detection of the International Space Station (ISS) [2] and other, smaller, objects in orbit [3]. However, for coherent processing, methods compensating for all aspects of motion migration are required in order to detect smaller satellites and space debris [4]. As passive radar systems have no control over the transmitter used for detection, improving processing gain through extended Coherent Processing Intervals (CPIs) is a method used to achieve the required sensitivity [5]. Orbital trajectories are ideal targets for such techniques, as stable and predictable relative motion allows for simpler measurement models. Such techniques have also been used with active radar, for improved sensitivity and processing gain [6] [7].

Consisting of 256 tiles spread across many square kilometres, the MWA’s sparse layout²²2At FM radio frequencies, even the compact configuration of MWA Phase II is sparse [8]. provides high angular resolution. Objects in orbit will therefore transit many beamwidths per second at the point of closest approach. Because of this, high angular resolution (normally a desirable attribute) can result in significant angular migration. Highly eccentric orbits will transit significantly faster. This is particularly challenging for the uncued detection of small objects, where longer integration times are needed to achieve sufficient sensitivity.

Individual radar detections consisting of a single measurement of range, Doppler, azimuth and elevation, only define a broad region of potential orbital parameters [9]. This region may be constrained by incorporating angular rates [10], and even further by including radial acceleration and jerk [7]. Usually, many radar detections are required to perform initial orbit determination. The mapping between radar measurement space and orbital parameters is an ongoing area of research [11].

A standard timeseries matched filter is a function to detect reflected copies of a reference signal $d(t)$ in the surveillance signal $s(t)$, specifically copies delayed by $\tau$ and frequency shifted by $f_{D}$:

This matched filter can be extended to more complicated motions by dechirping (or even applying higher order corrections to) the motion-induced frequency shift. For example, instead of matching to the radial velocity with a Doppler shift of $f_{D}$, higher order motions could be matched with a time varying frequency (that can be represented as a polynomial phase signal) given by $f_{D}+f_{C}t$, where $f_{C}$ is proportional to the radial acceleration. This can be extended to an arbitrary number of parameters at the cost of adding extra dimensions to the matched filter outputs. To account for any range migration, the delay term $\tau$ will also need to be a function of time to match the radial motion.

For a receiver array consisting of $N$ elements, the surveillance signal $s(t)$ can be formed by classical far-field beamforming in a direction of interest such that:

where $s_{n}(t)$ is the received signal at the $n^{\text{th}}$ antenna, $\boldsymbol{u}_{n}$ is the position of the $n^{\text{th}}$ antenna, and $\boldsymbol{k}(\theta,\phi)$ is the signal wavevector for azimuth $\theta$ and elevation $\phi$. Time varying adjustments can be made to every measurement parameter to create a filter, $\chi$, matched to the exact motion of an object with range $\rho(t)$ and slant range-rate $\dot{\rho}(t)$, in time-varying directions given by azimuth $\theta(t)$ and elevation $\phi(t)$:

where the delay to the target is now given by the total path distance scaled by $\frac{1}{c}$, and the Doppler shift is given by $\frac{2\dot{\rho}}{\lambda}$.

The most common elements used to parameterise an orbit are the Keplerian, or classical, orbital elements. These elements directly describe the size, shape, and orientation of an orbital ellipse (with one focus being at the centre of the central body), and the position of an object on this ellipse at some epoch, in the Earth-Centered Inertial (ECI) coordinate frame [12]. The ECI coordinate frame has its origin at the centre of the Earth, but it does not rotate with the Earth. It is also worth noting that a Keplerian orbit can, in fact, be any conic section. However, in this paper, it is assumed that orbits describe Earth-captured closed orbits.

The Keplerian orbital parameters are: the semi-major axis, $a$, and eccentricity, $e$, defining the size and shape of the ellipse; the right-ascension of the ascending node, $\Omega$, and inclination, $i$, which define the orientation of the elliptical plane to the Earth’s equatorial plane; the argument of periapsis, $w$, defining the orientation/rotation of the ellipse in the orbital plane; and finally, the true anomaly, $\nu$, defining the position of the object on the ellipse (refer to Figure 1).

It is also assumed that the only force acting on the object in orbit is due to the gravity of the dominant mass³³3Uniform acceleration does not take into account the ellipsoidal/oblate nature of the Earth or other forces, such as micro-atmospheric drag, solar weather, and gravity due to other celestial bodies. For the short duration of a single CPI, these factors are generally negligible., with the acceleration due to the Earth’s gravity $\boldsymbol{\ddot{r}}$, given by:

where $\mu$ is the standard gravitational parameter for the Earth.

Given the orbital parameters, and the acceleration due to the Earth’s gravity, the Cartesian position $\boldsymbol{r}$, and velocity $\boldsymbol{\dot{r}}$, for an object in Earth orbit is completely deterministic and is given by:

where $\boldsymbol{P}$ and $\boldsymbol{Q}$ represent axes of a coordinate system co-planar with the orbital plane in the Cartesian ECI coordinate frame (given by axes $\boldsymbol{I}$, $\boldsymbol{J}$, and $\boldsymbol{K}$). The third axis, $\boldsymbol{W}$, is perpendicular to the orbital plane [12]. These vectors are described by:

Note that a complicating factor with the ECI reference frame is that a nominally stationary position on the surface of the Earth, such as a fixed radar sensor, will have significant motion.

The methods described above enable coherent processing that matches orbital parameters; however, they are not suitable for searching to perform uncued detections. The parameter space is far too large to be practically searched, and the vast majority of orbits will not correspond to passes within a region of interest above the sensor. Although, as stated earlier in Section II-B, alternatives to the Keplerian parameter set are available. In fact, it is possible to parameterise a Keplerian orbit with the Cartesian position and velocity to constitute the six elements [12]. It is also possible to utilise combinations of both sets of elements in other formulations.

Instead of searching through classical orbital parameters, three parameters can be expressed as a hypothesised ECI position within a search volume of interest. This ensures any hypothesised orbit, determined from these initial parameters, will be within the search volume. Given this potential orbital position, $\boldsymbol{r}$, only three more additional parameters are needed to fully define an elliptical orbit. Although the three elements forming the orbital velocity could be treated as free variables, the majority of possible velocities would not correspond to valid Earth-captured orbits. Instead, given position $\boldsymbol{r}$ and semi-major axis $a$, the magnitude of the velocity of the corresponding orbit is given by the Vis-Viva equation [12]:

Furthermore, given position $\boldsymbol{r}$ and eccentricity $e$, the semi major axis length will itself be constrained between the potential limits of the orbit’s apogee and perigee ranges:

The semi-major axis is also constrained by realistic limits on an orbit’s range, as well as a sensor’s maximum detection range, represented by minimum and maximum allowable periapsides, ${rp}_{min}$ and ${rp}_{max}$:

Another constraint is the constant angular momentum of the orbit, $\boldsymbol{h}$. This vector is perpendicular to the orbital plane, parallel to $\boldsymbol{W}$, with a magnitude depending on the size and shape of the ellipse:

This cross-product may be rewritten to form an expression for the inner product between the position and velocity:

Combined with the magnitude of the velocity, from the Vis-Viva equation (24), as well as the magnitude of the constant angular momentum (27), an expression for this inner product can be formed which depends solely on the position $\boldsymbol{r}$ and the size and shape of the orbital ellipse:

Additionally, the specific relative angular momentum vector, $\boldsymbol{h}$, is perpendicular to both the orbital position $\boldsymbol{r}$ and orbital velocity $\boldsymbol{\dot{r}}$. This leads to the expressions $\boldsymbol{r}\cdot\boldsymbol{h}=0$ and $\boldsymbol{\dot{r}}\cdot\boldsymbol{h}=0$, which result in another constraint on the velocity, dependant on the right ascension of the ascending node, $\Omega$:

These expressions lead to a simple geometric solution for determining orbits when $\boldsymbol{r}$ (and other parameters) are known, and $\boldsymbol{\dot{r}}$ is unknown. For determining $\boldsymbol{\dot{r}}$, (24) defines a sphere of radius $\sqrt{\mu(\frac{2}{\lvert\boldsymbol{r}\rvert}-\frac{1}{a})}$, representing valid orbits in the velocity vector’s element space. Additionally, (29) defines two parallel planes of valid orbits, which intersect with (24) to define two circles. Finally, intersecting these two circles with the plane defined by the position and the right ascension of the ascending node, $\Omega$, (30) will result in a maximum of four intersection points, that is, four velocities, each corresponding to a valid orbit. An example diagram is shown in Figure 3. Although this means that a choice of six orbital parameters will result in up to four potential orbital matched filters, this approach will be far more efficient than methods outlined earlier in this section, as the orbit will be within the search volume, and each parameter choice restricts the range of subsequent parameters.

Therefore, given an orbital position, $\boldsymbol{r}$, a choice of eccentricity, $e$, semi-major axis, $a$, and right ascension of the ascending node, $\Omega$, four potential orbital velocities, $\boldsymbol{\dot{r}}$, are calculated, which leads to an expression for the complete matched filter:

The proposed method tests for only realistic orbits in a given search region. Also, given a set of orbit parameters, this matched filter should maximise a radar’s sensitivity to that orbit. Additionally, a detection in this matched filter corresponds to a detection in the orbital element space, providing initial orbit determination from a single detection.

This style of trajectory-match approach, has several advantages beyond just maximising sensitivity to motion models. Coupling measurement parameters together through a trajectory model can improve achievable resolution compared with using separate independent measurement parameters. As an example, a radar’s range resolution is determined solely by the signal bandwidth, but its Doppler and Doppler-rate resolution improve with the CPI length.

Through coupling the measurement parameters with the trajectory model, as a radar can resolve finer Doppler and Doppler rate measurements it can essentially resolve finer trajectory states. This can potentially improve target localisation as increasingly accurate state measurements could localise a target within a single range bin.

The flexibility of the geometric formulation in Section III-A allows radar parameters to be used alongside, and in place of, other orbital parameters to constrain the search space. A Doppler shift $f_{D}$ will define another plane in $\dot{\boldsymbol{r}}$ space, given by:

Equation (32) can be used to search for a particular Doppler shift instead of one of the orbital parameters. This is useful because it allows a blind search to constrain the search-space solely for objects in orbit at their point of closest approach to the sensor. As an object is passing overhead, its point of closest approach will correspond exactly with it being at zero Doppler, which is when it is most detectable⁴⁴4This may not necessarily hold in all instances, depending on particular beampattern and radar cross section factors.. If a radar is unable to detect an object at its point of closest approach, at its minimum range, there is little value trying to detect it as it moves further away, towards the horizon.

Another benefit to applying this constraint is that, as Doppler is proportional to the range-rate, this constraint will also restrict the orbit search-space to a point of minimal (or zero) range migration, which greatly simplifies matched-processing⁵⁵5Depending on the CPI length, it may be possible to make ${\rho}(t)\approx{\rho}$..

The vast majority of the objects in an Earth-captured orbit are in a circular, or near-circular, orbit. Searching solely for objects in a circular orbit greatly decreases the potential orbital search space. A circular orbit means the eccentricity of the orbital ellipse is zero, $e=0$, and so (25) becomes $a=\lvert\boldsymbol{r}\rvert$. In a circular orbit, the position and velocity vectors will always be perpendicular, so (29) simplifies to $\boldsymbol{r}\cdot\boldsymbol{\dot{r}}=0$, a single plane instead of two parallel planes. The result is that a three-parameter search, within a region of interest, provides sufficient information to match the closest approach of objects in a circular orbit. For a given position in a search-region, there will be at most two possible orbits to match against (determined from the intersection of (24), (29), and (32)). This type of search approach, attempting uncued detection of the most common types of orbit when they are most detectable, is a far more realisable and practical approach than a completely unbounded search through measurement parameters. Additionally, for an eccentric orbit, the orbital velocity and position are perpendicular at perigee [12]. For typical radar detection ranges, an object in a highly eccentric orbit is likely to be within a radar’s field of regard solely at, or near, perigee. Because of this, the same simplification of $\boldsymbol{r}\cdot\boldsymbol{\dot{r}}=0$ could be used to reduce the number of potential orbits.

Coupling together measurement parameters is not necessarily new; however, incorporating such techniques into the detection stage offers some significant advantages. By coupling together the measurement parameters using these ODBD methods, it is possible to apply this matched filtering to single beam radar systems. This could be a post-beamformed surveillance signal from an array or even a classic narrowbeam tracking radar. Because the trajectory model determines all measurement parameters, a particular polynomial phase signal which results in a detection is coupled to a particular location and orbit. This is shown in (31). The beamforming parameters do not determine the location; rather the (hypothesised) location determines the beamforming parameters. Removing the array processing, as in (33), does not remove the ability to localise a target using the algorithm.

In the case of a narrow beam radar, the pointing of the beam will be incorporated into the algorithm by determining the search region that is used. Because it handles sensor motion, this type of processing would be ideal for a satellite-based sensor, with the sensor location term $\boldsymbol{q}(t)$ (or its instantaneous components $\boldsymbol{q},\boldsymbol{\dot{q},\boldsymbol{\ddot{q}}}$, etc.) themselves determined by a known orbit rather than the motion of the Earth.