An ionospheric monitor

Abstract

The application belongs to the technical field of ionosphere monitoring, and relates to an ionosphere monitor arranged in a GNSS double-frequency receiver of an ionosphere monitoring station, comprising: a GNSS high-precision OEM board, which obtains GNSS double-frequency pseudo-range observations, carrier phase observations, receiver position information, satellite elevation angles, azimuth angles and satellite ephemeris data according to the length of a time window set, and respectively transmits the data to a TEC inversion calculation module and a storage module; the TEC inversion calculation module is realized based on an embedded processor, uses a sliding time window to pre-process the received data, establishes an ionosphere vertical TEC spherical harmonic function model, adopts Kalman filtering iteration to solve, obtains nonlinear model coefficients and hardware deviations of the receiver and the satellite, realizes quasi-real-time monitoring of the ionosphere vertical TEC value in the space above the monitoring station, and transmits the TEC value to the storage module; and the storage module stores GNSS measurement data and the ionosphere vertical TEC value.

Description

Technical Field

[0001] This invention belongs to the field of ionospheric monitoring technology, and particularly relates to an ionospheric monitoring instrument. Background Technology

[0002] The Earth's upper atmosphere absorbs solar extreme ultraviolet and X-rays, resulting in partial ionization and forming the ionosphere at an altitude of approximately 60–1000 kilometers above the Earth's surface. While the electron density of the ionosphere is typically three orders of magnitude lower than that of the neutral atmosphere, it is still significant enough to affect the propagation of radio waves used in systems such as radio astronomy, radar, satellite navigation and communication, and shortwave communication. Therefore, the ionosphere is an important component of the space (electromagnetic) environment and a crucial element for space weather forecasting.

[0003] The ionosphere exhibits regional characteristics, diurnal variations, seasonal variations, annual variations, and a fluctuating structure, primarily influenced by factors such as the sun, the geomagnetic field, and the neutral atmosphere. Consequently, complex ionospheric disturbances exist, affecting the performance of ground-to-air radio systems whose signals cross or are reflected from the ionosphere. Therefore, conducting ionospheric monitoring is crucial for implementing space weather early warning and forecasting.

[0004] For Global Navigation Satellite System (GNSS) positioning, the ionosphere is the largest source of error. By utilizing the influence of the ionosphere on GNSS measurements and retrieving the total electron content (TEC), an ionospheric correction model can be constructed and applied to GNSS to improve positioning accuracy. On the other hand, conducting all-weather, all-day wide-area ionospheric TEC monitoring using GNSS receiver observations has become an important means of studying the ionosphere and conducting ionospheric disturbance early warning and forecasting.

[0005] The total electron content (TEC) of the ionosphere is usually measured using GNSS dual-frequency signals, which can be used to carry out large-scale ionospheric observation and research around the clock and in all weather conditions. However, at present, GNSS dual-frequency signal measurement data is generally processed by computer on a daily basis, which cannot achieve real-time monitoring of the ionosphere.

[0006] For ionospheric TEC monitoring, two common methods are network monitoring and single-site monitoring. In areas where network deployment is difficult, data transmission is often limited. Even with single-site monitoring, current technologies typically employ post-processing, calculating GNSS satellite hardware deviations using the previous day's measurement data. This requires a large amount of data and lacks real-time performance. In scenarios with limited network bandwidth, such as maritime stations, transmitting data from dozens of satellites in real time results in an enormous data volume. If the TEC data over a station at a given moment is considered as one set, the corresponding navigation ephemeris, pseudorange, and carrier phase data would be two hundred times greater.

[0007] Furthermore, existing GNSS ionospheric monitors only provide slant STEC parameters along the line-of-sight paths between the receiver and each visible satellite. Slant STEC originates from three sources: 1. ionospheric delay, 2. satellite hardware delay, and 3. receiver hardware delay. On the other hand, due to the different elevation angles of different satellites, the STEC along the line-of-sight paths from different satellites to the receiver also varies. Therefore, STEC is a composite measurement and does not represent the actual amount of ionospheric STEC along the line-of-sight path between the receiver and the satellites. Truly meaningful GNSS ionospheric monitoring requires the monitor to provide the vertical TEC above the station, denoted as VTEC.

[0008] For single-station monitoring, there is an urgent need to develop a near real-time GNSS ionospheric TEC monitor. This monitor would process GNSS single-station data within the embedded processor chip of the GNSS receiver, solving the ionospheric monitoring problem in areas where networking and data transmission are difficult, and also addressing the limitations of equipment development environment and computing resources. It should be able to achieve VTEC monitoring with relatively small data volumes and maintain the required monitoring accuracy. Summary of the Invention

[0009] The purpose of this invention is to overcome the shortcomings of the prior art and to propose an ionospheric monitoring instrument.

