A Temporal Algorithm for Satellite Subset Selection in Multi-Constellation GNSS

GNSS receivers convert the measured pseudoranges from the visible GNSS satellites into an estimate of the position and clock offset of the receiver. For various reasons receivers might want to process only a subset of the visible satellites; it would be desired, of course, to use the best subset. In general, selecting the best subset is a combinatorics problem; selecting m objects from a choice of n allows for ( n m ) potential subsets. And since the typical performance criterion (e.g. Geometric Dilution of Precision) is nonlinear and non-separable in the satellites’ locations in the sky, finding the best subset is a brute force procedure; hence, a number of authors have described sub-optimal algorithms for choosing satellites. This paper revisits this problem, especially in the context of multiple GNSS constellations. The paper begins with a review of the existing subset selection algorithms. We note that all of these algorithms are what might be called “snapshot” in nature, selecting a subset for a single, fixed skyview of satellites. Through an example with the GPS constellation, we examine the time-sequential, or temporal, characteristics of the best subset selection noting: • That the best subset at a particular point (snapshot) in time is also the best subset for a significant time interval around that point (typically measured in minutes). • That the changes in the best subset over time are primarily, but not always, due to the loss or gain of a satellite crossing the horizon (or, more precisely, the receiver’s mask angle). Based upon these observations this paper develops several time-sequential, or temporal, algorithms that attempt to track the optimum subset of satellites over time at low computational cost. The accuracy and complexity of the algorithms are assessed with GPS constellation data. On a larger scale, these algorithms are then tested on combined GPS, GLONASS, and Galileo constellations with the resulting performance compared to optimal solutions found via exhaustive search. INTRODUCTION GNSS receivers convert the measured pseudoranges from the visible GNSS satellites into an estimate of the position and the clock offset of the receiver. The typical implementation of the solution algorithm is an iterative, linearized least squares method [1]. Assuming that pseudoranges from m satellites are measured, the “direction cosines” matrix is formed. Using an East, North, and Up coordinate frame this matrix is of the form G =  e1 n1 u1 1 e2 n2 u2 1 .. .. em nm um 1  (1) in which (ek, nk, uk) is the unit vector pointing toward the k th satellite from the assumed receiver position. This matrix is then used to form the pseudoinverse to solve the set of overdetermined, linearized pseudorange equations. Since the pseudoranges themselves are noisy, the resulting solution is also noisy. The accuracy of the resulting solution can be described statistically by its error covariance matrix, which is proportional to the inverse of GG. Rather than considering all of the individual elements of this covariance matrix, it is common to reduce it to a scalar parameter. The most used reduction is the GDOP (Geometric Dilution of Precision), the square root of the trace of this matrix GDOP = √ trace { (GTG) −1 } equivalently, proportional to the square root of the sum of the variances of the four estimates (three of position and one of time). For multiple constellations the definition of G can be extended by appending additional columns to account for different clock biases. Specifically, for L unsynchronized constellations, G is of the form G =  e1,1 n1,1 u1,1 1 0 0 . . . 0 .. .. .. e1,m1 n1,m1 u1,m1 1 0 0 . . . 0 e2,1 n2,1 u2,1 0 1 0 . . . 0 e2,2 n2,2 u2,2 0 1 0 . . . 0 .. .. .. eL,mL nL,mL uL,mL 0 0 0 . . . 1  (the first of the two subscripts is the constellation number, 1 to L; the second is the satellite number within the k constellation, 1 to mk) and the GDOP is still as written above, but now includes the variances of L + 3 variables, three for the receiver’s position and those of the L clock biases. In some instances a GNSS receiver cannot process all of the visible satellites. For example, the issue might be: • That the receiver physically cannot track all of the potential signals – this might be a hardware limitation (due to a fixed number of channels) or the desire to minimize power usage [2, 3]. • That the receiver is using corrections from some augmentation system and that the bandwidth of the correction channel is insufficient to provide information for all of the visible satellites (see, for example, [4]). In such a case the question arises: “If only m of n (m < n) visible satellites can be processed, which m should they be?” Consider the question for the GDOP criterion. Since the GDOP is nonlinear and non-separable in the satellites locations in the sky, finding the best subset is a combinatorics problem; selecting m objects from a choice of n allows for ( n m ) possible subsets. If the receiver limits its attention to the 15 or so visible GPS satellites then a brute force comparison of all subsets to find the optimal subset is possible (for n = 15, ( 15 m ) < 6, 500 for all 3 < m < 15, well within modern computational capability). The advent of other GNSS constellations exacerbates this problem. For example, desiring to process 16 of 35 visible satellites (such as frequently occurs with GPS, GLONASS, and Galileo) brute force comparison is no longer viable (e.g. with n = 35, ( 35 16 ) > 4 billion, impossible for real time usage). This question of selecting a subset of the possible satellites is not new in the navigation literature; multiple authors have described sub-optimal algorithms for choosing the satellite subset. (Note, however, that this question is still timely; all three of the papers referenced above, [2–4], are from 2016.) Most of these sub-optimum methods are greedy in nature (described below) and have been developed under the assumption of a single constellation. Further all of these algorithms are what might be called “snapshot” in nature, selecting a subset for a single, fixed skyview of satellites. This paper revisits this problem, especially in the context of multiple GNSS constellations. The next section briefly reviews the existing subset selection algorithms. This review is followed by a motivational example with the GPS constellation, examining the time-sequential, or temporal, characteristics of the best subset. Based upon these observations, this paper develops several time-sequential, or temporal, algorithms that attempt to track the optimum subset of satellites over time at low computational cost. The accuracy and complexity of the algorithms are assessed with GPS constellation data. Finally, these algorithms are tested on combined GPS, GLONASS, and Galileo constellations with the resulting performance compared to several optimal solutions found via full search (brute force) computation. PREVIOUS SUBSET SELECTION ALGORITHMS Multiple authors have presented sub-optimum satellite selection procedures; a number of these employ alternative performance measures beyond GDOP. These include volume of the polytope formed by the satellites [5, 6] and cosines of the angles between pairs of satellites [7–9], as well as combinations of vertical and horizontal protection limits [3]. The sub-optimum algorithms tend to be greedy algorithms, making optimum decisions on a sequence of smaller problems. For example, to generate a subset of size m the authors of [10] suggest starting with the full set of n satellites and iterating the following steps: 1. Assuming that the current subset consists of k satellites, compute k GDOPs, one for each proper subset of k − 1 satellites. 2. Of these k values identify and remove that satellite which results in the smallest increase in GDOP; the result is a subset of size k − 1. 3. Repeat steps 1 and 2 until k = m. This algorithm is greedy in that it makes an optimum choice at each step in the iteration although the final result might not be the global optimum. Specifically, a poor (but still locally optimum) choice at one step might lead to a globally sub-optimum solution for future iterations. It is noted in [10] that the loss to the global optimum for small constellations appears to be small; however, there is no guarantee that this is true for larger numbers of satellites. In a similar way it is possible to grow the subset greedily from the best set of 4 [11]. Another sub-optimum algorithm suggests starting with a subset of size m (and one could discuss how to select this initial subset!) and then iterate in a greedy fashion – growing the subset to m + 1 satellites by adding the most helpful (with respect to GDOP) of the unused satellites and then shrinking back down to m by removing the least helpful one, denoted a “revolving door” method [12] – until a equilibrium is reached. Simulated annealing has also been proposed as a technique to implement a GDOP-based satellite selection algorithm [13]. Finally, an algorithm could try to mimic the minimum GDOP constellation from the known results on lower bounds to GDOP. Specifically, in [14] it is shown that for m satellites the minimum GDOP is attained by a constellation consisting of approximately 30% of the satellites at zenith and the remaining 70% at the horizon. An algorithm could, then, choose high elevation satellites to match the number desired to be at zenith and then attempt to place the remaining satellites near the horizon following “balance” [15,16]. A MOTIVATING EXAMPLE To motivate the time sequential, or temporal, algorithms developed in this paper, consider the case of selecting a subset of satellites from the GPS constellation. The data for the experiment consists of azimuth and elevation angles of the GPS satellites visible at our location in New England for a 24 hour period; the data was decimated to one minute snapshots (a total of 1440 snapshots). Figure 1 shows the n