CN109541647A - GNSS multipath effect modification method based on hemisphere grid point model - Google Patents
GNSS multipath effect modification method based on hemisphere grid point model Download PDFInfo
- Publication number
- CN109541647A CN109541647A CN201811527763.6A CN201811527763A CN109541647A CN 109541647 A CN109541647 A CN 109541647A CN 201811527763 A CN201811527763 A CN 201811527763A CN 109541647 A CN109541647 A CN 109541647A
- Authority
- CN
- China
- Prior art keywords
- station
- satellite
- model
- grid
- grid point
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 230000000694 effects Effects 0.000 title claims abstract description 34
- 238000002715 modification method Methods 0.000 title abstract 2
- 238000000034 method Methods 0.000 claims description 41
- 238000004364 calculation method Methods 0.000 claims description 19
- 238000012545 processing Methods 0.000 claims description 14
- 238000012937 correction Methods 0.000 claims description 12
- 238000013507 mapping Methods 0.000 claims description 7
- 230000008859 change Effects 0.000 claims description 4
- 238000012360 testing method Methods 0.000 claims description 2
- 239000011159 matrix material Substances 0.000 description 18
- 238000010586 diagram Methods 0.000 description 9
- 230000003313 weakening effect Effects 0.000 description 6
- 230000008569 process Effects 0.000 description 5
- 238000005516 engineering process Methods 0.000 description 2
- 238000001914 filtration Methods 0.000 description 2
- 230000003068 static effect Effects 0.000 description 2
- 238000012935 Averaging Methods 0.000 description 1
- 238000007792 addition Methods 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 230000008030 elimination Effects 0.000 description 1
- 238000003379 elimination reaction Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000000737 periodic effect Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
- G01S19/01—Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
- G01S19/13—Receivers
- G01S19/22—Multipath-related issues
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
- G01S19/01—Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
- G01S19/13—Receivers
- G01S19/23—Testing, monitoring, correcting or calibrating of receiver elements
Landscapes
- Engineering & Computer Science (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Computer Networks & Wireless Communication (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Position Fixing By Use Of Radio Waves (AREA)
Abstract
The present invention provides a kind of GNSS multipath effect modification method based on hemisphere grid point model, double difference observation residual information including the fuzziness fixed period between survey station after resolving to multistation data sheet day is handled, and the parametrization of hemisphere grid point model indicates at each survey station；The residual error of record is mapped in the hemisphere grid on related survey station by successively utilizing fuzziness fixed period corresponding double difference observation residual information as observation, obtains first kind normal equation by the double difference observation residual error that will acquire；It is constrained based on model grid point value, constructs the second class normal equation；It is constrained based on changing value between model grid points, building third class normal equation is overlapped three classes normal equation to obtain final normal equation, solves each station hemisphere grid points model parameter, corrects multipath effect bring error in subsequent observation.The present invention can effectively weaken multipath effect error, not limited by conditions such as distance between sites.
Description
Technical Field
The invention belongs to the field of global satellite navigation systems, and particularly relates to a technology for weakening a multipath effect existing in GNSS precision data processing and improving the precision of a resolving result.
Background
At present, aiming at the influence of the multipath effect on the positioning precision, a modeling or signal processing technology is mainly used for analysis in the data processing process so as to achieve the effect of separating out the multipath error and weaken the influence of the multipath error on the GNSS precision calculation. Ding[1]And Seeber[2]The research of the people finds that the multipath error has repeatability in time and space, when static observation is carried out, the antenna of the receiver is unchanged relative to the surrounding environment, the relative motion of the satellite can cause the signal to generate multipath effect, and the motion trail of the satellite has certain periodic repetition relative to the antenna of the receiver and the surrounding environment. Ragheb[3]And Larson[4]The sidereal day filtering method proposed by the people is to correct multi-path errors based on a time domain, and perform filtering processing by using the sunday repeatability characteristic of a satellite constellation aiming at multi-path space-time repeatability under a static observation environment. But different satellite systems have different satellite types, and their altitude and orbit are different, such as for the BDS system, the orbit repetition period of GEO satellites and IGSO satellites is about one day, but the repetition period of MEO satellites is about seven days, which results in inconsistency of the satellite orbit repetition period,when the method is used for weakening the multipath effect, the problem of different repetition periods of different constellation satellites exists firstly; secondly, for observed quantities of two repetition periods, even if a small advance quantity existing in the satellite orbit repetition period is corrected, the satellite elevation angle and the satellite azimuth angle of the same satellite at the corresponding moment on the same observation station are not completely consistent, so that a certain deviation exists in the daytime repetition of the multipath effect, and the effectiveness of the method is caused to have a certain problem.
In recent years Dong[5]Et al also propose a method for multipath error cancellation based on the spatial domain. The method is based on an application mode of one antenna with multiple antennas, multiple antennas are required to share the same receiver clock, the multipath error of a single satellite is obtained through calculation, and the space is divided into a reticular hemisphere according to the altitude angle and the azimuth angle of the satellite. For a certain grid point, when a co-view satellite passes through the periphery of the grid point among the stations, averaging the residual errors of the observed values of the satellites after the single-difference ambiguity among the stations is fixed, weakening the noise influence, keeping the low-frequency multipath effect, and then utilizing the model to weaken the multipath effect error. The method is mainly limited by special hardware requirements, a receiver is required to support multiple antennas and can share a clock reference, the length of a base line is limited by hardware conditions such as signal attenuation during cable transmission, and the method cannot be applied to the base line with a longer distance (such as more than 500 meters).
Reference to the literature
[1]Ding,X.,Chen,Y.,Zhu,J.,&Huang,D.(1999).Surface deformationdetection using gps.Proceedings of International Technical Meeting of theSatellite Division of the Institute of Navigation,53-62.
[2]Seeber,G.,Menge,F.,C.,Wübbena,G.,&Schmitz,M.(1998).PreciseGPS Positioning Improvements by Reducing Antenna and Site DependentEffects.Advances in Positioning and Reference Frames.Springer BerlinHeidelberg.
[3]Ragheb,A.E.,Clarke,P.J.,&Edwards,S.J.(2007).Gps siderealfiltering:coordinate-and carrier-phase-level strategies.Journal of Geodesy,81(5),325-335.
[4]Larson,K.M.,Bilich,A.,&Axelrad,P.(2007).Improving the precision ofhigh-rate gps.Journal of Geophysical Research Solid Earth,112(B5).
[5]Dong,D.,Wang,M.,Chen,W.,Zeng,Z.,Song,L.,&Zhang,Q.,et al.(2016).Mitigation of multipath effect in gnss short baseline positioning by themultipath hemispherical map.Journal of Geodesy,90(3),255-262.
Disclosure of Invention
Aiming at the problem that the precision of GNSS precision data processing in a global satellite navigation system is influenced by a multipath effect, the invention provides a method for establishing a GNSS multipath model at each observation station by using a semi-celestial sphere grid point to weaken the multipath influence.
In order to solve the above technical problems, the present invention provides a GNSS multipath effect correction method based on a semi-celestial grid point model, comprising the following steps,
step 1, processing double-difference observation value residual information of a fixed period of ambiguity between stations after single-day solution of multi-station data, realizing the following,
step 1.1, setting an index list according to data containing ambiguity fixing information obtained by resolving, wherein each index relates to two survey stations and two satellites and is respectively set as a survey station m and n and a satellite j and k;
step 1.2, processing the resolved residual error information aiming at each index, acquiring a double-difference residual error formed by two co-view satellites on two stations related to each index, forming a data record, and simultaneously recording an epoch serial number and the altitude angle and the azimuth angle of the two co-view satellites on the two stations;
step 2, parameterization representation of the semi-celestial grid point model at each survey station is realized as follows,
step 2.1, grid point range value setting is carried out on the semi-celestial sphere model on each survey station, semi-celestial sphere grid points are divided according to the altitude angle as latitude, and the azimuth angle as longitude;
step 2.2, setting a half celestial sphere model grid point interval d on each test station;
step 2.3, taking the grid point value of the model on each measuring station as an unknown quantity to be solved;
step 3, using the obtained double-difference observation value residual error as an observation value, and constructing a method equation corresponding to a half-celestial sphere lattice point model parameter, wherein the method equation comprises the steps of mapping the recorded residual error to the half-celestial sphere lattice network on the related station according to the altitude angles and the azimuth angles of m and n satellites j and k of the station in record by sequentially using double-difference observation value residual error information corresponding to a ambiguity fixed time period to obtain a first type of method equation; constructing a second type of normal equation based on the constraint of the model grid point value; constructing a third-class equation based on the constraint of the change value between the points of the model grid;
and 4, superposing the three types of equation obtained in the step 3 to obtain a final equation, solving the half-celestial sphere grid point model parameters of each station, and correcting errors caused by multipath effects in subsequent observed values.
And the step of mapping the recorded residual errors into the semi-celestial grid network on the measuring station according to the recorded altitude angles and azimuth angles of the satellites j and k at the measuring stations m and n comprises the step of respectively solving model parameter calculation values between the station m and the satellite j, between the station m and the satellite k, between the station n and the satellite j and between the station n and the satellite k
Furthermore, the calculation value of the model parameter between the station a and the satellite b is obtaineda, m or n, b, j or k, as follows,
acquiring altitude angle and azimuth angle data of a satellite b on a station a, dividing according to grid points, indexing the position of a puncture point of the satellite b in a semispherical grid at the moment of recording the data, and setting longitude and latitude (x, y), wherein a model parameter calculation value in a semispherical grid point model isRepresenting a multipath error experienced by a non-differential observation between station a and satellite b;
when the multipath error between the station a and the satellite b relates to the values of four points on the grid model, extracting an interpolation coefficient by using a bilinear interpolation method; when the multipath error between the station a and the satellite b relates to the values of three points on the grid model, extracting an interpolation coefficient by using a triangular interpolation method;
expressing model parameter calculation values based on interpolation coefficients
And, in step 4, it is determined that the satellite b on the station a needs to be corrected for the multipath error, a grid point model of the station a is found, a grid where the satellite b is located is determined, since four grid point values or three grid point values in the grid where the satellite b is located are known, an error value of the multipath effect at the position where the satellite b is located is obtained according to a corresponding interpolation coefficient, and then the error value is obtainedThe observed value of the satellite b at the observation station a may be corrected as a correction number.
The invention can effectively weaken multipath effect error, and compared with the prior art, the invention has the advantages and technical effects that:
1. the invention can establish a unified model of a plurality of GNSS systems, and different systems and different satellite constellations of the same system can be corrected by using the same model, so that the establishment and the use of the model are simpler.
2. The model established by the invention can be used for superposition modeling by using multi-time-segment data, the model result is not influenced by the tiny lead of the satellite orbit repetition period, and the use effectiveness is further improved.
3. The invention can realize multi-station integral modeling during model building, integrally solves and obtains a model result, is not limited by an application mode, does not need the hardware condition of one machine with multiple antennas, and is not limited by the distance between stations.
Drawings
Fig. 1 is a schematic diagram of dividing semispherical grid points according to an embodiment of the invention.
Fig. 2 is a diagram illustrating interpolation of grid points according to an embodiment of the present invention, wherein fig. 2(a) is a diagram illustrating a case where a multipath error relates to values of four points on a grid model, and fig. 2(b) is a diagram illustrating a case where a multipath error relates to values of three points on the grid model.
Fig. 3 is a schematic diagram illustrating constraint of variation values between parameters of a mesh point according to an embodiment of the present invention, where fig. 3(a) is a schematic diagram illustrating constraints between a normal mesh point and four surrounding mesh points, fig. 3(b) is a schematic diagram illustrating constraints between a normal mesh point and three surrounding mesh points, and fig. 3(c) is a schematic diagram illustrating constraints between a normal mesh point and a plurality of surrounding mesh points.
FIG. 4 is a flow chart of an embodiment of the present invention.
Detailed Description
The technical solution is described in detail below with reference to the accompanying drawings and examples.
The method for weakening the multipath effect by the celestial grid point model effectively utilizes multipath information contained in double-difference observation value residual errors between stations (step 1), can flexibly set the starting and stopping range of grid points and the density of the grid points in the data processing process, and simultaneously utilizes a normal equation to superpose and solve the grid point values of the multipath model according to the practical station environment with reasonable constraint (step 2, step 3 and step 4.1). When the model is subsequently used for multi-path error correction, the realization method is simple and quick, and data of multiple periods and multiple time periods can be selected to participate in modeling together, so that the accuracy and the reliability of model construction are improved (step 4.2 and step 4.3).
Referring to fig. 4, the method for weakening the multipath effect by using the celestial sphere grid point model provided by the embodiment specifically includes the following steps:
step 1, performing single-day calculation on multi-station data, acquiring double-difference residual error information in a period of fixed ambiguity between stations, and processing the information, wherein the method specifically comprises the following steps:
step 1.1, an index list is set according to the data containing the ambiguity fixing information obtained through calculation, and each index record relates to two stations and two satellites.
In an embodiment, each index record includes the numbers of two stations, the PRN numbers (pseudorandom noise codes) of two satellites, and the number of start-stop epochs with fixed ambiguity.
Step 1.2, after the index records are obtained, processing the resolved residual error information aiming at each index, wherein the residual error information comprises satellite related information (specifically comprising epoch count, residual error, satellite altitude angle, azimuth angle and the like) participating in resolving on each observation station, obtaining double-difference residual errors formed by two co-view satellites on two observation stations related to each index, forming data records, and simultaneously recording the epoch sequence number and the altitude angle and the azimuth angle of the two co-view satellites on the two observation stations.
As shown in fig. 1, data of N periods are recorded, which are first period double-difference observation value residual data, second period double-difference observation value residual data, and … nth period double-difference observation value residual data.
Step 2, parameterization representation of the semi-celestial lattice point model at each survey station specifically comprises the following steps:
and 2.1, setting grid point range values of the semi-celestial sphere model on each survey station, wherein semi-celestial sphere grid points can be divided by taking the altitude angle as latitude and the azimuth angle as longitude. Setting the initial value of latitude value of grid point, i.e. the minimum value of altitude angle as B0. Because of poor observation quality of low altitude angle data, the GNSS precision data processing usually eliminates low altitude angles (smaller than B)0) The satellite data of (1) is generally set to be B less than or equal to 0 DEG0Is less than or equal to 10 degrees. For observation data with a large altitude angle, it is generally considered that the observation data is not affected by multipath effect or is less affected, so the maximum value of the grid point latitude circle, i.e. the maximum value of the altitude angle, is set as B1In general, B is 45 DEG or less1< 90 deg. The range value of the grid point longitude is set to be 0-360 degrees according to the range value of the azimuth. In addition, a grid point parameter is set at the zenith. The zenith is the maximum value of the altitude angle of 90 degrees, and the latitude circle at the moment is a point similar to the pole of the earth, so that the point only needs to be provided with a grid point parameter.
And 2.2, setting the space of the half celestial sphere model grid points on each survey station as d, wherein d is generally set as an integer for convenient division, such as 1 degree, 2 degrees, 5 degrees, 10 degrees and the like.
Step 2.3, after setting the semi-celestial sphere grid point range value and the grid point interval value, taking the grid point value of the model on each survey station as the unknown quantity to be solved, referring to fig. 1, and taking B as0＝0°，B1The form of dividing grid points on an observation station is given by an example of 80 ° and d 10 ° (black points are to-be-solved grid points):
step 3, using the obtained double-difference observation value residual error as an observation value, and constructing a normal equation corresponding to the half-celestial sphere lattice point model parameter, wherein the method comprises the steps of using the obtained double-difference observation value residual error as an observation value, constructing a normal equation corresponding to the half-celestial sphere lattice point model parameter, and mapping the recorded residual error into the half-celestial sphere lattice network on the related observation station according to the altitude angles and the azimuth angles of m and n satellites j and k in the observation station in record by sequentially using double-difference observation value residual error information corresponding to the ambiguity fixed time period to obtain a first kind of normal equation; constructing a second type of normal equation based on the constraint of the model grid point value; and constructing a third-class equation based on the constraint of the change value between the points of the model grid.
And mapping the recorded residual errors into the semi-celestial grid network on the related survey station according to the recorded altitude angles and azimuth angles of the satellites j and k at the survey stations m and n, wherein the model parameter calculation values between the station m and the satellite j, between the station m and the satellite k, between the station n and the satellite j and between the station n and the satellite k are respectively solved
The invention further provides that the calculation value of the model parameter between the station a and the satellite b is obtaineda, m or n, b, j or k, as follows,
acquiring altitude angle and azimuth angle data of a satellite b on a station a, dividing according to grid points, indexing the position of a puncture point of the satellite b in a semispherical grid at the moment of recording the data, and setting longitude and latitude (x, y), wherein a model parameter calculation value in a semispherical grid point model isRepresenting a multipath error experienced by a non-differential observation between station a and satellite b;
when the multipath error between the station a and the satellite b relates to the values of four points on the grid model, extracting an interpolation coefficient by using a bilinear interpolation method; when the multipath error between the station a and the satellite b relates to the values of three points on the grid model, extracting an interpolation coefficient by using a triangular interpolation method;
model parameter calculation based on interpolation coefficientsValue of
In an embodiment, step 3 specifically includes:
and 3.1, sequentially utilizing the double-difference observation value residual error information corresponding to the ambiguity fixed time period, and mapping the residual error of the data record into a semi-celestial grid network on the related corresponding observation station according to the altitude angles and azimuth angles of the satellites j and k at the m and n of the observation station in the record. The realization method comprises the following steps:
step 3.1.1, firstly, acquiring altitude angle and azimuth angle data of the satellite j on the station m, dividing according to grid points shown in fig. 1, indexing the data to record the position of the puncture point of the satellite j in the semispherical grid at the moment, and setting the longitude and latitude as (x, y), wherein the calculated value in the semispherical grid point model isNamely the multipath error borne by the non-differential observation value between the station m and the satellite j;
step 3.1.2 multipath error between station m and satellite j relates to the values of four points on the grid modelAnd(the altitude of satellite j is at B0And B1In between) or three pointsAnd(the altitude of satellite j is greater than B1Time), two cases are shown in fig. 2(a) and 2(b) (squares represent grid point parameters, circles represent signal propagation from satellite jThe positions of puncture points in the semi-celestial sphere model on the station m also correspond to model parameter calculation values of the semi-celestial sphere grid point model)；
Step 3.1.3, if the position relation between the puncture point in the satellite j signal propagation direction and the semi-celestial grid point parameter is as shown in fig. 2(a), then a bilinear interpolation method is used, and 4 grid point parameters are assumedThe corresponding longitude and latitude coordinates are respectively (x)1,y1)、(x2,y1)、(x1,y2)、(x2,y2) Then, thenThe model parameter values of (a) may be expressed in the form:
wherein the interpolation coefficientsThe following is obtained by taking the following equation,
the interpolation mode is bilinear interpolation, and the four grid point parameters respectively correspond to an interpolation coefficient.
If the position relationship between the puncture point in the signal propagation direction of the satellite j and the parameters of the semi-celestial grid point is as shown in fig. 2(b), a plane interpolation method is adopted, and the three grid point parameters are assumed to beThe corresponding longitude and latitude coordinates are respectively (x)1,y1)、(x2,y2)、(x3,y3) Then, thenThe model parameter values of (a) may be expressed in the form:
wherein the interpolation coefficients u, v are obtained as follows,
the interpolation is a triangular interpolation, in which the coefficients of one lattice point can be expressed by the coefficients of the other two lattice points.
Step 3.1.4, similar to steps 3.1.2 and 3.1.3, can calculate the calculated values of the model parameters between station m and satellite k, between station n and satellite j, and between station n and satellite k
Solving the model parameter calculation values between the station m and the satellite k, between the station n and the satellite j, and between the station n and the satellite kThen, step 3.1.3 may be repeated, and the strategy in step 3.1.3 may be followed when 4 or 3 mesh points are involved.
At this time, if the two-difference observed value residual value in the data record is s:
step 3.2, constructing a first type of normal equation: and (3) constructing a method equation according to the formulas (1) to (3), wherein a matrix formed by grid point parameters related to the observation equation (3) is a grid point parameter vector X to be solved, s is an observed value, b is a coefficient matrix, and a grid point parameter vector containing the upper celestial sphere model of each station to be solved is X.
The calculated values of the model parameters obtained by the formula (1) are given belowAnd (3) the corresponding observation equation expansion equation:
wherein:
wherein,in order to interpolate the coefficients of the image,the grid point parameters.
When the calculated value of the model parameter is obtained by the formula (2)When the observation equation expansion corresponding to the formula (2) is expressed in the form of the formula (4), the corresponding terms are respectively:
wherein u and v are interpolation coefficients related to formula (2), the upper corner mark represents a survey station, and the lower corner mark represents a satellite.
There are also model parameter calculationsAnd simultaneously, relating to two conditions of the formulas (1) and (2), repeating the steps, and developing according to the formula (3) to obtain b and x corresponding to the formula (4). The weight corresponding to the observation equation (4) is p, and the variance D of the double-difference residual value is set0＝a0In meters, a0The value of (b) is generally set to 0.001m to 0.005m, and the weight is given by p being 1m/D0. One double-difference residual error observed value relates to four non-difference multipath, each non-difference multipath also relates to four or three model lattice points, namely a coefficient matrix b in an expression (4) contains 12-16 elements and corresponds to 12-16 elements in X, and after an index value (namely row and column numbers) of the group of elements in the X is obtained, an n matrix and a w matrix corresponding to an equation of the following formula can be used:
wherein n is a normal equation coefficient array, and w is a normal equation constant array.
3.3, constructing a second type of normal equation based on the constraint of the model grid point value: according to the normal equation obtained in step 3.2, if data records are insufficient, and some grid point parameters to be estimated may have situations which are not related to the observation equation in step 3.2, a rank deficiency situation is bound to occur in solving the solution equation, in order to solve the problem and ensure the rationality of solving the grid point parameters, additional constraint conditions need to be added to the grid point parameters, and considering that the numerical values of multipath have a certain range, the size of the grid point parameters can be constrained first. Then for any grid point parameter Q, the observation equation at this time can be obtained according to equation (4), where:
b＝1 x＝Q s＝0 (6)
for the additional constraint, the additional constraint is used as a virtual observation value to participate in the solution, at this time, the weight p of the constraint is determined, and the variance D of the grid point parameter is set1＝a1，a1The value of (A) is generally set to 0.05m to 0.10m, and p is 1m/D1. At this time, one data record constrains one grid point parameter, so the coefficient matrix b only relates to one element, all data records are gradually traversed, constraints are added to all grid points, and an n matrix and a w matrix corresponding to a normal equation are obtained by using a formula (5).
Step 3.4, constructing a third-class equation based on the constraint of the change value between the model grid points: when the size of the lattice point parameters is constrained, the multi-path semi-celestial lattice point model established in the same environment is considered, and the lattice point parameters are not excessively mutated theoretically, so that the variation values among the lattice point parameters are constrained (including the constraint in the longitude direction and the latitude direction). When such constraints are added, there are several cases:
FIG. 3(a) shows a normal grid point (with grid point parameters of) And the surrounding four grid points (the grid point parameter is set as) The observation equation at this time can be obtained from equation (4), where:
FIG. 3(b) shows the bottom edge of the mesh model, a mesh point (with the mesh point parameters as) Three grid points related to the surroundings (set grid point parameters as) With the constraints in between, equation (4) may result in the observation equation at this time, where:
FIG. 3(c) is a diagram of constraints existing between vertices and surrounding points of the mesh point model, where the observation equation can be obtained according to equation (4), where:
wherein,and (3) expressing the nth grid point parameter on the observation station m, wherein the value of n is determined according to specific grid division.
Similar to step 3.3, the weights p of such constraints corresponding to the observation equations need to be determined. Variance D of grid point parameter variation value2＝d×a2D is the grid point value spacing (in degrees) in step 2.2, a2The value of (A) is generally set to 0.01 m/degree to 0.05 m/degree, and the weight is given as p being 1m/D2At this time, one data record constrains two grid point values, so the coefficient matrix b relates to two elements, all the data records are gradually traversed, constraints are added among all the grid points, and an n matrix and a w matrix corresponding to a normal equation are obtained by using a formula (5).
And 4, superposing the normal equations obtained in the step 3 to obtain a final normal equation, solving the parameters of the half-celestial sphere grid point model of each station, and correcting errors caused by the multipath effect in the subsequent observed values by using the model.
The invention further provides that the multi-path error of the satellite b on the station a needs to be corrected, a grid point model of the station a is found, a grid where the satellite b is located is determined, four grid point values or three grid point values in the grid where the satellite b is located are known, the error value of the multi-path effect of the position where the satellite b is located is obtained according to a corresponding interpolation coefficient, and then the satellite b is subjected to multi-path error correctionThe observed value of the satellite b at the observation station a may be corrected as a correction number.
Step 4 of the embodiment specifically includes:
step 4.1, step 3.2 to step 3.4, respectively obtaining corresponding N and W matrixes, wherein the N and W matrixes are matrixes related to a normal equation when the X parameter is solved, assuming that the matrixes corresponding to the final parameter vector X are respectively N and W, positions of the N and W matrixes in the N matrix and the W matrix can be recorded according to index values corresponding to the X parameter when the step 3.2 to step 3.4 is carried out, the final N matrix and the final W matrix are obtained by overlapping and updating the N and W matrix elements, namely the final normal equation coefficient matrix and the final normal equation constant matrix,
utilizing the following formula:
X＝N-1W (10)
x, namely the parameters of the semi-celestial sphere lattice point model of each station can be obtained through solving.
Step 4.2, when a semi-celestial grid point model is used for multi-path error elimination, the implementation steps are similar to those of step 3.1, after the relevant information of each satellite of each survey station is obtained, indexes can be carried out according to the grid point model of different stations according to the altitude angle and the azimuth angle, and the hypothesis that the satellite on the m stations needs to be subjected to satellite positioning is thatAnd (4) the satellite j corrects the multipath error, finds the grid point model of the m stations, and determines the grid where the satellite j is positioned, wherein the four grid point values in the grid where the j is positionedAndor three grid point valuesAndis known as the error value of multipath effect at the location of satellite jCan be obtained according to the formula (1) or (2), and thenAnd correcting the observed value of the satellite j at the observation station m as a correction number.
And 4.3, correcting the multipath effect of the observed value by adopting the method for all the visible satellites of all the stations, namely repeating the step 4.2 for all the visible satellites on all the stations, obtaining the multipath effect error value as a correction number for correction, and effectively weakening the multipath error influence in the subsequent GNSS data processing process.
In specific implementation, the above process can adopt computer software to realize automatic operation process.
The specific embodiments described herein are merely illustrative of the spirit of the invention. Various modifications or additions may be made to the described embodiments or alternatives may be employed by those skilled in the art without departing from the spirit or ambit of the invention as defined in the appended claims.
Claims (4)
1. A GNSS multi-path effect correction method based on a semi-celestial sphere grid point model is characterized in that: comprises the following steps of (a) carrying out,
step 1, processing double-difference observation value residual information of a fixed period of ambiguity between stations after single-day solution of multi-station data, realizing the following,
step 1.1, setting an index list according to data containing ambiguity fixing information obtained by resolving, wherein each index relates to two survey stations and two satellites and is respectively set as a survey station m and n and a satellite j and k;
step 1.2, processing the resolved residual error information aiming at each index, acquiring a double-difference residual error formed by two co-view satellites on two stations related to each index, forming a data record, and simultaneously recording an epoch serial number and the altitude angle and the azimuth angle of the two co-view satellites on the two stations;
step 2, parameterization representation of the semi-celestial grid point model at each survey station is realized as follows,
step 2.1, grid point range value setting is carried out on the semi-celestial sphere model on each survey station, semi-celestial sphere grid points are divided according to the altitude angle as latitude, and the azimuth angle as longitude;
step 2.2, setting a half celestial sphere model grid point interval d on each test station;
step 2.3, taking the grid point value of the model on each measuring station as an unknown quantity to be solved;
step 3, using the obtained double-difference observation value residual error as an observation value, and constructing a method equation corresponding to a half-celestial sphere lattice point model parameter, wherein the method equation comprises the steps of mapping the recorded residual error to the half-celestial sphere lattice network on the related station according to the altitude angles and the azimuth angles of m and n satellites j and k of the station in record by sequentially using double-difference observation value residual error information corresponding to a ambiguity fixed time period to obtain a first type of method equation; constructing a second type of normal equation based on the constraint of the model grid point value; constructing a third-class equation based on the constraint of the change value between the points of the model grid;
and 4, superposing the three types of equation obtained in the step 3 to obtain a final equation, solving the half-celestial sphere grid point model parameters of each station, and correcting errors caused by multipath effects in subsequent observed values.
2. The GNSS multipath effect correcting method based on the semi-celestial grid point model of claim 1, wherein: and mapping the recorded residual errors into the semi-celestial grid network on the related survey station according to the recorded altitude angles and azimuth angles of the satellites j and k at the survey stations m and n, wherein the model parameter calculation values between the station m and the satellite j, between the station m and the satellite k, between the station n and the satellite j and between the station n and the satellite k are respectively solved
3. The GNSS multipath effect correcting method based on the semi-celestial grid point model as claimed in claim 2, wherein: calculating model parameter calculation value between station a and satellite ba, m or n, b, j or k, as follows,
acquiring altitude angle and azimuth angle data of a satellite b on a station a, dividing according to grid points, indexing the position of a puncture point of the satellite b in a semispherical grid at the moment of recording the data, and setting longitude and latitude (x, y), wherein a model parameter calculation value in a semispherical grid point model isRepresenting a multipath error experienced by a non-differential observation between station a and satellite b;
when the multipath error between the station a and the satellite b relates to the values of four points on the grid model, extracting an interpolation coefficient by using a bilinear interpolation method; when the multipath error between the station a and the satellite b relates to the values of three points on the grid model, extracting an interpolation coefficient by using a triangular interpolation method;
expressing model parameter calculation values based on interpolation coefficients
4. The GNSS multi-path effect correction method based on semi-celestial grid point model as claimed in claim 3, wherein: step 4, it is required to correct the multipath error of the satellite b on the station a, find the grid point model of the station a, determine the grid where the satellite b is located, because the four grid point values or three grid point values in the grid where the satellite b is located are known, the error value of the multipath effect of the position where the satellite b is located is obtained according to the corresponding interpolation coefficient, and then the error value of the multipath effect of the satellite b is obtainedThe observed value of the satellite b at the observation station a may be corrected as a correction number.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811527763.6A CN109541647B (en) | 2018-12-13 | 2018-12-13 | GNSS multi-path effect correction method based on semi-celestial sphere grid point model |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811527763.6A CN109541647B (en) | 2018-12-13 | 2018-12-13 | GNSS multi-path effect correction method based on semi-celestial sphere grid point model |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109541647A true CN109541647A (en) | 2019-03-29 |
| CN109541647B CN109541647B (en) | 2019-12-10 |
Family
ID=65855088
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201811527763.6A Active CN109541647B (en) | 2018-12-13 | 2018-12-13 | GNSS multi-path effect correction method based on semi-celestial sphere grid point model |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109541647B (en) |
Cited By (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111596321A (en) * | 2020-05-29 | 2020-08-28 | 武汉大学 | Multi-GNSS multi-path error star day filtering method and system using non-difference correction |
| CN112433240A (en) * | 2020-10-13 | 2021-03-02 | 武汉理工大学 | Phase multipath extraction and correction method based on non-differential non-combination PPP model |
| CN112612039A (en) * | 2020-12-23 | 2021-04-06 | 武汉大学 | GNSS indirect signal detection and elimination method and system for static survey station |
| CN114488228A (en) * | 2022-04-11 | 2022-05-13 | 南京北斗创新应用科技研究院有限公司 | GNSS multi-path error weakening method suitable for dynamic carrier platform |
| CN115343734A (en) * | 2022-10-13 | 2022-11-15 | 武汉地铁集团有限公司 | GNSS deformation monitoring method based on bilinear interpolation hemisphere model |
| CN118311615A (en) * | 2024-03-29 | 2024-07-09 | 武汉大学 | GNSS troposphere and multipath joint modeling correction method and system |
Citations (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US5600329A (en) * | 1995-06-30 | 1997-02-04 | Honeywell Inc. | Differential satellite positioning system ground station with integrity monitoring |
| US6594582B1 (en) * | 1999-05-14 | 2003-07-15 | The United States Of America As Represented By The Administrator Of The National Aeronautics And Space Administration | GPS compound eye attitude and navigation sensor and method |
| JP2005037241A (en) * | 2003-07-14 | 2005-02-10 | Sony Corp | Receiving apparatus of gps |
| US20060140254A1 (en) * | 2004-12-29 | 2006-06-29 | Nokia Corporation | Multi-path detection method for CDMA receivers |
| US20070211793A1 (en) * | 2006-03-09 | 2007-09-13 | Shaowei Han | Multipath error estimation in satellite navigation receivers |
| CN101903794A (en) * | 2007-11-12 | 2010-12-01 | 高通股份有限公司 | Suppression of multipath effects for received SPS signal |
| CN104122566A (en) * | 2014-07-01 | 2014-10-29 | 华东师范大学 | Multi-path error removing method of navigation satellite system and multi-path hemisphere model |
| CN106643709A (en) * | 2016-10-10 | 2017-05-10 | 东南大学 | Combined navigation method and device for offshore carrier |
| CN106646538A (en) * | 2016-10-31 | 2017-05-10 | 东南大学 | Single-difference filtering-based deformation monitoring GNSS (global navigation satellite system) signal multi-path correction method |
| CN107193016A (en) * | 2017-04-06 | 2017-09-22 | 广州中硕云空间信息技术有限公司 | A kind of method and system of city GNSS navigation quality evaluations and prediction |
-
2018
- 2018-12-13 CN CN201811527763.6A patent/CN109541647B/en active Active
Patent Citations (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US5600329A (en) * | 1995-06-30 | 1997-02-04 | Honeywell Inc. | Differential satellite positioning system ground station with integrity monitoring |
| US6594582B1 (en) * | 1999-05-14 | 2003-07-15 | The United States Of America As Represented By The Administrator Of The National Aeronautics And Space Administration | GPS compound eye attitude and navigation sensor and method |
| JP2005037241A (en) * | 2003-07-14 | 2005-02-10 | Sony Corp | Receiving apparatus of gps |
| US20060140254A1 (en) * | 2004-12-29 | 2006-06-29 | Nokia Corporation | Multi-path detection method for CDMA receivers |
| US20070211793A1 (en) * | 2006-03-09 | 2007-09-13 | Shaowei Han | Multipath error estimation in satellite navigation receivers |
| CN101903794A (en) * | 2007-11-12 | 2010-12-01 | 高通股份有限公司 | Suppression of multipath effects for received SPS signal |
| CN104122566A (en) * | 2014-07-01 | 2014-10-29 | 华东师范大学 | Multi-path error removing method of navigation satellite system and multi-path hemisphere model |
| CN106643709A (en) * | 2016-10-10 | 2017-05-10 | 东南大学 | Combined navigation method and device for offshore carrier |
| CN106646538A (en) * | 2016-10-31 | 2017-05-10 | 东南大学 | Single-difference filtering-based deformation monitoring GNSS (global navigation satellite system) signal multi-path correction method |
| CN107193016A (en) * | 2017-04-06 | 2017-09-22 | 广州中硕云空间信息技术有限公司 | A kind of method and system of city GNSS navigation quality evaluations and prediction |
Non-Patent Citations (2)
| Title |
|---|
| XUAN ZOU等: "Instantaneous BDS + GPS undifferenced NRTK positioning with dynamic atmospheric constraints", 《GPS SOLUTIONS》 * |
| 李星星: "GNSS精密单点定位及非差模糊度快速确定方法研究", 《中国博士学位论文全文数据库基础科学辑》 * |
Cited By (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111596321A (en) * | 2020-05-29 | 2020-08-28 | 武汉大学 | Multi-GNSS multi-path error star day filtering method and system using non-difference correction |
| CN111596321B (en) * | 2020-05-29 | 2022-04-01 | 武汉大学 | Multi-GNSS multi-path error star day filtering method and system using non-difference correction |
| CN112433240A (en) * | 2020-10-13 | 2021-03-02 | 武汉理工大学 | Phase multipath extraction and correction method based on non-differential non-combination PPP model |
| CN112612039A (en) * | 2020-12-23 | 2021-04-06 | 武汉大学 | GNSS indirect signal detection and elimination method and system for static survey station |
| CN112612039B (en) * | 2020-12-23 | 2023-08-15 | 武汉大学 | GNSS non-direct signal detection and elimination method and system for static station |
| CN114488228A (en) * | 2022-04-11 | 2022-05-13 | 南京北斗创新应用科技研究院有限公司 | GNSS multi-path error weakening method suitable for dynamic carrier platform |
| CN114488228B (en) * | 2022-04-11 | 2022-07-01 | 南京北斗创新应用科技研究院有限公司 | GNSS multi-path error weakening method suitable for dynamic carrier platform |
| WO2023197714A1 (en) * | 2022-04-11 | 2023-10-19 | 南京北斗创新应用科技研究院有限公司 | Gnss multi-path error reducing method suitable for dynamic carrier platform |
| CN115343734A (en) * | 2022-10-13 | 2022-11-15 | 武汉地铁集团有限公司 | GNSS deformation monitoring method based on bilinear interpolation hemisphere model |
| CN118311615A (en) * | 2024-03-29 | 2024-07-09 | 武汉大学 | GNSS troposphere and multipath joint modeling correction method and system |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109541647B (en) | 2019-12-10 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109541647B (en) | GNSS multi-path effect correction method based on semi-celestial sphere grid point model | |
| CN110058236B (en) | InSAR and GNSS weighting method oriented to three-dimensional surface deformation estimation | |
| US10883829B2 (en) | Systems and methods for GNSS SNR probabilistic localization and 3-D mapping | |
| US20240159915A1 (en) | System and method for determining gnss positioning corrections | |
| Bisnath et al. | Precise point positioning | |
| CN105929424A (en) | BDS/GPS high-accuracy positioning method | |
| CN111596321B (en) | Multi-GNSS multi-path error star day filtering method and system using non-difference correction | |
| CN101743453A (en) | The post-mission high accuracy position and azimuth determining system | |
| CN114488228B (en) | GNSS multi-path error weakening method suitable for dynamic carrier platform | |
| WO2017160617A1 (en) | Navigation satellite wide-lane bias determination system and method | |
| CN116391138A (en) | Positioning method, device, equipment, system and storage medium | |
| CN114646992A (en) | Positioning method, positioning device, computer equipment, storage medium and computer program product | |
| CN115343734B (en) | GNSS deformation monitoring method based on bilinear interpolation hemisphere model | |
| Brack et al. | Operational multi-GNSS global ionosphere maps at GFZ derived from uncombined code and phase observations | |
| CN112068173A (en) | Collaborative navigation method based on loop and product data association algorithm | |
| CN105527639A (en) | Satellite positioning method based on smoothness and extrapolation | |
| CN116299574B (en) | GLONASS occultation corresponding reference star PRN correction method based on altitude angle | |
| Chu et al. | A new approach to modernized GPS phase-only ambiguity resolution over long baselines | |
| CN113777633A (en) | Positioning method, electronic device and computer storage medium | |
| CN116931007B (en) | Ionosphere delay processing method, ionosphere delay processing device, ionosphere delay processing equipment and storage medium | |
| CN115407367B (en) | Method for estimating navigation positioning precision attenuation factor of mixed constellation satellite | |
| Noinak et al. | Testing Horizontal Coordinate Correction Model Used for Transformation from PPP GNSS Technique to Thai GNSS CORS Network Based on ITRF2014 | |
| Kelecy | Precise relative sea surface positioning of a floating platform using the Global Positioning System (GPS) | |
| Alojaiman | Tropospheric delay modeling using gnss observations from continuously operating reference stations (cors) | |
| Park | Determination of glacial isostatic adjustment parameters based on precise point positioning using GPS |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |