A robust method for GNSS regional ionospheric TEC modeling and differential code bias estimation

Abstract

This invention relates to the field of ionospheric remote sensing and space environment monitoring technology for Global Navigation Satellite Systems (GNSS), specifically to a robust estimation method for GNSS regional ionospheric TEC modeling and differential code bias. By constructing an ionospheric geometry-independent observation model, an adaptive variable altitude mapping function and spherical cap harmonics are introduced to refine the regional ionospheric VTEC model. Furthermore, QR orthogonal transformation is used to decouple the estimation of ionospheric parameters and differential code bias parameters. Simultaneously, Tikhonov regularization constraints are introduced during the ionospheric parameter solution stage, thereby achieving stable, high-precision, and continuous reconstruction of the regional TEC. This invention reduces mapping errors through an adaptive mapping function, eliminates boundary effects using spherical cap harmonics, improves the reliability of DCB estimation through QR decomposition, and finally ensures spatiotemporal continuity through Tikhonov regularization. The corresponding structure and algorithm are coupled and synergistic, simultaneously solving key bottleneck problems in traditional methods from both physical modeling and numerical solution perspectives.

Description

Technical Field

[0001] This invention relates to the field of ionospheric remote sensing and space environment monitoring technology for global navigation satellite systems, specifically to a robust method for GNSS regional ionospheric TEC modeling and differential code bias estimation. Background Technology

[0002] The ionosphere is a region of charged plasma in the Earth's upper atmosphere formed by ionization from solar radiation, and its electron density distribution exhibits significant variations in both time and space. Free electrons in the ionosphere cause frequency-dependent dispersion delays in the propagation of Global Navigation Satellite System (GNSS) signals, directly affecting the accuracy and reliability of GNSS positioning, navigation, and timing (PNT). Therefore, retrieving the total electron content (TEC) of the ionosphere using dual-frequency GNSS observation data and establishing a high-precision ionospheric model has become an important technical means for space weather monitoring, ionospheric dynamics research, and GNSS error correction.

[0003] Currently, GNSS observation-based ionospheric modeling methods typically employ a single-layer thin-shell model (SLM). This model converts the slant total electron content (STEC) along the satellite-to-receiver line-of-sight (STS) into the vertical total electron content (VTEC) using a mapping function, and then spatially parameterizes it using spherical harmonics (SH) or polynomial functions. However, traditional methods still have several shortcomings in regional-scale ionospheric modeling. First, the fixed-height thin-shell assumption usually sets the equivalent ionospheric height to approximately 450 km, while the actual peak height of the F2 layer (hmF2) varies significantly under different local times, seasons, and solar activity conditions, leading to substantial mapping errors in low-elevation-angle observations. Second, standard spherical harmonics are primarily suitable for global-scale modeling. When applied to regional CORS networks, due to the limited observation coverage, they are prone to significant oscillations at regional boundaries (i.e., the "Runge phenomenon"), thus reducing the stability and accuracy of regional TEC reconstruction.

[0004] Furthermore, in the GNSS ionospheric inversion process, the observation equations not only include ionospheric TEC parameters but also receiver-satellite hardware delay errors, i.e., Differential Code Bias (DCB). Ionospheric parameters and DCB parameters are highly correlated in the observation model; without proper handling, this often leads to severe rank deficiency in the normal equations, thus affecting the stability and reliability of parameter estimation. Traditional methods typically address this problem by setting reference satellites or constraints, but under complex ionospheric disturbances or poor observation conditions, unstable solutions or amplified errors are still prone to occur. Summary of the Invention

[0005] The purpose of this invention is to provide a robust estimation method for GNSS regional ionospheric TEC modeling and differential code bias, which is used to reduce ionospheric mapping errors, suppress regional modeling boundary effects, and achieve stable decoupling of ionospheric spatial parameters and GNSS hardware delay parameters in regional ionospheric modeling. At the same time, it improves the robustness of differential code bias estimation and the numerical stability of TEC inversion results under complex space environment conditions, thereby achieving high-precision continuous reconstruction of regional ionospheric TEC.