[0010] To achieve the above objectives, the present invention provides an ionospheric monitor, deployed in a GNSS dual-frequency receiver at an ionospheric monitoring station, the monitor comprising:

[0011] The GNSS high-precision OEM board is used to receive GNSS dual-frequency pseudorange observations, carrier phase observations, receiver position information, satellite elevation angle, azimuth angle and satellite ephemeris data according to the time window set, and transmit them to the TEC inversion calculation module and storage module respectively.

[0012] The TEC inversion calculation module, implemented based on an embedded processor, is used to preprocess the received data using a sliding time window, establish a spherical harmonic function model of the ionospheric vertical TEC, and solve iteratively using Kalman filtering to obtain the nonlinear model coefficients and the hardware deviation between the receiver and the satellite. This enables near real-time monitoring of the ionospheric vertical TEC value above the monitoring station and transmits the TEC value to the storage module.

[0013] The storage module is used to store GNSS measurement data received by the GNSS high-precision OEM board, and also to store the ionospheric vertical TEC value output by the TEC inversion calculation module.

[0014] Preferably, the monitoring device further includes:

[0015] A GNSS antenna is used to receive GNSS satellite signals and transmit the received signals to a high-precision GNSS OEM board; and

[0016] The wireless transmission module is used to periodically send data files from the storage module to the server.

[0017] Preferably, the processing of the GNSS high-precision OEM board includes:

[0018] According to the time window set, receive and decode ranging information from multiple GNSS satellites using dual-frequency pseudorange and carrier phase measurements.

[0019] Receive elevation and azimuth information for each satellite;

[0020] The ranging information for each satellite, including its number, observation time, elevation angle, azimuth angle, GNSS dual-frequency pseudorange, and carrier phase measurement, is cached separately.

[0021] Preferably, the TEC inversion calculation module includes a preprocessing unit, and the processing includes:

[0022] Step S1) Obtain the oblique TEC value of satellite j based on dual-frequency pseudorange observations within the time window set by the cached data. and the oblique TEC value of satellite j based on dual-frequency carrier phase observations

[0023]

[0024] Where f1 and f2 are two frequencies of the GNSS satellite signal, These are the code pseudorange observations of satellite j at two frequencies f1 and f2, respectively. λ1 and λ2 are the carrier phase observations of the j-th satellite at two frequencies f1 and f2, respectively, where λ1 and λ2 are the wavelengths corresponding to the two frequencies f1 and f2, respectively.

[0025] Step S2) Calculate the value of satellite j according to the following formula. and Baseline values ​​between

[0026]

[0027] Where i represents the measurement time of the i-th data in the data received within the time window set duration, and Num represents the total number of times corresponding to the data within the time window set duration; This represents the elevation angle of satellite j relative to the receiver at the i-th data measurement time.

[0028] Step S3) Calculate the oblique TEC value (TEC) of the j-th satellite after carrier phase smoothing pseudorange according to the following formula. s,sys j :

[0029]

[0030] Wherein, the subscript s represents a satellite, and the subscript sys represents a GNSS satellite system, where sys=1 is the GPS system and sys=2 is the BeiDou system;

[0031] Step S4) Calculate the vertical TEC value (TEC) directly above the ionospheric monitoring station using the following formula:

[0032]

[0033] Among them, b r,sys To accommodate receiver hardware discrepancies across various GNSS satellite systems, χ represents the satellite hardware deviation of satellite j. j The zenith angle of the line-of-sight path between satellite j and the receiver at the ionospheric puncture point.

[0034] Preferably, the TEC inversion calculation module includes a model building unit, and the spherical harmonic function model it builds satisfies the following equation:

[0035]

[0036] Where TEC represents the vertical TEC value directly above the ionospheric monitoring station, M and N represent the order and series of the spherical harmonic function, respectively, m and n represent the m-th and n-th series, respectively, and θ j Let φ be the geographical latitude of the point where the signal from satellite j penetrates the ionosphere. j The solar hour angle at the point where the signal of satellite j is punctured. For the normalized Legendre function, A nm and B nm Let be the spherical harmonic coefficients to be determined.

[0037] Preferably, the TEC inversion calculation module includes an iterative solution unit, and the processing includes:

[0038] Step T1) Construct the state equation and measurement equation for the Kalman filter;

[0039] Step T2) Calculate the zenith angle χ at the puncture point of the ionosphere for the signal of satellite j. j latitude θ j and longitude φ j ;

[0040] Step T3) Determine if this is the first Kalman filter calculation. If yes, initialize the parameter x0 to be estimated with 0 and the cross-covariance P0 of the parameter to be estimated with an initial value. If no, inherit the parameter x0 from the previous calculation. k-1 The cross-covariance P of the parameter to be estimated k-1 ;

[0041] Step T4) Calculate the estimated value of the parameter to be estimated. The covariance matrix P k ;

[0042] Step T5) Cache the estimated value of the parameter to be estimated. The covariance matrix P k ; and P k If the hardware deviation of satellite j is relatively and P k-1 middle If there are changes, update; if satellite j appears for the first time, then... and P k Add hardware deviations to satellite j. variable;

[0043] Step T6) The time sliding window slides backward by a fixed length, and the next round of Kalman filtering iteration begins based on the next set of preprocessed observations.

[0044] Preferably, the state equation for the Kalman filter in step T1) is:

[0045] x k =Φ k,k-1 x k-1 +W k

[0046] Where, Φ k,k-1 Let W be the state transition matrix from time (k-1) to time k. k Let x be the system noise matrix at time k. k Let the vector of parameters to be estimated at time k satisfy the following equation:

[0047]

[0048] Where the subscript k represents time; A 00,k A 10,k ,…A 33,k B 33,k b represents the 3rd order spherical harmonic coefficients to be calculated at time k; r，1,k b represents the receiver hardware deviation corresponding to the GPS signal at time k. r，2,kThis represents the receiver hardware deviation corresponding to the BeiDou system signal at time k. This represents the satellite hardware deviation numbered 1, 2, ..., J at time k, and the superscript T on the matrix indicates matrix transpose.

[0049] Preferably, the measurement equation for the Kalman filter in step T1) is:

[0050] z k =H k x k +V k

[0051] Among them, z k The oblique TEC measurement value after carrier phase smoothing pseudorange at time k. The set of j∈[1,J],

[0052]

[0053] H k The measurement matrix satisfies the following equation:

[0054]

[0055] V k For each satellite at time k Measurement noise.

[0056] Preferably, step T4) includes:

[0057] Step T4-1) Fuse the predicted values ​​at time k and measurement value z k The result of the parameter to be estimated at time k is obtained.

[0058]

[0059] in, To predict the estimate of the parameter vector at time k using the state equation, the following equation must be satisfied:

[0060]

[0061] in, The Kalman filter output of the parameter vector to be estimated at time k-1;

[0062] K k The Kalman gain at time k satisfies the following equation:

[0063]

[0064] R k =Cov(V k)

[0065]

[0066] Among them, R k To measure noise V k The covariance matrix, Let P be the estimation error covariance matrix of the parameter to be estimated at time k. k-1 Let Q be the error covariance matrix of the parameter to be estimated at time k-1. k The system noise W at time k k The covariance matrix satisfies the following equation:

[0067] Q k =Cov(W k )

[0068] Step T4-2) Update the error covariance matrix P at time k according to the following formula. k Used for estimating the error covariance matrix of the parameter to be estimated at time k+1:

[0069]

[0070] Where I is the identity matrix.

[0071] Preferably, the embedded processor is an ARM chip.

[0072] Compared with the prior art, the advantages of the present invention are:

[0073] 1. The ionospheric TEC monitoring instrument disclosed in this invention processes GNSS single-station data in the embedded processor chip of the GNSS receiver, solving the problem of ionospheric monitoring in areas where networking and data transmission are difficult.

[0074] 2. The ionospheric TEC monitor of the present invention is deployed in the GNSS dual-frequency receiver of the ionospheric monitoring station and is implemented based on the ARM chip, which solves the problem of limited equipment development environment and computing resources;

[0075] 3. By using a sliding time window data preprocessing method, spherical harmonic function modeling, and the proposed Klaman filter single-station quasi-real-time ionospheric vertical TEC inversion algorithm, this invention enables VTEC monitoring with a relatively small amount of data and achieves the required monitoring accuracy. Attached Figure Description

[0076] Figure 1 This is a block diagram showing the structural composition of the ionosphere monitoring instrument of the present invention;

[0077] Figure 2 Flowchart of the sliding time window data preprocessing method in the TEC inversion calculation module;

[0078] Figure 3 This is a flowchart of the Kalman filter TEC inversion calculation process in the TEC inversion calculation module;

[0079] Figure 4 The output results are from the ionosphere monitor using the present invention. Detailed Implementation

[0080] This invention discloses an ionospheric monitor, deployed within a GNSS dual-frequency receiver at an ionospheric monitoring station. It can be used for land and sea observations, and is particularly suitable for scenarios where communication bandwidth is limited, preventing the transmission of large amounts of raw data to a data center for processing. The ionospheric monitor includes:

[0081] GNSS antenna, used to receive GNSS satellite signals and transmit the received signals to the GNSS high-precision OEM board;

[0082] The GNSS high-precision OEM board is used to receive GNSS dual-frequency pseudorange observations, carrier phase observations, receiver position information, satellite elevation angle, azimuth angle and satellite ephemeris data according to the time window set, and transmit them to the TEC inversion calculation module and storage module respectively.

[0083] The TEC inversion calculation module, implemented based on an embedded processor, is used to preprocess the received data using a sliding time window, establish a spherical harmonic function model of the ionospheric vertical TEC, and solve iteratively using Kalman filtering to obtain the nonlinear model coefficients and the hardware deviation between the receiver and the satellite. This enables near real-time monitoring of the ionospheric vertical TEC value and transmits the TEC value to the storage module.

[0084] The storage module is used to store GNSS measurement data received by the GNSS high-precision OEM board, and also to store TEC values ​​output by the TEC inversion calculation module;

[0085] The wireless transmission module is used to periodically send data files from the storage module to the server.

[0086] The technical solution of the present invention will be described in detail below with reference to the accompanying drawings and embodiments.

[0087] Example 1

[0088] like Figure 1 As shown, the present invention discloses an ionospheric monitor, comprising: a TEC inversion calculation module, a GNSS high-precision OEM board, a storage module, a wireless transmission module, and a GNSS antenna.

[0089] The GNSS high-precision OEM board is used to measure the full-system, full-frequency signals of BDS, GPS, GLONASS, Galileo and QZSS. It outputs GNSS dual-frequency carrier phase and pseudorange observations, carrier-to-noise ratio, satellite elevation and azimuth data, and satellite ephemeris data, which are respectively provided to the TEC inversion calculation module and the TEC and GNSS measurement data file storage module.

[0090] The TEC inversion calculation module is used to calculate TEC results in near real-time based on measurement data from a high-precision GNSS OEM board. The data used includes high-precision dual-frequency pseudorange observations, carrier phase observations, and satellite ephemeris data from GNSS. The TEC inversion calculation module is implemented using an embedded processor; in one embodiment, an ARM chip is employed. A sliding time window data preprocessing method, nonlinear function TEC modeling, GNSS satellite hardware bias information inheritance, and a Kalman filtering algorithm are used to invert and calculate the vertical TEC of the ionosphere above the monitor.

[0091] The storage module stores TEC calculation result files and GNSS measurement data files. The data files in the storage module are periodically sent to the server.

[0092] The wireless transmission module is used to send TEC or GNSS measurement data files to the server. Each time a TEC result is calculated, it is sent to the data server via the wireless transmission module.

[0093] GNSS antennas are used to receive GNSS satellite signals and transmit the received signals to GNSS high-precision OEM boards.

[0094] When acquiring ionospheric TEC using observations of dual-frequency GNSS signals, both pseudorange-based and carrier phase-based observations can obtain oblique TEC along the line-of-sight path from the satellite to the receiver. Pseudorange-based observations have high noise, resulting in poor TEC accuracy; carrier phase-based observations have low noise, resulting in high TEC accuracy, but due to integer ambiguity in the carrier phase, the obtained TEC has an unknown bias. Typically, to obtain higher-accuracy TEC free of integer ambiguity, carrier phase smoothing pseudorange is used.

[0095] In the post-processing of ionospheric TEC data retrieved using GNSS observations, carrier phase smoothing pseudorange is performed using satellite observation information across the entire visible arc. This yields more reliable results. However, this method is not suitable for real-time ionospheric TEC monitoring; ionospheric TEC monitoring requires more real-time vertical TEC results.

[0096] Implementing ionospheric vertical TEC monitoring on an ARM processor requires consideration of ARM's storage space to conserve hardware resources as much as possible, while simultaneously obtaining reliable vertical TEC results. Therefore, a 20-minute data length is chosen for data preprocessing to obtain the oblique TEC after carrier phase smoothing pseudorange. However, this duration is not mandatory.

[0097] The core component, the TEC inversion calculation module, includes a preprocessing unit, a model building unit, and an iterative solution unit. These will be discussed in detail below.

[0098] Preprocessing unit:

[0099] like Figure 2 As shown, the sliding time window data preprocessing method includes the following steps:

[0100] Step 1: Receive and decode the ranging information from the full system dual-frequency pseudorange and carrier phase measurement output by the GNSS high-precision OEM board.

[0101] Step 2: Receive the elevation and azimuth information of GNSS system satellites output from the GNSS high-precision OEM board;

[0102] Step 3: Cache the observation time, satellite system, satellite number, elevation angle, azimuth angle, code pseudorange 1, code pseudorange 2, carrier phase measurement 1, and carrier phase measurement 2 within 20 minutes;

[0103] Step 4: Preprocess the 20 minutes of buffered data from Step 3, and calculate the oblique TEC based on dual-frequency pseudorange observations and the oblique TEC based on dual-frequency carrier phase observations for each satellite using the 20 minutes of observations:

[0104]

[0105] In the formula, the superscript j represents the satellite number; For the oblique TEC of satellite j based on dual-frequency pseudorange observations, The oblique TEC of satellite j is based on dual-frequency carrier phase observations; f1 and f2 are two frequencies of the GNSS satellite signal, respectively. These are the code pseudorange observations of the j-th satellite at two frequencies, f1 and f2, respectively. λ1 and λ2 are the carrier phase observations of the j-th satellite at two frequencies f1 and f2, respectively; λ1 and λ2 are the wavelengths corresponding to the two frequencies f1 and f2, respectively.

[0106] Step 5, calculate for satellite j. and Baseline values ​​between

[0107]

[0108] in For satellite j calculated using 20 minutes of data and The baseline value between; the superscript j indicates the satellite number, This represents the summation of 20 minutes of data, where i represents the measurement time of the i-th data point in the 20-minute data set, and Num represents the total number of times corresponding to the 20-minute data set. Represents the time of the i-th data measurement. Represents the time of the i-th data measurement. This represents the elevation angle of satellite j relative to the receiver at the measurement time of the i-th data point in the 20-minute data set.

[0109] Step 6: Calculate the slant TEC after carrier phase smoothing pseudorange, denoted as TEC. s,sys j :

[0110]

[0111] The superscript j represents the satellite number, and the subscript sys represents the GNSS satellite system. sys=1 represents the GPS system, and sys=2 represents the BeiDou system.

[0112] The oblique TEC (TEC) obtained after smoothing the pseudorange of the carrier phase of satellite j in step 6 s,sys j This includes hardware discrepancies between the receiver and the satellite. TEC s,sys j The vertical TEC (denoted as TEC) at the puncture point of the ionosphere by the satellite signal has the following relationship:

[0113]

[0114] That is, it can also be expressed as:

[0115]

[0116] Among them, b r,sys To accommodate receiver hardware deviations across various GNSS satellite systems, such as b r,1 For the receiver hardware deviation of the GPS system, b r,2 This is a hardware deviation in the BeiDou system receiver. χ represents the satellite hardware deviation of satellite j. j The zenith angle of satellite j is the angle between the line-of-sight path of the satellite and the receiver at the ionospheric puncture point and the zenith direction.

[0117] Model building unit:

[0118] In low-latitude regions, particularly in the Equatorial Ionization Anomaly (EIA), the vertical TEC of the ionosphere exhibits a nonlinear variation with latitude, forming a TEC maximum (the Hump of the EIA) near 20°N, and then rapidly decreasing north of the Hump. Therefore, using a nonlinear function to model TEC yields more accurate TEC inversion calculation results. When a spherical harmonic function is used for TEC modeling, the specific expression for the vertical TEC of the ionosphere above the station is as follows:

[0119]

[0120] Where TEC represents the vertical TEC value directly above the ionospheric monitoring station, M and N represent the order and series of the spherical harmonic function, respectively, m and n represent the m-th and n-th series, respectively, and θ j φ represents the geographical latitude of the point where the signal from satellite j penetrates the ionosphere, in radians. j The solar hour angle at the point where the signal of satellite j is punctured. For the normalized Legendre function, A nm and B nm Let be the spherical harmonic coefficients to be determined.

[0121] The main algorithms for single-station ionospheric inversion include the grid method, Kalman filtering method, and polynomial method. In the polynomial method, the polynomial can be a quadratic polynomial, a trigonometric series, or a spherical harmonic function. Generally, quadratic polynomials and trigonometric series are only suitable for inversions on a relatively small scale, while spherical harmonic functions, by changing their order and degree, can better describe the variations of the ionosphere at different spatial scales.

[0122] It should be noted that this embodiment uses a spherical harmonic function as an example, but it is not limited to this; it can also be a polynomial function, etc. For a spherical harmonic function, when a third-order spherical harmonic function is used, M = N = 3, but it is not limited to the third order; it can also be of other orders.

[0123] A nm and B nm The spherical harmonic coefficients can be calculated using the Kalman filtering method. The number of spherical harmonic coefficients is:

[0124] N SHC = (N+1) 2 -(NM)(N-M+1) (7)

[0125] In the Kalman filter calculation process, the parameter vector to be estimated at time k is defined as x. k for:

[0126]

[0127] Where the subscript k represents time; A 00,k A 10,k ,…A 33,k B 33,k b represents the 3rd order spherical harmonic coefficients to be calculated at time k; r，1,k b represents the receiver hardware deviation corresponding to the GPS signal at time k. r，2,k This represents the receiver hardware deviation corresponding to the BeiDou system signal at time k. This represents the satellite hardware bias numbered 1, 2, ..., J at time k. Kalman filtering is used to calculate x. k The state equation is:

[0128] x k =Φ k,k-1 x k-1 +W k (9)

[0129] Where, Φ k,k-1 The state transition matrix from time (k-1) to time k is shown below:

[0130]

[0131] W k Here is the system noise matrix at time k:

[0132]

[0133] W k Each value in the vector represents the corresponding parameter to be estimated, x. k The noise variance at time k.

[0134] According to formulas (5) and (6), the Kalman filter calculates x. k The measurement equations during the process are expressed as follows:

[0135] z k =H k x k +V k (12)

[0136] Among them, z k The oblique TEC measurement value after carrier phase smoothing pseudorange at time k. According to formula (5):

[0137]

[0138] According to formulas (5) and (6), H in the measurement equation k as follows:

[0139]

[0140] In the formula, V k For each satellite at time k Measurement noise.

[0141] Assuming measurement noise V k and system noise W k All are zero-mean white noise sequences. V k The covariance matrix is ​​taken as:

[0142] R k =Cov(V k ) = 0.01I J (15)

[0143] Among them, I J W is a J×J identity matrix. k The covariance matrix is ​​taken as:

[0144]

[0145] in For (N) SHC +2+J)×(N SHC +2+J) is the identity matrix; is the total number of spherical harmonic coefficients of M and N order spherical harmonic functions, refer to formula (7); 2 represents the use of satellite data from both GPS and Beidou systems; J is the total number of satellites used.

[0146] Iterative solution unit:

[0147] The recursive process of Kalman filtering calculation can be divided into two steps: the first step is prediction, and the second step is to correct the prediction result using measurement information.

[0148] The first step of prediction involves using the state equation to predict the estimate of the parameter vector at time k, denoted as .

[0149]

[0150] Where, Φ k,k-1 For the state transition matrix, refer to formulas (9) and (10); This is the Kalman filter output of the parameter vector to be estimated at time k-1. It also includes the estimation error covariance matrix for predicting time k.

[0151]

[0152] P in the formula k-1 Let Q be the estimation error covariance matrix at time k-1. k Refer to formula (16).

[0153] The second step involves using measurement information to correct the prediction results, including using the measurement matrix, the measurement covariance matrix, and the Kalman gain K at time k defined by equation (18). k :

[0154]

[0155] H in the formula k Referring to formula (14), R k Refer to formula (15).

[0156] Fusion of predicted values ​​at time k and measurement value z k The results of the parameters to be estimated at time k are obtained:

[0157]

[0158] Update the estimated error covariance matrix for the next estimation, i.e.:

[0159]

[0160] In the formula, I is the identity matrix.

[0161] If this is the first Kalman filter calculation, initialize the parameter to be estimated to 0; or if the statistical characteristics of the initial state x0 are known, initialize the parameter to be estimated as follows:

[0162]

[0163] Using equations (17) to (22), the Kalman filtering process can converge quickly, thereby estimating the vertical TEC and hardware bias.

[0164] like Figure 3 As shown, the GNSS satellite hardware bias information inheritance and Kalman filtering algorithm of the iterative solution unit includes the following steps:

[0165] Step 1: Construct the state equation and measurement equation for the Kalman filter;

[0166] Step 2: Calculate the zenith angle χ at the ionospheric puncture point. j latitude θ j and longitude φ j ,j=1,2,…,J is the satellite number at the current epoch.

[0167] Step 3: Determine if this is the first time using the Kalman filter.

[0168] Step 4: If yes, initialize the parameter to be estimated x0 to 0 and the cross-covariance P0 of the parameter to be estimated to 0; otherwise, inherit the parameter to be estimated from the previous step. k-1The cross-covariance P of the parameter to be estimated k-1 ;

[0169] Step 5: Calculate the covariance matrix R of the measurement noise and the system noise. k and Q k ;

[0170] Step 6, calculate the state transition matrix Φ k,k-1 and measurement matrix H k ;

[0171] Step 7: Calculate the covariance matrix of the prediction error. and Kalman filter gain matrix K k

[0172] Step 8: Calculate the estimated value of the parameter to be estimated. The covariance matrix P k ;

[0173] Step 9, cache the estimated value The covariance matrix P k ; and P k If the hardware deviation of satellite j is relatively and P k-1 middle If there are changes, update; if satellite j appears for the first time, then... and P k Add hardware deviations to satellite j. variable;

[0174] Step 10: When the time sliding window slides backward by 5 minutes, re-enter Step 2 based on the next set of preprocessed observations to start the next round of Kalman filtering iteration.

[0175] Due to the movement of GNSS satellites, received GNSS satellite signals may fade out of the receiver's range, and satellites may re-enter the receiver's range. Satellite hardware bias is among the parameters to be estimated. The number of satellites varies with the number of visible satellites and is not a fixed vector. The parameter x to be estimated in the Kalman filtering process... k The intermediate calculation matrix must also change along with the changes in the observed values.

