A method for estimating positioning performance of a mega-satellite network

Abstract

This invention discloses a method for estimating the positioning performance of a mega-satellite network. Under TOA, TDOA, and FOA measurements, this invention uses the azimuth and elevation angles of the satellites relative to the receiver to represent items in the Fisher information matrix. After substituting relevant satellite constellation parameters, the average CRLB closed-form solution of the considered constellation is obtained by weighting the joint probability density of elevation and azimuth angles. This average CRLB can theoretically reflect the average positioning performance of a constellation. This invention is of great significance for exploring the positioning potential of mega-satellite constellations and guiding constellation construction.

Description

Technical Field

[0001] This invention relates to the field of satellite positioning technology, and more specifically to a method for estimating the positioning performance of a giant satellite network. Background Technology

[0002] Satellite positioning has been widely used for a long time. With the construction of numerous giant low-Earth orbit satellite constellations, the number of satellites in orbit has grown exponentially. Evaluating the positioning performance of a large number of satellites is of great significance for tapping the positioning potential of giant satellite constellations and guiding constellation construction.

[0003] The Cramer-Rao lower bound (CRLB) is the theoretical lower bound of the error for an unbiased estimator in parameter estimation problems. In positioning problems, the CRLB coefficients are calculated based on the relative geometry of the receiver and satellites. Multiplying these CRLB coefficients by the measurement noise variance effectively characterizes the theoretical lower bound of the positioning error (CRLB). Existing research often involves providing the satellite distribution at a specific time, substituting the known satellite position information into the Fisher Information Matrix (FIM) J, finding the trace of the inverse of the FIM, and then multiplying the result by the measurement error variance to obtain the lower bound of the positioning error variance (CRLB) for that specific time. However, this method can only calculate the lower bound of the positioning error at a single time, characterizing the accuracy of the algorithm at that moment. Since the geometry of the satellite network changes over time, the above method cannot characterize the average positioning performance of a satellite constellation.

[0004] There is currently limited research on methods for estimating the positioning performance of satellite constellations, and no reports have been found both domestically and internationally. Summary of the Invention

[0005] In view of this, the present invention provides a method for estimating the positioning performance of a mega-satellite network, which can calculate the average positioning CRLB of a mega-satellite constellation based on receiver position and constellation parameters under the conditions of time of arrival (TOA), time difference of arrival (TDOA), and frequency of arrival (FOA), thereby characterizing the average positioning performance of the satellite constellation.

[0006] The method for estimating the positioning performance of a giant satellite network according to the present invention includes:

[0007] Step 1: Obtain the time of arrival (TOA), time difference of arrival (TDOA), or frequency of arrival (FOA) of the satellite signal to the receiver, and calculate its error variance. and in, and The unit is m 2 , Convert the units to (m / s) 2 ;

[0008] Step 2, calculate the average positioning CRLB coefficient of the satellite constellation, specifically including:

[0009] S21, based on TOA, TDOA, or FOA, transform the variables of each item in the Fisher information matrix J into the azimuth angle θ and elevation angle of the satellite relative to the receiver.

[0010] S22, based on the latitude and longitude information of the receiver and the orbital parameters of the satellite constellation, calculate the angular division corresponding to the time interval for the receiver to acquire observations, and then calculate the azimuth angle θ and elevation angle of the satellite constellation with respect to the receiver's position. joint probability density function And for the joint probability density function Perform normalization processing;

[0011] S23, perform a process based on each item of the FIM. Weighted average; This is the normalized joint probability density function;

[0012] S24, based on S23 Substituting the weighted FIM terms into the CRLB coefficient calculation formula yields the CRLB coefficient c. CRLB ;

[0013] Step 3: Combining the observation time interval Δt, simulate the average number of visible satellites N of the satellite constellation; calculate the average positioning CRLB of the satellite constellation:

[0014]

[0015] Step 4: Based on the satellite constellation average positioning CRLB obtained in Step 3, complete the positioning performance estimation of the mega-satellite network.

[0016] Preferably, in step 1, the measurement error variance is calculated according to formula (1). Convert to (m / s) 2 as a unit