[0006] To achieve the above objectives, this invention provides a method for GNSS region ionospheric TEC modeling and robust estimation of differential code bias, comprising the following steps:

[0007] Step 1: Construct a GNSS ionospheric geometry-independent observation model;

[0008] Step 2: Establish an adaptive mapping function that takes into account the variation of the ionospheric peak height;

[0009] Step 3: Construct a regional ionospheric VTEC model based on spherical cap harmonic analysis;

[0010] Step 4: Construct the network observation equations and perform linearization;

[0011] Step 5: Use QR orthogonal transformation to achieve parameter decoupling;

[0012] Step 6: Perform a global robust estimation of the differential code bias parameters;

[0013] Step 7: Perform regularized reconstruction of ionospheric parameters;

[0014] Step 8: Output the regional ionosphere TEC reconstruction results.

[0015] Optionally, in step 1, a geometrically independent combination is constructed using GNSS dual-frequency pseudorange observations and carrier phase observations, and observation noise and multipath effects are reduced by carrier phase smoothing pseudorange technology to obtain high-precision geometrically independent ionospheric observations, thereby establishing an ionospheric observation equation that includes STEC and the differential code deviation between the receiver and the satellite.

[0016] Optionally, the adaptive mapping function converts STEC to VTEC based on the dynamic change of the peak height of the F2 layer of the ionosphere, wherein the shell height is provided with a priori constraints by the international reference ionospheric model and adaptively changes with the geomagnetic latitude of the ionospheric puncture point, the fixed solar longitude, and the observation epoch.

[0017] Optionally, in step 3, VTEC is expanded into a series by constructing a generalized connected Legendre function basis that satisfies the boundary conditions within the spherical cap region.

[0018] Optionally, in step 4, the adaptive mapping function and the regional ionospheric VTEC model based on spherical cap harmonic analysis are substituted into the geometry-independent observation model to uniformly model all receiver and visible satellite observation data in the regional GNSS network, forming a batch observation equation that includes ionospheric cap harmonic analysis coefficients and differential code deviation parameters, and is expressed in matrix form as a network-level linear observation model.

[0019] Optionally, in step 5, the network-level linear observation model is decomposed into two independent subsystems by QR orthogonal transformation, one of which contains only the differential code bias parameter.

[0020] Optionally, in step 6, the differential code deviation subsystems of multiple time batches are accumulated to construct a global normal equation. The parameter rank deficiency problem is solved by introducing a zero-sum constraint condition for satellite differential code deviation. A robust estimate of the receiver-satellite differential code deviation is obtained by using the constrained least squares method.

[0021] Optionally, step 7 introduces Tikhonov regularization constraints during the ionospheric parameter solution process, describes the continuous change of ionospheric state between adjacent time batches through a smooth difference operator, and adaptively determines the regularization parameters using the L-curve method.

[0022] Optionally, in step 8, the VTEC distribution at any ionospheric puncture point within the study area is calculated based on the obtained SCHA coefficients.

[0023] This invention provides a robust method for modeling regional ionospheric TEC (Vibrational Temperature Coefficient) and differential code bias estimation in GNSS regions. By constructing an ionospheric geometry-independent observation model, an adaptive variable altitude mapping function and spherical cap harmonics are introduced to refine the modeling of regional ionospheric VTEC. Furthermore, QR orthogonal transformation is used to decouple the estimation of ionospheric parameters and differential code bias parameters. Simultaneously, Tikhonov regularization constraints are introduced during the ionospheric parameter solution stage, thereby achieving stable, high-precision, and continuous reconstruction of regional TEC. This invention reduces mapping errors through an adaptive mapping function, eliminates boundary effects using spherical cap harmonics, improves the reliability of DCB (Discretionary Code Bias) estimation through QR decomposition, and finally ensures spatiotemporal continuity through Tikhonov regularization. The corresponding structure and algorithm are coupled and work synergistically, simultaneously solving key bottleneck problems in traditional methods from both physical modeling and numerical solution perspectives. Attached Figure Description

[0024] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0025] Figure 1 This is a schematic diagram of the overall process for fine modeling of the regional ionosphere in an embodiment of the present invention.