[0176] Because of parameters to be estimated This refers to hardware deviations in the satellite, which generally change very little. After the satellite fades out of the observation area and then re-enters the observation area, the changes in its hardware deviations become more pronounced. The initial values ​​are inherited from the results retained before the fade-out. Similarly, the cross-correlation matrix of the observations also inherits the initial values ​​retained before the fade-out, serving as the initial values ​​for this satellite's entry. This avoids the instability caused by multiple initializations during the Kalman filtering process, ensuring the reliability of the results.

[0177] like Figure 4 This paper presents the results of a single-station quasi-real-time vertical TEC inversion algorithm using Klaman filtering for TEC modeling with spherical harmonics, as implemented in this invention. The results are compared with those of a single-station algorithm using Klaman filtering for TEC modeling with a vertical TEC reference value and a linear function modeling TEC. In the figure, the TEC reference value is the result of inverting TEC using spherical harmonics on a full day's measurement data, while Klaman+linear represents the result of a single-station algorithm using Klaman filtering for TEC modeling with a linear function. The results show that the results obtained using this invention are closer to the reference results.

[0178] In summary, the ionospheric monitor provided by this invention achieves near real-time single-station inversion calculation of ionospheric TEC by using a sliding time window data preprocessing method, nonlinear function TEC modeling, GNSS satellite hardware deviation information inheritance, and Kalman filtering algorithm.

[0179] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit it. Although the present invention has been described in detail with reference to the embodiments, those skilled in the art should understand that modifications or equivalent substitutions to the technical solutions of the present invention do not depart from the spirit and scope of the technical solutions of the present invention, and all such modifications or substitutions should be covered within the scope of the claims of the present invention.

Claims

1. An ionospheric monitor, characterized in that, The GNSS dual-frequency receiver deployed at the ionospheric monitoring station includes: The GNSS high-precision OEM board is used to receive GNSS dual-frequency pseudorange observations, carrier phase observations, receiver position information, satellite elevation angle, azimuth angle and satellite ephemeris data according to the time window set, and transmit them to the TEC inversion calculation module and storage module respectively. The TEC inversion calculation module, implemented using an embedded processor, preprocesses the received single-station data using a sliding time window, establishes a spherical harmonic function model or polynomial model of the ionospheric vertical TEC, and iteratively solves it using Kalman filtering to obtain the nonlinear model coefficients and the hardware deviations between the receiver and the satellite. This enables near real-time monitoring of the ionospheric vertical TEC values ​​above the monitoring station, and the TEC values ​​are then transmitted to the storage module. The storage module is used to store GNSS measurement data received by the GNSS high-precision OEM board, and also to store the ionospheric vertical TEC value output by the TEC inversion calculation module. The TEC inversion calculation module includes a preprocessing unit, and the processing includes: Step S1) Obtain the oblique TEC value of satellite j based on dual-frequency pseudorange observations within the time window set by the cached data. and the oblique TEC value of satellite j based on dual-frequency carrier phase observations Step S2) Calculate the value of satellite j according to the following formula. and Baseline values ​​between Where i represents the measurement time of the i-th data in the data received within the time window set duration, and Num represents the total number of times corresponding to the data within the time window set duration; This represents the elevation angle of satellite j relative to the receiver at the i-th data measurement time. Step S3) According to and Calculate the oblique TEC value (TEC) of satellite j after carrier phase smoothing pseudorange. s,sys j ; Step S4) According to TEC s,sys j By combining hardware deviation and zenith angle, the vertical TEC value TEC above the ionospheric monitoring station is calculated. The TEC inversion calculation module includes an iterative solution unit, and the processing procedure includes: Step T1) Construct the state equation and measurement equation for the Kalman filter; Step T2) Calculate the zenith angle χ at the puncture point of the ionosphere for the signal of satellite j. j latitude θ j and longitude φ j ; Step T3) Determine if this is the first Kalman filter calculation. If yes, initialize the parameter x0 to be estimated with 0 and the cross-covariance P0 of the parameter to be estimated with an initial value. If no, inherit the parameter x0 from the previous calculation. k-1 The cross-covariance P of the parameter to be estimated k-1 ; Step T4) Calculate the estimated value of the parameter to be estimated. The covariance matrix P k ; Step T5) Cache the estimated value of the parameter to be estimated. The covariance matrix P k ; and P k If the hardware deviation of satellite j is relatively and P k-1 middle If there are changes, update; if satellite j appears for the first time, then... and P k Add hardware deviations to satellite j. variable; Step T6) The time sliding window slides backward by a fixed length, and the next round of Kalman filtering iteration begins based on the next set of preprocessed observations.