[0017]

[0018] Where c is the speed of light; f c This refers to the carrier frequency of the satellite downlink signal.

[0019] Preferably, in step S21,

[0020] For the TOA measurement method, the Fisher information matrix J contains... and They are respectively:

[0021]

[0022] Among them, P i =[P i,x ,P i,y ,P i,z [p] represents the three-dimensional position of satellite i, where p = [p] x ,p y ,p z ] represents the three-dimensional position of the receiver, and || represents the L2 norm.

[0023] For the TDOA measurement method, the Fisher information matrix J contains... and They are respectively:

[0024]

[0025] Where P i-1 =[P i-1,x ,P i-1,y ,P i-1,z [] indicates the three-dimensional position of satellite i-1;

[0026] For the FOA measurement method, the Fisher information matrix J contains... and They are respectively:

[0027]

[0028] Among them, V i =[V i,x V i,y V i,z ] and v = [v x ,v y ,v z [ ] represents the three-dimensional velocities of satellite i and receiver i, respectively; and θ v,i These represent the elevation and azimuth angles of satellite i relative to the receiver velocity, respectively; a i This represents the ratio of the magnitude of satellite i's velocity relative to the receiver to the magnitude of its distance.

[0029] Preferably, in step S22, the joint probability density function is... Each element in Normalize using the following method:

[0030]

[0031] The subscript “′” indicates the normalized probability density distribution.

[0032] Preferably, in step S23,

[0033] For TOA or TDOA measurement methods, based on The weighted FIM values ​​are as follows:

[0034]

[0035] Where θ i =Δθ,2Δθ,...,2π / Δθ

[0036] For FOA measurement methods, based on The weighted FIM values ​​are as follows:

[0037]

[0038] Beneficial effects:

[0039] This invention enables the use of azimuth and elevation angles of satellites relative to the receiver to represent items in the Fisher information matrix under TOA, TDOA, and FOA measurements. After substituting relevant satellite constellation parameters, the average CRLB closed-form solution for the considered constellation is obtained by weighting the joint probability density of elevation and azimuth angles. This average CRLB theoretically reflects the average positioning performance of a constellation. This invention is of great significance for exploring the positioning potential of mega-satellite constellations and guiding constellation construction. Attached Figure Description

[0040] Figure 1 This is a schematic diagram of the positioning performance estimation method for giant satellite networks according to the present invention. Detailed Implementation

[0041] The present invention will now be described in detail with reference to the accompanying drawings and embodiments.

[0042] This invention provides a method for estimating the positioning performance of a mega-satellite network. Considering commonly used receiver parameters such as Time of Arrival (TOA), Time Difference of Arrival (TDOA), and Frequency of Arrival (FOA), which are related to the positions of the satellite and receiver, these parameters can be converted into the azimuth and elevation angles of the satellite relative to the receiver. Then, using the variance of these measurement errors, the parameters in the Functional Information Model (FIM) are revised, the CRLB coefficient is calculated, and the average positioning performance of the satellite constellation is evaluated. The flowchart of this invention is shown below. Figure 1 As shown, the specific steps include the following:

[0043] Step 1: Obtain the time of arrival (TOA), time difference of arrival (TDOA), or frequency of arrival (FOA) of the satellite signal to the receiver, and calculate its error variance σ. 2 And the error variance σ 2 The units used should be standardized. Specifically, the unit of measurement error variance for TOA and TDOA should be standardized to m. 2The unit for the measurement error variance of FOA is uniformly set to (m / s). 2 .

[0044] Specifically,

[0045] If TOA measurement is selected, the variance of the pseudorange measurement error in meters is calculated. but

[0046] If TDOA measurement is selected, calculate the variance of the pseudorange measurement error in meters. Because the TDOA differential process doubles the measurement error variance, we get... but

[0047] If FOA measurement is selected, the variance of the Doppler measurement error is calculated in Hertz (Hz). Then, according to formula (1), the measurement error variance is converted into meters per second (m / s).

[0048]

[0049] Where c is the speed of light, taken as 3*10 8 m / s; f c This refers to the carrier frequency of the satellite downlink signal.

[0050] but

[0051] Step 2: Calculate the average positioning CRLB coefficient of the satellite constellation to be evaluated.

[0052] Each item in the Fisher information matrix can be represented by azimuth and elevation angles. This invention leverages this characteristic of the Fisher information matrix, obtaining the average positioning CRLB coefficient by averaging these two variables. First, TOA, TDOA, and FOA measurements are converted into azimuth and elevation angles of the satellite relative to the receiver. Then, the joint probability density function of the satellite constellation with respect to azimuth and elevation angles is calculated. Next, the joint probability density function is weighted for each item in the Fisher information matrix with azimuth and elevation angles as variables. Finally, the CRLB coefficient is calculated based on the weighted Fisher information matrix. This CRLB coefficient obtained using joint probability density weighting is the average CRLB, which theoretically reflects the average positioning performance of a constellation.

[0053] Assume the receiver's observations of the i-th satellite are d. i For the Fisher information matrix J, its structure is as follows:

[0054]

[0055] With j xx and j yz For example, let's introduce its calculation method:

[0056]

[0057] As long as θ and express and The CRLB coefficients can be obtained by finding the trace of the inverse of J.

[0058] For TOA and TDOA measurements, the variables can be converted into the azimuth angle θ and elevation angle of the satellite relative to the receiver. Then, the sampling angle interval Δθ is calculated based on the actual sampling time interval Δt, and the satellite constellation is calculated with respect to... The joint probability density function of θ with intervals of Δθ Finally, it will be with Each item in the Fisher information matrix with θ as variables is used The CRLB coefficient can be calculated by weighting (this coefficient is dimensionless);

[0059] For FOA measurements, since frequency measurements depend not only on the positions of the satellite and receiver but also on their relative velocities, the values ​​at θ and... Based on this, the azimuth angle θ of the relative velocity between the satellite and the receiver is introduced. v and elevation angle Assuming a uniform distribution, we introduce the ratio 'a' between the satellite's relative velocity and the receiver's distance, and use the edge distribution of the satellite's elevation angle relative to the receiver. Calculate the mean of a, E[a]. Then, using E[a], we can... Each item in the Fisher information matrix with θ as variables is used The CRLB coefficient can be calculated by weighting (this coefficient has dimensions due to the presence of 'a', and is s). 2 The unit in FOA in step (1) is (m / s). 2 The product of the measurement error variances yields a result in units of m. 2 CRLB)

[0060] Specifically, it includes the following sub-steps:

[0061] Step 201: Convert the variables measured by TOA, TDOA, and FOA into the azimuth angle θ and elevation angle of the satellite relative to the receiver.

[0062] In the TOA measurement method, the observed quantity is the pseudorange: d i =||P i -p||, where P i =[P i,x ,Pi,y ,P i,z ] represents the three-dimensional position of satellite i, p = [p x ,p y ,p z ] represents the receiver's three-dimensional position (the variable to be determined), and || represents the L2 norm. The J matrix in the TOA method... and They are respectively:

[0063]

[0064] For the TDOA measurement method, the observed quantity is a pseudo-distance difference component: d i =||P i -p||-||P i-1 -p||, where P i-1 =[P i-1,x ,P i-1,y ,P i-1,z [] represents the three-dimensional position of satellite i-1, assuming satellite positions are independent. The J matrix in the TDOA method... and They are respectively:

[0065]

[0066] Since the satellite position is independent, substituting equation (4) into equation (2) to calculate the CRLB coefficient, the CRLB coefficient of TDOA is half of the CRLB coefficient of TOA. Therefore, we no longer need to calculate each term in J for TDOA separately, but only need to halve the coefficient of TOA when calculating the CRLB coefficient at the end.

[0067] For the FOA measurement method, the observations are Doppler observations in meters per second: Where V i =[V i,x V i,y V i,z ] and v = [v x ,v y ,v z [ ] represents the three-dimensional velocities of satellite i and receiver i, respectively. The J matrix in the FOA method... and They are respectively:

[0068]

[0069] FOA measurement in and θ i In addition, it also introduced θ v,i and a i ,in and θv,i Let these represent the elevation angle and azimuth angle of satellite i relative to the receiver velocity, respectively.

[0070]

[0071] This represents the ratio of the magnitude of satellite i's velocity relative to the receiver to the magnitude of its range. In the subsequent derivation, we assume... and θ v,i It follows a uniform distribution and requires a i Find the mean E[a].

[0072] Step 202, according to Kepler's third law, calculate the angular velocity ω of the satellite orbiting the Earth:

[0073]

[0074] Among them, R e Let G be the equivalent radius of Earth, taken as 6371 km; let G be the gravitational constant, taken as 6.67 × 10⁻⁶. -11 m 3 / (kg·s 2 M is the mass of the Earth, taken as 5.965 × 10⁻⁶. 24 kg;

[0075] Step 203: Calculate based on the satellite orbital altitude h, ω, and the time interval Δt between the receiver acquiring the observations. of The angular division Δφ with θ:

[0076]

[0077] Step 204, based on the receiver's longitude λ0 and latitude The orbital inclination i and orbital altitude h of the satellite constellation are used to calculate the azimuth angle of the satellite constellation with respect to the receiver position at intervals of Δφ. The joint probability density function of the elevation angle θ The calculation of the joint probability density can be found in J. Khalife and ZZMKassas, "Performance-Driven Design of Carrier Phase Differential Navigation Frameworks With Megaconstellation LEO Satellites," in IEEE Transactions on Aerospace and Electronic Systems, vol. 59, no. 3, pp. 2947-2966, June 2023, doi:10.1109 / TAES.2023.3234521.

[0078] Next, regarding Normalize: For Each element in Normalize using the following method:

[0079]

[0080] The subscript “′” indicates the normalized probability density distribution.

[0081] Step 205, perform a process based on each item of the FIM. Weighted average;

[0082] If TOA or TDOA measurement is selected, proceed to step 20501; if FOA measurement is selected, proceed to step 20502.

[0083] Step 20501, calculate each term of the FIM under the TOA or TDOA method:

[0084]

[0085] Where θ i =Δθ,2Δθ,...,2π / Δθ

[0086] Step 20502: Since the numerator and denominator of 'a' are independent, we can calculate the mean of both the numerator and denominator, and then substitute them back into the formula.

[0087]

[0088] To obtain the value of E[a], where E||V i -v|| and E||P i -p|| represents the mean values ​​of the relative velocity norm and the distance norm, respectively:

[0089] Step 2050201, for the molecule of E[a], i.e. the mean of the relative velocity norm, E||V i -v||Estimation:

[0090] Since the receiver's velocity is much smaller than the satellite's velocity, the receiver can be approximated as stationary. The relative velocity is determined solely by the satellite's velocity, and its norm is a constant value that can be obtained using Kepler's third law.

[0091]

[0092] Step 2050202, the denominator of E[a] is the mean of the distance E||P i -p||Estimation:

[0093] First, based on the satellite's elevation angle relative to the receiver and satellite elevation angle relative to Earth The conversion method calculates each The corresponding k = 1, 2, ..., π / 2Δθ k = 1, 2, ..., π / 2Δθ, the conversion formula is as follows:

[0094]

[0095] Then, the joint probability distribution Integrating the azimuth angle yields the edge distribution with respect to the elevation angle. for Each item in The marginal distribution is calculated using the following method:

[0096]

[0097] Finally, equation (8) is used to estimate the mean of the denominator of a, i.e., the distance ||Pp||:

[0098]

[0099] Step 2050203: Substitute the results of steps 2050201 and 2050202 into formula (5) to obtain the mean of a.

[0100] Step 2050204, substitute equation (5) and E[a] into equation (2), assuming and θ v,i Following a uniform distribution, each term of the FIM under the FOA method is calculated using the following formula:

[0101]

[0102] Step 206: Substitute the terms in the obtained FIM into the CRLB coefficient calculation formula to obtain the CRLB coefficient c. CRLB :

[0103]

[0104] It is worth noting that, as described in step 201, if the TDOA method is selected, the obtained CRLB coefficient c CRLB It needs to be halved.

[0105] Step 3: Combine the observation time interval Δt to simulate and obtain the average number of visible satellites N of the satellite constellation.

[0106] Step 4: Calculate the satellite constellation average positioning CRLB based on the results of the first three steps:

[0107]

[0108] The obtained average positioning CRLB is the average positioning accuracy calculated using this constellation, which can be used to evaluate the average positioning accuracy of the constellation.

[0109] In summary, the above are merely preferred embodiments of the present invention and are not intended to limit the scope of protection of the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A method for estimating the positioning performance of a giant satellite network, characterized in that, include: Step 1: Obtain the time of arrival (TOA), time difference of arrival (TDOA), or frequency of arrival (FOA) of the satellite signal to the receiver, and calculate its error variance. and in, and The unit is m 2 , Convert the units to (m / s) 2 ; Step 2, calculate the mean positioning Cramer-Rao lower bound (CRLB) coefficient of the satellite constellation, specifically including: S21, based on TOA, TDOA, or FOA, transform the variables of each item in the Fisher information matrix J into the azimuth angle θ and elevation angle of the satellite relative to the receiver. S22, based on the latitude and longitude information of the receiver and the orbital parameters of the satellite constellation, calculate the angular division corresponding to the time interval for the receiver to acquire observations, and then calculate the azimuth angle θ and elevation angle of the satellite constellation with respect to the receiver's position. joint probability density function And for the joint probability density function Perform normalization processing; S23, perform a process based on each item of the Fisher information matrix J. Weighted average; This is the normalized joint probability density function; S24, based on S23 Substituting each term of the weighted Fisher information matrix J into the CRLB coefficient calculation formula yields the CRLB coefficient c. CRLB ; Step 3: Combining the observation time interval Δt, simulate the average number of visible satellites N of the satellite constellation; calculate the average positioning CRLB of the satellite constellation: Where, σ 2 This represents the error variance; if TOA measurement is selected, then... If TDOA measurement is selected, then If FOA measurement is selected, the measurement error variance will be... Convert to (m / s) 2 as a unit Step 4: Based on the satellite constellation average positioning CRLB obtained in Step 3, complete the positioning performance estimation of the mega-satellite network.

2. The method as described in claim 1, characterized in that, In step 1, the measurement error variance is calculated according to formula (1). Convert to (m / s) 2 as a unit Where c is the speed of light; f c This refers to the carrier frequency of the satellite downlink signal.

3. The method as described in claim 1, characterized in that, In step S21 For the TOA measurement method, the Fisher information matrix J contains... and They are respectively: Where, d i For the receiver's observations of the i-th satellite, P i =[P i,x ,P i,y ,P i,z [p] represents the three-dimensional position of satellite i, where p = [p] x ,p y ,p z ] represents the three-dimensional position of the receiver, and |||| represents the calculation of the L2 norm; For the TDOA measurement method, the Fisher information matrix J contains... and They are respectively: Where P i-1 =[P i-1,x ,P i-1,y ,P i-1,z [] indicates the three-dimensional position of satellite i-1; For the FOA measurement method, the Fisher information matrix J contains... and They are respectively: Among them, V i =[V i,x V i,y V i,z ] and v = [v x ,v y ,v z [ ] represents the three-dimensional velocities of satellite i and receiver i, respectively; and θ v,i These represent the elevation and azimuth angles of satellite i relative to the receiver velocity, respectively; a i This represents the ratio of the magnitude of satellite i's velocity relative to the receiver to the magnitude of its distance.

4. The method as described in claim 1 or 3, characterized in that, In step S22, the joint probability density function Each element in Normalize using the following method: The subscript "′" indicates the normalized probability density distribution.

5. The method as described in claim 3, characterized in that, In step S23 For TOA or TDOA measurement methods, based on The weighted Fisher information matrix J consists of the following terms: among them i =Δθ,2Δθ,...,2π / Δθ, For FOA measurement methods, based on The weighted Fisher information matrix J consists of the following terms: Where E represents the mean.

Images