[0026] Figure 2 This is a schematic diagram of the VTEC extraction experiment results in an embodiment of the present invention.

[0027] Figure 3 This is a schematic diagram of the zero-baseline experimental results in an embodiment of the present invention.

[0028] Figure 4 This is a schematic diagram of the satellite DCB estimation experiment results in an embodiment of the present invention.

[0029] Figure 5 This is a schematic diagram of the experimental results of satellite DCB stability analysis in an embodiment of the present invention.

[0030] Figure 6 This is a schematic diagram of the experimental results of DCB stability analysis at the station in an embodiment of the present invention.

[0031] Figure 7 This is a schematic diagram of the VTEC modeling experiment results in an embodiment of the present invention.

[0032] Figure 8 This is a schematic diagram of the VTEC error analysis experimental results in an embodiment of the present invention.

[0033] Figure 9This is a schematic diagram of a comparative experiment between the method described in this invention and a traditional method in an embodiment of the invention. Detailed Implementation

[0034] Embodiments of the present invention are described in detail below, examples of which are illustrated in the accompanying drawings, wherein the same or similar reference numerals denote the same or similar elements or elements having the same or similar functions throughout. The embodiments described below with reference to the accompanying drawings are exemplary and intended to explain the present invention, and should not be construed as limiting the present invention.

[0035] This invention provides a robust method for ionospheric TEC modeling and differential code bias estimation in GNSS regions, comprising the following steps:

[0036] Step 1: Construct a GNSS ionospheric geometry-independent observation model;

[0037] Step 2: Establish an adaptive mapping function that takes into account the variation of the ionospheric peak height;

[0038] Step 3: Construct a regional ionospheric VTEC model based on spherical cap harmonic analysis;

[0039] Step 4: Construct the network observation equations and perform linearization;

[0040] Step 5: Use QR orthogonal transformation to achieve parameter decoupling;

[0041] Step 6: Perform a global robust estimation of the differential code bias parameters;

[0042] Step 7: Perform regularized reconstruction of ionospheric parameters;

[0043] Step 8: Output the regional ionosphere TEC reconstruction results.

[0044] This invention introduces an adaptive variable-height mapping function (VHMF) and spherical cap harmonic analysis (SCHA) to perform refined modeling of the vertical total electron content (VTEC) of the regional ionosphere. The following description, combined with the execution process and specific embodiments, provides further explanation:

[0045] Step 1: Construct a GNSS ionospheric geometry-independent observation model

[0046] When GNSS signals pass through the ionosphere, their propagation delay exhibits significant frequency dispersion. By linearly combining dual-frequency pseudorange and carrier phase observations and employing carrier-to-code smoothing pseudorange technology, most of the non-dispersion errors related to the geometry path can be eliminated, extracting high-precision geometry-independent smoothed code observations.

[0047] (1)

[0048] in, This represents geometrically independent observations after phase smoothing; For satellite To the receiver Slant TEC (Slant Total Electron Content) is the integral along the line of sight. and These represent the differential code offset (DCB) between the receiver and the satellite, respectively. Residual multipath effects and observation noise; The ionospheric dispersion proportionality coefficient is defined as follows: ,in and This corresponds to the carrier frequency.

[0049] Step 2: Establish an adaptive mapping function that takes into account the variation of ionospheric peak height.