2. The ionospheric monitor according to claim 1, characterized in that, Also includes: GNSS antenna, used to receive GNSS satellite signals and transmit the received signals to the GNSS high-precision OEM board; and The wireless transmission module is used to periodically send data files from the storage module to the server.

3. The ionospheric monitor according to claim 1, characterized in that, The processing of the GNSS high-precision OEM board includes: According to the time window set, receive and decode ranging information from multiple GNSS satellites using dual-frequency pseudorange and carrier phase measurements. Receive elevation and azimuth information for each satellite; The ranging information for each satellite, including its number, observation time, elevation angle, azimuth angle, GNSS dual-frequency pseudorange, and carrier phase measurement, is cached separately.

4. The ionospheric monitor according to claim 1, characterized in that, In step S1), the oblique TEC value of satellite j based on dual-frequency pseudorange observations. and the oblique TEC value of satellite j based on dual-frequency carrier phase observations We obtain the following formula: Where f1 and f2 are two frequencies of the GNSS satellite signal, These are the code pseudorange observations of satellite j at two frequencies f1 and f2, respectively. λ1 and λ2 are the carrier phase observations of the j-th satellite at two frequencies f1 and f2, respectively, where λ1 and λ2 are the wavelengths corresponding to the two frequencies f1 and f2, respectively. The oblique TEC value (TEC) after smoothing the pseudorange of the carrier phase of the j-th satellite in step S3) s,sys j We obtain the following formula: Wherein, the subscript s represents a satellite, and the subscript sys represents a GNSS satellite system, where sys=1 is the GPS system and sys=2 is the BeiDou system; In step S4), the vertical TEC value (TEC) above the ionospheric monitoring station is obtained according to the following formula: Among them, b r,sys To accommodate receiver hardware discrepancies across various GNSS satellite systems, χ represents the satellite hardware deviation of satellite j. j To calculate the zenith angle of the signal from satellite j at the puncture point in the ionosphere.

5. The ionospheric monitor according to claim 1, characterized in that, The TEC inversion calculation module includes a model building unit, which, when building a spherical harmonic function model, satisfies the following equation: Where TEC represents the vertical TEC value above the ionospheric monitoring station, M and N represent the order and series of the spherical harmonic function, respectively, m and n represent the m-th and n-th series, respectively, and θ j Let φ be the geographical latitude of the point where the signal from satellite j penetrates the ionosphere. j The solar hour angle at the point where the signal of satellite j is punctured. For the normalized Legendre function, A nm and B nm Let be the spherical harmonic coefficients to be determined.

6. The ionospheric monitor according to claim 1, characterized in that, The state equation for the Kalman filter in step T1) is: x k =Φ k,k-1 x k-1 +W k Where, Φ k,k-1 Let W be the state transition matrix from time (k-1) to time k. k Let x be the system noise matrix at time k. k Let the vector of parameters to be estimated at time k satisfy the following equation: Where the subscript k represents time; A 00,k A 10,k ,…A 33,k B 33,k b represents the 3rd order spherical harmonic coefficients to be calculated at time k; r，1,k b represents the receiver hardware deviation corresponding to the GPS signal at time k. r，2,k This represents the receiver hardware deviation corresponding to the BeiDou system signal at time k. This represents the satellite hardware deviation numbered 1, 2, ..., J at time k, and the superscript T on the matrix indicates matrix transpose.

7. The ionospheric monitor according to claim 6, characterized in that, The measurement equation for the Kalman filter in step T1) is: z k =H k x k +V k Among them, z k The oblique TEC measurement value after carrier phase smoothing pseudorange at time k. The set of j∈[1,J], H k The measurement matrix satisfies the following equation: V k For each satellite at time k Measurement noise.

8. The ionospheric monitor according to claim 7, characterized in that, Step T4) includes: Step T4-1) Fuse the predicted values ​​at time k and measurement value z k The result of the parameter to be estimated at time k is obtained. in, To predict the estimate of the parameter vector at time k using the state equation, the following equation must be satisfied: in, The Kalman filter output of the parameter vector to be estimated at time k-1; K k The Kalman gain at time k satisfies the following equation: R k =Cov(V k ) Among them, R k To measure noise V k The covariance matrix, Let P be the estimation error covariance matrix of the parameter to be estimated at time k. k-1 Let Q be the error covariance matrix of the parameter to be estimated at time k-1. k The system noise W at time k k The covariance matrix satisfies the following equation: Q k =Cov(W k ) Step T4-2) Update the error covariance matrix P at time k according to the following formula. k Used for estimating the error covariance matrix of the parameter to be estimated at time k+1: Where I is the identity matrix.

9. The ionospheric monitor according to claim 1, characterized in that, The embedded processor is an ARM chip.

Images