A GNSS baseline joint solution method of a multi-reference station network and a computer readable medium
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- WUHAN UNIV
- Filing Date
- 2023-06-14
- Publication Date
- 2026-07-21
Smart Images
Figure CN116882134B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to high-precision satellite navigation and positioning data processing technology, and more particularly to a method for joint calculation of GNSS baselines of a multi-reference station network and a computer-readable medium. Background Technology
[0002] GNSS technology, with its advantages of automation, real-time performance, high precision, and all-weather capability, has become a common method for monitoring surface deformation. For monitoring scenarios with relatively slow deformation, such as dams and bridges, a multi-static baseline network monitoring method is typically used. Static deformation monitoring mainly involves periodically and repeatedly observing deformation patterns to understand their characteristics, thereby providing a data source for geological disaster safety monitoring and early warning.
[0003] In GNSS static network monitoring, baseline resolution technology, as a crucial component of high-precision data post-processing, is key to ensuring millimeter-level accuracy in monitoring results. Baseline resolution technology utilizes the correlation between multiple stations observing multiple satellites and employs different linear combinations of carrier phase observations to calculate the relative vectors between stations. Regarding research findings on baseline resolution methods, Zhou Letao et al. analyzed different baseline resolution modes, using a star network composed of multiple baselines as the basic resolution unit for ambiguity resolution; Tang Weiming et al. studied ambiguity fixing and baseline resolution for a single epoch of the BeiDou Navigation Satellite System and conducted a comparative accuracy analysis; Wang Zhenfang et al. conducted a detailed analysis of the correlation of observations under GPS multi-baseline mode; Chu Lixin introduced the concept of equivalent equations in multi-baseline resolution of GPS / GLONASS combinations to construct multi-baseline resolution equations; Zhu Huizhong et al. studied a high-precision static positioning method for long-distance baselines in BDS, which, based on the resolution of wide-lane integer ambiguity, estimates the zenith tropospheric delay error and simultaneously performs carrier phase integer ambiguity resolution and positioning calculation.
[0004] The aforementioned studies primarily employed single-baseline or multi-baseline holistic solutions in the baseline calculation process; however, both approaches have limitations. The single-baseline approach, in solving independent baselines, fails to fully utilize the correlation between stations, resulting in a theoretically less rigorous approach. While the multi-baseline holistic solution approach considers the interrelationships between parameters and has a more rigorous mathematical model, it becomes more complex when dealing with a large number of stations, leading to a greater number of parameters and a more intricate mathematical model and solution process, thus limiting the applicability of this method to some extent. Summary of the Invention
[0005] To address the aforementioned technical problems, this invention proposes a joint solution method for GNSS baselines in a multi-reference station network and a computer-readable medium.
[0006] The technical solution of this invention is a joint solution method for GNSS baselines of a multi-reference station network, comprising the following steps:
[0007] Step 1: Obtain the initial 3D coordinates of each base station, transform the initial 3D coordinates of each base station into the coordinates of each base station in the Gaussian plane through the Gaussian projection algorithm, construct the base station plane coordinate point set through the coordinates of multiple base stations in the Gaussian plane, and obtain multiple baselines to be solved by the triangular mesh growth algorithm of the base station plane coordinate point set.
[0008] Step 2: Solve the double-difference integer ambiguity of each baseline to be solved;
[0009] Step 3: Construct the error equation for each baseline to be adjusted, the error equation for the multi-baseline triangulation network, the weight matrix of the double-difference observations of the multi-baseline triangulation network, and the adjustment mathematical model for multiple reference stations in sequence. Combine the adjustment mathematical model of multiple reference stations to calculate the corrected three-dimensional coordinates of each reference station.
[0010] Preferably, the set of plane coordinate points of the reference station mentioned in step 1 is defined as follows:
[0011] P = {P} i |i = 1, 2, ..., k}
[0012] Among them, P i Let k represent the coordinates of the i-th reference station in the Gaussian plane, and k represent the number of reference stations.
[0013] Step 1 involves using a triangular mesh growth algorithm to obtain multiple baselines to be solved from the set of plane coordinate points of the base station. The specific steps are as follows:
[0014] Step 1.1: Arbitrarily select a point in the point set P as the first vertex of the initial triangle;
[0015] Step 1.2: Search for the point in the point set P that is closest to the first vertex, and use it as the second vertex of the initial triangle. Connect the two points to generate the initial edge.
[0016] Step 1.3: Search in the point set P for the point that is closest to the midpoint of the initial edge and is not on a straight line with the first vertex and the second vertex of the initial triangle. Use this point as the third vertex of the initial triangle. Connect the three vertices to form a triangle, which is the initial triangle.
[0017] Step 1.4: Using the two edges determined in Step 1.3 as a base, repeat Step 1.3 until all edges are constructed. At this point, all edges can form a planar triangulation.
[0018] In the planar triangulation, the vector formed by the two reference stations corresponding to each edge in the three-dimensional space is the baseline to be solved.
[0019] Preferably, step 2 involves calculating the double-difference integer ambiguity of each baseline to be solved, as follows:
[0020] Step 2.1: Perform data preprocessing on the GNSS carrier phase observations of each baseline to be solved to obtain the preprocessed GNSS carrier phase observation data of the two reference stations;
[0021] Each baseline to be solved is selected. Assuming that the reference stations at both ends of the baseline are the i-th reference station and the j-th reference station, the GNSS carrier phase observation data of these two reference stations are subjected to gross error detection and cycle slip detection respectively. The GNSS carrier phase observation data containing gross errors are subjected to gross error removal processing, and the GNSS carrier phase observation data containing cycle slips are subjected to cycle slip repair processing, so as to obtain the GNSS carrier phase observation data of the two reference stations that make up the baseline after preprocessing.
[0022] Step 2.2: Determine the weight matrix of the double-difference observations using the law of error propagation;
[0023] Using the preprocessed GNSS carrier phase observation data of the i-th and j-th reference stations, a linearized carrier phase double-difference observation equation is constructed for each baseline to be solved; and based on the weight matrix of the original observation values determined by the satellite elevation angle weighting model, the weight matrix P of the double-difference observation values is determined by the error propagation law.
[0024] Step 2.3: Solve the double-difference integer ambiguity of each baseline to be solved;
[0025] The linearized carrier phase double-difference observation equation for the baseline formed by reference stations i and j is transformed into an error equation, as follows:
[0026] V = Aa + Bb - l
[0027] In the formula, V is the residual vector of double-difference observations, A is a diagonal matrix with diagonal elements of λ, a is the integer ambiguity vector of double-difference observations, B is a coefficient matrix including direction cosines and double-difference tropospheric projection functions, and b is the coordinate correction dX from the reference station j. j The relative tropospheric wet delay dT between two reference stations in the zenith direction w The parameter vector formed; l is the difference between the carrier phase double-difference observation value and its corresponding calculated value;
[0028] Based on the above error equation, ignoring the integer property of double-difference integer ambiguities, least squares adjustment is used to solve for the estimated double-difference ambiguity vector. Baseline and tropospheric parameter vectors Its cofactor matrix Q is expressed by the following formula:
[0029]
[0030]
[0031] In the formula, The cofactor matrix of the double-difference ambiguity. The cofactor matrix of baseline and tropospheric parameters, and They are transpose matrices and are both cofactor matrices between the ambiguity parameters and the baseline and tropospheric parameters;
[0032] Estimation using double-difference ambiguity vector and its cofactor matrix By combining the least squares ambiguity decorrelation adjustment algorithm, the optimal double-difference integer ambiguity vector is obtained through search. The ratio R of the suboptimal ambiguity to the quadratic form of the optimal ambiguity residual and the ambiguity resolution success rate Ps;
[0033] If both the ratio R and the success rate Ps are greater than the set threshold, then the optimal double-difference integer ambiguity vector is the final double-difference fixed value of the selected baseline, and the double-difference integer ambiguity solution of the baseline is completed.
[0034] Repeat steps 2.1-2.3 until the double-difference integer ambiguity of all baselines to be solved is completed.
[0035] Preferably, step 3 involves constructing the error equation for each baseline to be adjusted, as follows:
[0036] Each baseline to be solved is treated as a baseline to be adjusted;
[0037] Select the baseline to be adjusted from the plane triangulation network, which consists of the i-th and j-th reference stations. Substitute the corresponding fixed ambiguity values back into the carrier phase double-difference observation equation, and linearize the coordinates of the two stations constituting the baseline. After simplification, the following error equation can be obtained:
[0038]
[0039] i, j ∈ {1, 2, ..., k} and i ≠ j
[0040] In the formula, Let be the residual vector of double-difference observations between the i-th and j-th reference stations after ambiguity recovery, and let be the coefficient matrix. Let be the coefficient matrix of the j-th reference station, which includes direction cosines and tropospheric projection functions. Let be the coefficient matrix of the i-th reference station, which includes direction cosines and tropospheric projection functions. Let dT be the coordinate correction for the j-th reference station and the tropospheric wet delay in the zenith direction.w，j The parameter vector formed Let dT be the coordinate correction for the i-th reference station and the tropospheric wet delay in the zenith direction. w，i The constructed parameter vector, constant term The vector is the double-difference observation between the i-th and j-th base stations after ambiguity recovery, minus the calculated initial value, where k represents the number of base stations.
[0041] The error equation for constructing the multi-baseline triangulation network in step 3 is as follows:
[0042] By superimposing the error equations of all baselines to be adjusted in the triangulation network, the error equations of the multi-baseline triangulation network observations are constructed, as follows:
[0043] V o =B o XL o
[0044] In the formula, V o B is the residual vector of all double-difference observations of the baselines to be adjusted after ambiguity recovery. o The coefficient matrix of the parameter X to be estimated is formed by superimposing the error equations of all baselines to be adjusted. That is, a parameter vector consisting of the coordinate corrections of each reference station and the tropospheric wet delay in the zenith direction. Composition, i = 1, 2, ..., k, where k represents the number of base stations, constant term L o This is the vector obtained by subtracting the initial values from the double-difference observations of all baselines to be adjusted after ambiguity recovery.
[0045] Step 3 involves constructing the double-difference observation weight matrix of the multi-baseline triangulation network, as detailed below:
[0046] For GNSS carrier phase observation data of all baselines to be adjusted after ambiguity recovery, a satellite elevation angle weighting model is used to determine the variance matrix of the original observations. Then, the variance-covariance matrix of the double-difference observations is determined using the error propagation law. Taking the inverse of the matrix yields the observation weight matrix P. o .
[0047] Step 3, which involves constructing the mathematical model for the adjustment of the benchmark station with coordinate datum constraints, is as follows:
[0048] If coordinate datum constraints need to be introduced during triangulation adjustment, then the constraint equations are added:
[0049] V c =TX
[0050] In the formula, V cLet T be the residual of the coordinate datum constraint, and let T be the coefficient matrix of 4c×4k, where c is the number of constrained datums and c < k, and k is the number of reference stations. The matrix T is defined as follows: Assume that the m-th reference station, m = 1, 2, ..., k, is the reference station constrained by the coordinate datum, and its n-th reference station, n = 1, 2, ..., c, is sorted among all the reference stations constrained by the coordinate datum. Then the elements in the (4n-3)-th row and (4m-3)-th column, the (4n-2)-th row and (4m-2)-th column, and the (4n-1)-th row and (4m-1)-th column of matrix T are 1, and all other elements are 0.
[0051] Based on the observation weight matrix, the weight matrix P of the constraint equation is added in the form of an additional block diagonal matrix. c Its dimension is 4c×4c;
[0052] Matrix P c The definition is as follows: The (4n-3), (4n-2), and (4n-1)th diagonal elements of the matrix, n = 1, 2, ..., c, take the maximum value M, and all other elements take the value 0. Then, the adjustment mathematical model for multiple reference stations is:
[0053]
[0054] Step 3 involves calculating the corrected three-dimensional coordinates of each reference station using the adjustment mathematical model of multiple reference stations, as detailed below:
[0055] The adjustment mathematical models of multiple reference stations are used as the observation equations for Kalman filtering. A random walk process is used to describe the state changes of the station coordinates and the zenith tropospheric delay parameter, and the state equations of the parameters to be estimated are constructed. The coordinate corrections and the filtered solutions of the tropospheric parameters of each reference station are obtained through the Kalman filtering algorithm. The initial three-dimensional coordinates of each reference station are superimposed with the coordinate corrections of each reference station to obtain the corrected three-dimensional coordinates of each reference station.
[0056] The present invention also provides a computer-readable medium storing a computer program executed by an electronic device, wherein when the computer program is run on the electronic device, it executes the steps of the GNSS baseline joint solution method of the multi-reference station network.
[0057] The beneficial effects of this invention are:
[0058] This invention, based on the traditional single-baseline solution for obtaining station coordinates and ambiguity parameters, reconstructs the observation equations of each baseline in the triangulation network during network adjustment to uniformly solve for the station coordinate parameters. By considering the correlation between observed values and various parameters, the strength of the adjustment model is increased, improving the accuracy and reliability of the adjustment results.
[0059] This invention, based on the recovery of integer ambiguities for each satellite pair using a single baseline model, utilizes a multi-baseline overall solution model to jointly estimate coordinates and zenith tropospheric parameters. Through the ambiguity resolution in the first step, the number of parameters to be estimated during baseline network adjustment is significantly reduced, thus improving the speed of triangulation network adjustment. Attached Figure Description
[0060] Figure 1 : A schematic diagram of the point distribution in an embodiment of the present invention;
[0061] Figure 2 : Flowchart of the method according to an embodiment of the present invention. Detailed Implementation
[0062] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.
[0063] In specific implementation, the method proposed in the technical solution of this invention can be automatically executed by those skilled in the art using computer software technology. System devices for implementing the method, such as computer-readable storage media storing the corresponding computer program of the technical solution of this invention and computer equipment including the computer program running the corresponding computer program, should also be within the protection scope of this invention.
[0064] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to embodiments using four base stations. A schematic diagram of the location distribution of the four base stations is shown below. Figure 1 As shown, stations 1 and 3 are known and have fixed coordinates, and will be used as constraint references in subsequent adjustment processing. It should be understood that the specific embodiments described herein are for illustrative purposes only and are not intended to limit the scope of the invention.
[0065] The following is combined Figure 1-2 The technical solution of this invention is a joint solution method for GNSS baselines in a multi-reference station network, as detailed below:
[0066] like Figure 2 The diagram shown is a flowchart of a method according to an embodiment of the present invention.
[0067] Step 1: Obtain the initial 3D coordinates of each base station, transform the initial 3D coordinates of each base station into the coordinates of each base station in the Gaussian plane through the Gaussian projection algorithm, construct the base station plane coordinate point set through the coordinates of multiple base stations in the Gaussian plane, and obtain multiple baselines to be solved by the triangular mesh growth algorithm of the base station plane coordinate point set.
[0068] Step 2: Solve the double-difference integer ambiguity of each baseline to be solved;
[0069] Step 3: Construct the error equation for each baseline to be adjusted, the error equation for the multi-baseline triangulation network, the weight matrix of the double-difference observations of the multi-baseline triangulation network, and the adjustment mathematical model for multiple reference stations in sequence. Combine the adjustment mathematical model of multiple reference stations to calculate the corrected three-dimensional coordinates of each reference station.
[0070] The set of plane coordinate points of the reference station mentioned in step 1 is defined as follows:
[0071] P = {P} i |i = 1, 2, ..., k}
[0072] Among them, P i Let k represent the coordinates of the i-th reference station in the Gaussian plane, and k = 4 represent the number of reference stations.
[0073] Step 1 involves using a triangular mesh growth algorithm to obtain multiple baselines to be solved from the set of plane coordinate points of the base station. The specific steps are as follows:
[0074] Step 1.1: Arbitrarily select point 1 in the point set P as the first vertex of the initial triangle;
[0075] Step 1.2: Search for point 2 in the point set P that is closest to point 1, and use it as the second vertex of the initial triangle. Connect the two points to generate the initial edge L12.
[0076] Step 1.3: Search for point 3 in point set P that is closest to a point in the initial edge L12 and is not on the same straight line as points 1 and 2. Use this point 3 as the third vertex of the initial triangle. Connect the three vertices to form a triangle, which is the initial triangle.
[0077] Step 1.4: Based on the two edges L23 and L13 determined in Step 1.3, repeat Step 1.3 until all edges are constructed. At this point, all edges can form a planar triangulation.
[0078] Then, the vector formed by the two reference stations corresponding to each edge in the planar triangulation network in three-dimensional space is the baseline to be solved;
[0079] Step 2 involves calculating the double-difference integer ambiguity of each baseline to be solved, as detailed below:
[0080] Step 2.1: Perform data preprocessing on the GNSS carrier phase observations of each baseline to be solved to obtain the preprocessed GNSS carrier phase observation data of the two reference stations;
[0081] Each baseline to be solved is selected. Assuming that the reference stations at both ends of the baseline are reference station 1 and reference station 2, gross error detection and cycle slip detection are performed on the GNSS carrier phase observation data of the two reference stations respectively. The GNSS carrier phase observation data containing gross errors are processed to remove gross errors, and the GNSS carrier phase observation data containing cycle slips are processed to repair cycle slips, so as to obtain the GNSS carrier phase observation data of the two reference stations after preprocessing the baseline.
[0082] Step 2.2: Determine the weight matrix of the double-difference observations using the law of error propagation;
[0083] Using the preprocessed GNSS carrier phase observation data from reference station 1 and reference station 2, a linearized carrier phase double-difference observation equation is constructed for each baseline to be solved; and based on the weight matrix of the original observation values determined by the satellite elevation angle weighting model, the weight matrix P of the double-difference observation values is determined by the error propagation law.
[0084] Step 2.3: Solve the double-difference integer ambiguity of each baseline to be solved;
[0085] The linearized carrier phase double-difference observation equation for the baseline formed by reference station 1 and reference station 2 is transformed into an error equation, as follows:
[0086] V = Aa + Bb - l
[0087] In the formula, V is the residual vector of double-difference observations, A is a diagonal matrix with diagonal elements of λ, a is the integer ambiguity vector of double-difference observations, B is a coefficient matrix including direction cosines and double-difference tropospheric projection functions, and b is the coordinate correction dX2 of reference station 2 and the relative tropospheric wet delay dT between the two reference stations in the zenith direction. w The parameter vector formed; l is the difference between the carrier phase double-difference observation value and its corresponding calculated value;
[0088] Based on the above error equation, ignoring the integer property of double-difference integer ambiguities, least squares adjustment is used to solve for the estimated double-difference ambiguity vector. Baseline and tropospheric parameter vectors Its cofactor matrix Q is expressed by the following formula:
[0089]
[0090]
[0091] In the formula, The cofactor matrix of the double-difference ambiguity. The cofactor matrix of baseline and tropospheric parameters, and They are transpose matrices and are both cofactor matrices between the ambiguity parameters and the baseline and tropospheric parameters;
[0092] Estimation using double-difference ambiguity vector and its cofactor matrix By combining the least squares ambiguity decorrelation adjustment algorithm, the optimal double-difference integer ambiguity vector is obtained through search. The ratio R of the suboptimal ambiguity to the quadratic form of the optimal ambiguity residual and the ambiguity resolution success rate Ps. If both the ratio R and the success rate Ps are greater than the set threshold, then the optimal double-difference integer ambiguity vector is the final double-difference fixed value of the selected baseline, and at this time the double-difference integer ambiguity resolution of the baseline is completed;
[0093] Repeat steps 2.1-2.3 until the double-difference integer ambiguity of all baselines to be resolved is completed;
[0094] Step 3 involves constructing the error equation for each baseline to be adjusted, as detailed below:
[0095] Each baseline to be solved is treated as a baseline to be adjusted;
[0096] Select the baseline to be adjusted from the plane triangulation network, which consists of the i-th and j-th reference stations. Substitute the corresponding fixed ambiguity values back into the carrier phase double-difference observation equation, and linearize the coordinates of the two stations constituting the baseline. After simplification, the following error equation can be obtained:
[0097]
[0098] i, j ∈ {1, 2, ..., k} and i ≠ j
[0099] In the formula, Let be the residual vector of double-difference observations between the i-th and j-th reference stations after ambiguity recovery, and let be the coefficient matrix. Let be the coefficient matrix of the j-th reference station, which includes direction cosines and tropospheric projection functions. Let be the coefficient matrix of the i-th reference station, which includes direction cosines and tropospheric projection functions. Let dT be the coordinate correction for the j-th reference station and the tropospheric wet delay in the zenith direction. w，j The parameter vector formed Let dT be the coordinate correction for the i-th reference station and the tropospheric wet delay in the zenith direction. w，i The constructed parameter vector, constant term The vector is the double difference observation between the i-th and j-th base stations after ambiguity recovery, minus the calculated initial value. k = 4 represents the number of base stations. For this example, (i, j) can be (1, 2), (1, 3), (1, 4), (2, 3), (3, 4).
[0100] The error equation for constructing the multi-baseline triangulation network in step 3 is as follows:
[0101] By superimposing the error equations of all baselines to be adjusted in the triangulation network, the error equations of the multi-baseline triangulation network observations are constructed, as follows:
[0102] V o =B o XL o
[0103] In the formula, V o B is the residual vector of all double-difference observations of the baselines to be adjusted after ambiguity recovery. o The coefficient matrix of the parameter X to be estimated is formed by superimposing the error equations of all baselines to be adjusted. That is, a parameter vector consisting of the coordinate corrections of each reference station and the tropospheric wet delay in the zenith direction. Composition, i = 1, 2, ..., k, k = 4 represents the number of base stations, constant term L o This is the vector obtained by subtracting the initial values from the double-difference observations of all baselines to be adjusted after ambiguity recovery; in this example,
[0104]
[0105]
[0106]
[0107]
[0108] Step 3 involves constructing the double-difference observation weight matrix of the multi-baseline triangulation network, as detailed below:
[0109] For GNSS carrier phase observation data of all baselines to be adjusted after ambiguity recovery, a satellite elevation angle weighting model is used to determine the variance matrix of the original observations. Then, the variance-covariance matrix of the double-difference observations is determined using the error propagation law. Taking the inverse of the matrix yields the observation weight matrix P. o ;
[0110] Step 3, which involves constructing adjustment mathematical models for multiple reference stations, is detailed below:
[0111] If coordinate datum constraints need to be introduced during triangulation adjustment, then the constraint equations are added:
[0112] Vc =TX
[0113] In the formula, V c Let T be the residual of the coordinate datum constraint, a 4c×4k coefficient matrix, where c=2 is the number of constrained datums and k=4 is the number of reference stations. The matrix T is defined as follows: Assume the m-th reference station (m=1,2,...,k) is a coordinate-constrained reference station, and its order among all coordinate-constrained reference stations is n-th (n=1,2,...,c). Then, the elements in matrix T at the (4n-3)-th row and (4m-3)-th column, the (4n-2)-th row and (4m-2)-th column, and the (4n-1)-th row and (4m-1)-th column are 1, and all other elements are 0. For this example,
[0114]
[0115] In the formula, 0 4×4 This represents a 4×4 zero-value matrix. This represents a 4×4 diagonal matrix where the first three diagonal elements are 1.
[0116]
[0117] Based on the observation weight matrix, the weight matrix P of the constraint equation is added in the form of an additional block diagonal matrix. c Its dimension is 4c×4c, where c=2 is the number of constraint datums;
[0118] Matrix P c The definition is as follows: the (4n-3), (4n-2), and (4n-1)th diagonal elements of this matrix, n = 1, 2, ..., c, take the maximum value M = 100000, and all other elements take the value 0. For this example:
[0119]
[0120] The adjustment mathematical model for multiple reference stations is as follows:
[0121]
[0122] Step 3 involves calculating the corrected three-dimensional coordinates of each reference station using the adjustment mathematical model of multiple reference stations, as detailed below:
[0123] The adjustment mathematical models of multiple reference stations are used as the observation equations for Kalman filtering. A random walk process is used to describe the state changes of the station coordinates and the zenith tropospheric delay parameter, and the state equations of the parameters to be estimated are constructed. The coordinate corrections and the filtered solutions of the tropospheric parameters of each reference station are obtained through the Kalman filtering algorithm. The initial three-dimensional coordinates of each reference station are superimposed with the coordinate corrections of each reference station to obtain the corrected three-dimensional coordinates of each reference station.
[0124] A specific embodiment of the present invention also provides a computer-readable medium.
[0125] The computer-readable medium is a server workstation;
[0126] The server workstation stores the computer program executed by the electronic device. When the computer program runs on the electronic device, it causes the electronic device to execute the steps of the GNSS baseline joint solution method for multi-reference station network according to the present invention.
[0127] It should be understood that any parts not described in detail in this specification belong to the prior art.
[0128] It should be understood that the above description of the preferred embodiments is quite detailed, but it should not be considered as a limitation on the scope of protection of this invention. Those skilled in the art, under the guidance of this invention, can make substitutions or modifications without departing from the scope of protection of the claims of this invention, and all such substitutions or modifications fall within the scope of protection of this invention. The scope of protection of this invention should be determined by the appended claims.
Claims
1. A method for joint calculation of GNSS baselines in a multi-reference station network, characterized in that, Includes the following steps: Step 1: Obtain the initial 3D coordinates of each base station, transform the initial 3D coordinates of each base station into the coordinates of each base station in the Gaussian plane through the Gaussian projection algorithm, construct the base station plane coordinate point set through the coordinates of multiple base stations in the Gaussian plane, and obtain multiple baselines to be solved by the triangular mesh growth algorithm of the base station plane coordinate point set. Step 2: Solve the double-difference integer ambiguity of each baseline to be solved; Step 3: Construct the error equation for each baseline to be adjusted, the error equation for the multi-baseline triangulation network, the weight matrix of the double-difference observations of the multi-baseline triangulation network, and the adjustment mathematical model for multiple reference stations in sequence. Combine the adjustment mathematical model of multiple reference stations to calculate the corrected three-dimensional coordinates of each reference station. The error equation for constructing the multi-baseline triangulation network in step 3 is as follows: By superimposing the error equations of all baselines to be adjusted in the triangulation network, the error equations of the multi-baseline triangulation network observations are constructed, as follows: In the formula, This represents the residual vector of all double-difference observations of the baselines to be adjusted after ambiguity recovery. The parameters to be estimated are the sum of the error equations of all the baselines to be adjusted. The coefficient matrix, the parameters to be estimated That is, a parameter vector consisting of the coordinate corrections of each reference station and the tropospheric wet delay in the zenith direction. composition, Indicates the number of reference stations, constant term The vector is the result of subtracting the initial values from the double-difference observations of all baselines to be adjusted after ambiguity recovery. Step 3 involves constructing the double-difference observation weight matrix of the multi-baseline triangulation network, as detailed below: For GNSS carrier phase observation data of all baselines to be adjusted after ambiguity recovery, a satellite elevation angle weighting model is used to determine the variance matrix of the original observations. Then, the error propagation law is used to determine the variance-covariance matrix of the double-difference observations. Taking the inverse of the matrix gives the observation weight matrix. ; Step 3, which involves constructing adjustment mathematical models for multiple reference stations, is detailed below: If coordinate datum constraints need to be introduced during triangulation adjustment, then constraint equations are added to the error equations for constructing multi-baseline triangulation networks: In the formula, The residuals are the coordinate datum constraints. for The coefficient matrix, c The number of constraint datums and c <k,k Number of base stations ， matrix The definition is as follows: Assume the first m, m = 1, 2, …, k The first base station is a coordinate-constrained base station, and its ranking among all coordinate-constrained base stations is [number]. n, n = 1, 2, …, c Then the matrix The Middle (4n-3) Line number (4m-3) Column, No. (4n-2) Line number (4m-2) List and No. (4n-1) Line number (4m-1) The first element in the column is 1, and all other elements are 0. Based on the observation weight matrix, the weight matrix of the constraint equations is increased by adding a block diagonal matrix. Its dimension is ; matrix The definition is as follows: the first element in this matrix... (4n-3) , No. (4n-2) Passing the exam (4n-1) One diagonal element, n = 1, 2, …, c Take the maximum value M All other elements take values. 0 The adjustment mathematical model for multiple reference stations is as follows: , 。 2. The GNSS baseline joint solution method for a multi-reference station network according to claim 1, characterized in that: The set of plane coordinate points of the reference station mentioned in step 1 is defined as follows: PointSet ={ Point i |i =1,2,…, k} in, Point i Indicates the first i The coordinates of a reference station in the Gaussian plane k Indicates the number of base stations; Step 1 involves using a triangular mesh growth algorithm to obtain multiple baselines to be solved from the set of plane coordinate points of the base station. The specific steps are as follows: Step 1.1: In the point set PointSet Choose any point in the triangle as the first vertex of the initial triangle; Step 1.2: Search for the point set PointSet The point closest to the first vertex is chosen as the second vertex of the initial triangle, and the two points are connected to form the initial edge. Step 1.3: In the point set PointSet The point closest to the midpoint of the initial edge and not on the same straight line as the first and second vertices of the initial triangle is selected as the third vertex of the initial triangle. The three vertices are connected to form a triangle, which is the initial triangle. Step 1.4: Using the two edges determined in Step 1.3 as a base, repeat Step 1.3 until all edges are constructed. At this point, all edges can form a planar triangulation. In the planar triangulation, the vector formed by the two reference stations corresponding to each edge in the three-dimensional space is the baseline to be solved.
3. The GNSS baseline joint solution method for a multi-reference station network according to claim 2, characterized in that: The calculation of the double-difference integer ambiguity of each baseline to be solved is as follows: Step 2.1: Preprocess the GNSS carrier phase observation data of each baseline to be solved to obtain the preprocessed GNSS carrier phase observation data of the two reference stations; For each baseline to be solved, assume that the reference stations at both ends of the baseline are the 1st, 2nd, and 3rd base stations respectively. i The first base station and the first j For each of the two reference stations, gross error detection and cycle slip detection are performed on the GNSS carrier phase observation data. Gross error removal is performed on the GNSS carrier phase observation data containing gross errors, and cycle slip repair is performed on the GNSS carrier phase observation data containing cycle slips. This results in the GNSS carrier phase observation data of the two reference stations after baseline preprocessing. Step 2.2: Determine the weight matrix of the double-difference observations using the law of error propagation; After preprocessing, the first i The first base station and the first j The GNSS carrier phase observation data of each reference station are used to construct the linearized carrier phase double-difference observation equation for each baseline to be solved. Based on the weight matrix of the original observations determined by the satellite elevation angle weighting model, the weight matrix of the double-difference observations is then determined using the error propagation law. P ; Step 2.3: Solve the double-difference integer ambiguity of each baseline to be solved; Base station i and j The linearized carrier phase double-difference observation equations for the baseline are transformed into error equations, as follows: In the formula, The residual vector of double-difference observations. diagonal elements are diagonal matrix, This is a double-difference integer ambiguity vector; The coefficient matrix includes direction cosines and double-difference tropospheric projection functions. For the base station j coordinate corrections Relative tropospheric wet delay between two reference stations in the zenith direction The parameter vector formed; This is the difference between the carrier phase double-difference observation value and its corresponding calculated value; Based on the above error equation, ignoring the integer property of double-difference integer ambiguities, least squares adjustment is used to solve for the estimated double-difference ambiguity vector. Baseline and tropospheric parameter vectors and cofactor matrix The formula is expressed as follows: In the formula, The cofactor matrix of the double-difference ambiguity. The cofactor matrix of baseline and tropospheric parameters, and They are transpose matrices and are both cofactor matrices between the ambiguity parameters and the baseline and tropospheric parameters; Estimation using double-difference ambiguity vector and cofactor matrix By combining the least squares ambiguity decorrelation adjustment algorithm, the optimal double-difference integer ambiguity vector is obtained through search. The ratio of the suboptimal fuzziness to the quadratic residual of the optimal fuzziness R And ambiguity resolution success rate Ps ; If the ratio R and success rate Ps If all values are greater than the set threshold, then the optimal double-difference integer ambiguity vector is the final double-difference fixed value of the selected baseline, and the double-difference integer ambiguity solution of the baseline is completed. Repeat steps 2.1-2.3 until the double-difference integer ambiguity of all baselines to be solved is completed.
4. The GNSS baseline joint solution method for a multi-reference station network according to claim 3, characterized in that: Step 3 involves constructing the error equation for each baseline to be adjusted, as detailed below: Each baseline to be solved is treated as a baseline to be adjusted; Select the first triangulation in the planar triangulation. i The first base station and the first j The baseline to be adjusted is composed of several reference stations. The corresponding fixed ambiguity values are resubstituted into the carrier phase double-difference observation equation, and the coordinates of the two stations constituting the baseline are linearized. After simplification, the following error equation can be obtained: In the formula, For the fuzziness recovery of the first i The first base station and the first j Residual vector of double-difference observations between reference stations, coefficient matrix For the first j Each reference station contains a coefficient matrix including direction cosines and tropospheric projection functions. For the first i Each reference station contains a coefficient matrix including direction cosines and tropospheric projection functions. For the first j Coordinate corrections for each reference station and tropospheric wet delay in the zenith direction The parameter vector formed For the first i Coordinate corrections for each reference station and tropospheric wet delay in the zenith direction The constructed parameter vector, constant term For the fuzziness recovery of the first i The first base station and the first j The vector obtained by subtracting the initial values from the double-difference observations between the reference stations k Indicates the number of base stations.
5. The GNSS baseline joint solution method for a multi-reference station network according to claim 4, characterized in that: Step 3 involves calculating the corrected three-dimensional coordinates of each reference station using the adjustment mathematical model of multiple reference stations, as detailed below: The adjustment mathematical models of multiple reference stations are used as the observation equations for Kalman filtering. A random walk process is used to describe the state changes of the station coordinates and the zenith tropospheric delay parameter, and the state equations of the parameters to be estimated are constructed. The coordinate corrections and the filtered solutions of the tropospheric parameters of each reference station are obtained through the Kalman filtering algorithm. The initial three-dimensional coordinates of each reference station are superimposed with the coordinate corrections of each reference station to obtain the corrected three-dimensional coordinates of each reference station.
6. A computer-readable medium, characterized in that, It stores a computer program executed by an electronic device, which, when run on the electronic device, causes the electronic device to perform the steps of the method as described in any one of claims 1-5.