[0050] To normalize STEC along different geometric paths in space to the vertical total electron content (VTEC), traditional models often use a fixed height (e.g., The assumption of a single-layer thin shell. However, the peak height of the F2 ionosphere ( The altitude fluctuates dramatically depending on local time and solar activity conditions. To mitigate the mapping error introduced by a fixed altitude in low-altitude observations, this invention constructs an adaptive variable-height mapping function (VHMF):

[0051] (2)

[0052] (3)

[0053] In the formula, The elevation angle of the satellite at the station; The average radius of the Earth; An empirical correction factor. Dynamic thin-shell height. Prior constraints are provided by the International Reference Ionospheric Model (IRI), depending on the geomagnetic latitude of the ionospheric puncture point (IPP). Sun's fixed longitude and observation epochs Adaptive evolution.

[0054] Step 3: Construct a regional ionospheric VTEC model based on spherical cap harmonic analysis

[0055] Due to the spatially limited nature of regional CORS network observation data, the standard spherical harmonic function (SH) is highly susceptible to ill-conditioned normal equations and the "Runge phenomenon" in boundary regions. This invention employs Spherical Cap Harmonic Analysis (SCHA) to rigorously parameterize the regional VTEC. The study area is limited to a half-angle of... Within the spherical cap, VTEC operates in the local spherical cap coordinate system. Expand below:

[0056] (4)

[0057] in, The maximum truncation order; and The SCHA coefficient to be estimated; This is the generalized normalized associated Legendre function. Unlike the standard SH function, its order is... For non-integer values, it is necessary to make the basis functions on the spherical cap boundary ( The characteristic equations satisfying specific Dirichlet or Neumann boundary conditions are obtained by solving them, thereby ensuring the orthogonality and high spatial resolution of the region modeling.

[0058] Step 4: Construct the network observation equations and perform linearization.

[0059] Substituting equations (2) and (4) into equation (1), a rigorous nonlinear ionospheric observation model for a single line-of-sight vector can be established. In the first... Within a time batch (usually set to 1 to 2 hours to ensure DCB stability), its linearization is simplified as follows:

[0060] (5)

[0061] in, To include all SCHA coefficients to be estimated The state vector; The row vectors of the design matrix are calculated based on the mapping function and the SCHA basis.

[0062] Traversing the first A batch-level linear observation model is constructed by using observation records from all receivers and visible satellites within a given time batch in the regional network.

[0063] (6)

[0064] in, This is the phase-smoothed pseudorange observation vector; Dense design matrix for SCHA parameters of the ionosphere; This is a global hardware delay parameter vector that includes all receivers and satellite DCBs in the entire observation network (assuming it remains constant over several consecutive days). This is the corresponding DCB mapping matrix (most elements are 0 or 1); This is the observed noise vector.

[0065] Step 5: Use QR orthogonal transformation to achieve parameter decoupling

[0066] To solve the problem of ionospheric spatial coefficient With global hardware latency parameter vector The strong correlation between them, and the large design matrix in this invention Perform QR orthogonal decomposition:

[0067] (7)

[0068] in, It is an orthogonal matrix, divided into blocks. and ; It is a non-singular upper triangular matrix. According to the null space property of QR decomposition, it must satisfy:

[0069] (8)

[0070] orthogonal transformation matrix Multiplying both ends of equation (6) on the left, the original highly coupled system can be rigorously transformed into two independent subsystems:

[0071] (9)

[0072] In the formula, ; ; Thus, a pure DCB observation subsystem was created, eliminating the influence of ionospheric parameters. They were successfully isolated.

[0073] Step 6: Perform a global robust estimation of the differential code bias parameter.

[0074] Due to the separated subsystems It no longer includes highly dynamic ionospheric parameters, and can be batched across all time periods over multiple consecutive days (e.g., 3-5 days). Accumulate data to construct a global normal equation:

[0075] (10)

[0076] To address the inherent rank deficiency problem in GNSS hardware delay calculation, a baseline constraint is introduced that the sum of the DCBs of all visible satellites is zero. Thus, the optimal estimate of constrained least squares is obtained:

[0077] (11)

[0078] Step 7: Perform regularized reconstruction of ionospheric parameters

[0079] Accurately obtaining network DCB parameters Afterwards, the traditional method directly used back substitution The ionospheric coefficients are solved independently for each epoch. However, during nighttime hours or in sparsely populated edge regions of the station, The ill-conditioned nature of the solution can lead to divergence. Considering the continuity of ionospheric state evolution, this paper introduces a Tikhonov regularization penalty function in the back-substitution step:

[0080] (12)

[0081] in, The prior state coefficients of the previous time batch; It is a smoothing difference operator matrix; The regularization smoothing parameter is dynamically determined by the L-curve method based on the weights between the observation residuals and the smoothing constraints. This regularization mechanism can effectively absorb numerical oscillations caused by observation blind spots, ultimately achieving robust and high-precision continuous reconstruction of regional TEC.

[0082] For further information, please refer to [link / reference]. Figures 1 to 9 The present invention illustrates and verifies the key steps of the method through four embodiments:

[0083] Example 1: Ionospheric observation modeling based on geometrically independent observations and VHMF mapping

[0084] like Figure 1As shown, this embodiment first constructs a geometry-independent GNSS ionospheric observation model. By linearly combining the dual-frequency pseudorange and carrier phase observations and introducing carrier phase smoothing pseudorange technology, multipath effects and thermal noise are effectively suppressed, thereby obtaining high-precision oblique total electron content (STEC) observations. This observation model utilizes the dispersion characteristics of the ionosphere to signals of different frequencies from a physical mechanism perspective, separating the ionospheric delay from the geometric term to achieve direct and sensitive observation of electron content.

[0085] Building upon this, an adaptive variable height mapping function (VHMF) is introduced. By using the peak height (hmF_2) of the F2 layer corresponding to the ionospheric puncture point (IPP) as a dynamic parameter, and incorporating prior information provided by the international reference ionospheric model, the mapping function dynamically varies with geomagnetic latitude, local solar time, and observation epoch. Compared to traditional fixed-height thin-shell models, this method significantly reduces mapping error under low-height-angle observation conditions. Figure 9 .

[0086] From a working principle perspective, traditional single-layer models assume that the electron density of the ionosphere is concentrated at a fixed height, leading to a systematic deviation in the STEC to VTEC conversion under long slant-range paths (especially at low elevation angles). This invention, by introducing dynamics (hmF2), allows the equivalent ionospheric shell height to accurately reflect the evolution of ionospheric thermodynamic and photochemical processes, thus making the mapped path closer to the actual propagation path. This structural improvement reduces the STEC to VTEC conversion error by approximately 20%–35%, significantly improving modeling accuracy. Figure 3 .

[0087] Example 2: A Refined Modeling Method for Regional VTEC Based on Spherical Harmonic Analysis (SCHA)

[0088] like Figure 7 , 8 As shown, this embodiment addresses the issues of uneven distribution of observation stations and boundary effects in regional ionospheric modeling by employing a spherical cap harmonic analysis method to parametrically represent VTEC. By restricting the study area to a specific spherical cap and constructing an orthogonal basis function system that satisfies the boundary conditions, high-resolution reconstruction of the electron content of the regional ionosphere is achieved.

[0089] Specifically, this method solves for non-integer order associated Legendre functions to satisfy Dirichlet or Neumann conditions at the spherical cap boundary, thereby ensuring the stability and orthogonality of the basis functions at the region boundary. Compared to traditional spherical harmonic function expansions, it avoids the "Runge phenomenon" caused by higher-order expansions and effectively suppresses ill-conditioned problems in the normal equations.

[0090] From a working mechanism perspective, after regional truncation, the basis functions of traditional global spherical harmonic functions are no longer orthogonal at the boundaries, leading to fitting oscillations and extrapolation errors. SCHA, however, reconstructs local orthogonal bases, ensuring the model relies solely on observational data within the region, thus avoiding the introduction of external spurious information. Therefore, this method reduces VTEC reconstruction errors by over 30% at regional edges (such as land-sea boundaries or sparsely populated areas), significantly improving spatial consistency and physical plausibility.

[0091] Example 3: Decoupling of Ionospheric Parameters and DCB Based on QR Orthogonal Decomposition

[0092] like Figure 4 , 5 As shown in Figure 6, this embodiment proposes to use QR orthogonal decomposition to perform structural transformation on the observation equations, achieving strict decoupling between ionospheric spatial parameters and hardware delay parameters (DCB). By performing orthogonal decomposition on the design matrix, the original highly coupled observation system is transformed into two independent subsystems, one containing only DCB parameters and the other containing ionospheric parameters.

[0093] Its core principle lies in utilizing the null space properties after QR decomposition to construct an orthogonal complement space, thereby eliminating the projection of ionospheric parameters and obtaining a pure DCB observation subsystem. Essentially, this process involves subspace decomposition of observational information, enabling orthogonal separation of parameters from different physical sources in mathematical space.

[0094] This method addresses the rank deficiency and solution instability issues caused by strong parameter correlation in traditional joint estimation. It constructs a global normal equation by accumulating data from multiple days and introduces a zero-mean constraint, effectively restoring the identifiability of DCB parameters. Experimental results show that, compared to traditional least-squares joint estimation methods, this method improves DCB estimation stability by approximately 40% and significantly reduces the impact of solar ionization disturbances on hardware delay estimation.

[0095] Example 4: Robust Reconstruction of Ionospheric Parameters with Tikhonov Regularization

[0096] like Figure 7 As shown, after completing the DCB solution, this embodiment further proposes to introduce Tikhonov regularization constraints in the ionospheric parameter back-substitution solution stage. By constructing a time-difference smoothing term, the ionospheric states of adjacent time batches are linked, thereby suppressing solution oscillations caused by observation sparsity or geometric weakening.

[0097] Its working principle is as follows: the spatiotemporal evolution of ionospheric electron density has continuity and gradual variation. By introducing a smoothing operator matrix, a physically reasonable constraint is imposed on the solution space, so that the solution results maintain temporal continuity while satisfying the observation fit. The regularization parameter is adaptively selected through the L-curve method to achieve the optimal balance between observation error and model smoothing.

[0098] This method effectively solves the ill-conditioned problem of the normal equations caused by insufficient observational information under low electron density conditions in the ionosphere at night. In simulation experiments, when the number of observed satellites is reduced by 30%, the traditional method's VTEC solution error can reach 5–8 TECU, while this method can control it within 2–3 TECU, significantly improving the robustness of the solution under extreme conditions.

[0099] In summary, compared with the prior art, the present invention has the following beneficial effects:

[0100] 1. Improved accuracy of ionospheric projection modeling. By introducing an adaptive variable height mapping function (VHMF) that takes into account the dynamic changes in the peak height of the F2 layer (hmF2), the STEC–VTEC projection error generated by the traditional fixed thin-shell height model under low elevation angle observation conditions is effectively overcome. This method can adaptively adjust the thin-shell height according to the geomagnetic latitude, fixed solar longitude, and local time variation of the ionospheric puncture point (IPP), thereby significantly improving the accuracy of regional vertical total electron content (VTEC) inversion.

[0101] 2. Significantly reduces boundary effects in regional modeling. Spherical cap harmonic analysis (SCHA) is used to parameterize the regional ionosphere, constructing an orthogonal basis function system that satisfies boundary conditions within a limited spatial range. This effectively avoids Runge oscillations and ill-conditioned problems in the normal equations caused by traditional spherical harmonic functions at the regional boundary, thereby improving the spatial resolution and modeling stability of the regional TEC.

[0102] 3. Achieve strict decoupling between ionospheric parameters and hardware delay parameters. By implementing QR orthogonal decomposition on the design matrix, a pure differential code bias observation subsystem is constructed, achieving mathematical decoupling between the ionospheric spatial coefficients and the receiver / satellite differential code bias (DCB) parameters. This strategy effectively reduces the strong correlation between parameters, making DCB estimation no longer affected by rapid spatiotemporal changes in the ionosphere, thereby significantly improving the robustness and reliability of GNSS hardware delay estimation.

[0103] 4. Improve numerical stability under extreme ionospheric conditions. Tikhonov regularization constraints are introduced in the back-substitution solution stage of ionospheric parameters, and the smoothing parameters are adaptively determined using the L-curve method. This enables the model to maintain stable solutions even under conditions such as low electron density at night, sparse observations, or uneven distribution of regional edge stations, thereby achieving continuous spatiotemporal reconstruction of regional TEC.

[0104] 5. Achieve high-precision continuous inversion of the regional ionosphere. By constructing a two-step solution framework of DCB global robust estimation + ionospheric parameter regularization reconstruction, while ensuring stable estimation of hardware delay, high temporal and spatial resolution modeling of the regional TEC field is achieved, providing more reliable ionospheric information products for ionospheric monitoring, space weather research, and GNSS precise positioning error correction.

[0105] The above description discloses only one or more preferred embodiments of the present invention, and should not be construed as limiting the scope of the present invention. Those skilled in the art will understand that all or part of the processes of the above embodiments can be implemented, and equivalent changes made in accordance with the claims of the present invention are still within the scope of the invention.

Claims

1. A robust method for GNSS regional ionospheric TEC modeling and differential code bias estimation, characterized in that, Includes the following steps: Step 1: Construct a GNSS ionospheric geometry-independent observation model; Step 2: Establish an adaptive mapping function that takes into account the variation of the ionospheric peak height; Step 3: Construct a regional ionospheric VTEC model based on spherical cap harmonic analysis; Step 4: Construct the network observation equations and perform linearization; Step 5: Use QR orthogonal transformation to achieve parameter decoupling; Step 6: Perform a global robust estimation of the differential code bias parameters; Step 7: Perform regularized reconstruction of ionospheric parameters; Step 8: Output the regional ionosphere TEC reconstruction results.

2. The GNSS regional ionospheric TEC modeling and differential code bias robust estimation method as described in claim 1, characterized in that, In step 1, a geometrically independent combination is constructed using GNSS dual-frequency pseudorange observations and carrier phase observations. The observation noise and multipath effect are reduced by carrier phase smoothing pseudorange technology to obtain high-precision geometrically independent ionospheric observations, thereby establishing an ionospheric observation equation that includes STEC and the differential code deviation between the receiver and the satellite.

3. The GNSS regional ionospheric TEC modeling and differential code bias robust estimation method as described in claim 2, characterized in that, The adaptive mapping function converts STEC to VTEC based on the dynamic change of the peak height of the F2 layer of the ionosphere. The shell height is provided with a priori constraints by the international reference ionospheric model and changes adaptively with the geomagnetic latitude of the ionospheric puncture point, the fixed solar longitude, and the observation epoch.

4. The GNSS regional ionospheric TEC modeling and differential code bias robust estimation method as described in claim 3, characterized in that, In step 3, a series expansion of VTEC is performed by constructing a generalized connected Legendre function basis that satisfies the boundary conditions within the spherical cap region.

5. The GNSS regional ionospheric TEC modeling and differential code bias robust estimation method as described in claim 4, characterized in that, In step 4, the adaptive mapping function and the regional ionospheric VTEC model based on spherical cap harmonic analysis are substituted into the geometry-independent observation model to uniformly model all receiver and visible satellite observation data in the regional GNSS network, forming a batch observation equation that includes ionospheric cap harmonic analysis coefficients and differential code deviation parameters, and is expressed in matrix form as a network-level linear observation model.

6. The GNSS regional ionospheric TEC modeling and differential code bias robust estimation method as described in claim 5, characterized in that, In step 5, the network-level linear observation model is decomposed into two independent subsystems by QR orthogonal transformation, one of which contains only the differential code bias parameter.

7. The GNSS regional ionospheric TEC modeling and differential code bias robust estimation method as described in claim 6, characterized in that, In step 6, the differential code deviation subsystems of multiple time batches are accumulated to construct a global normal equation. The parameter rank deficiency problem is solved by introducing a zero-sum constraint condition for satellite differential code deviation. The robust estimate of the receiver-satellite differential code deviation is obtained by using the constrained least squares method.

8. The GNSS regional ionospheric TEC modeling and differential code bias robust estimation method as described in claim 7, characterized in that, Step 7 introduces Tikhonov regularization constraints during the ionospheric parameter solution process. The smooth difference operator describes the continuous change of the ionospheric state between adjacent time batches, and the L-curve method is used to adaptively determine the regularization parameters.

9. The GNSS regional ionospheric TEC modeling and differential code bias robust estimation method as described in claim 8, characterized in that, In step 8, the VTEC distribution at any ionospheric puncture point within the study area is calculated based on the obtained SCHA coefficients.

Images