Beidou-based engineering measurement positioning method and system

Abstract

The invention discloses an engineering surveying and positioning method and system based on Beidou. The method comprises the following steps: S1, establishing an oil field operation environment three-dimensional model, and performing satellite visibility and multipath prediction; according to the invention, by constructing a multi-path geometric model and a dual-antenna array discrimination and weighted suppression algorithm, reflection interference caused by a metal facility dense environment in a well site and a station site is effectively weakened; by fusing a local DEM terrain model and earth surface reflection characteristics and introducing vertical Kalman filtering, the stability and precision of Beidou elevation calculation are significantly improved; under the scene of insufficient visibility such as shielding of a drilling machine and large equipment, continuous and high-reliability positioning capability is kept through a tight coupling fusion strategy of IMU and vision assistance, so that high-robustness and high-precision Beidou engineering measurement positioning can be realized in a complex environment of an oil field.

Description

Technical Field

[0001] This invention relates to the field of engineering surveying and positioning technology, and in particular to an engineering surveying and positioning method and system based on BeiDou. Background Technology

[0002] Oilfield surface engineering typically includes well site layout, station process construction, gathering and transmission pipeline laying, road construction, drainage and ancillary facility construction, etc. These projects are large-scale and complex, requiring extremely high accuracy and stability in measurement and positioning. While existing engineering measurement and positioning methods based on BeiDou (BDS) can achieve centimeter-level positioning, they still have the following prominent problems in the oilfield environment: 1. Dense metal facilities at well sites and stations lead to severe multipath effects: There are a large number of metal pipelines, wellhead devices, tanks, separators, and station steel structures in oilfield well sites. Beidou signals are easily superimposed on metal reflections, causing carrier phase jumps and increased multipath errors, resulting in unstable RTK fixed solutions and drifting positioning coordinates. 2. Oilfields are mostly characterized by undulating terrain, desert and Gobi environments, resulting in large BeiDou elevation errors: In large terrain scenarios with poor surface reflection of sand and gravel and drastic changes in satellite elevation angle, the vertical error of BeiDou is more obvious than the horizontal error, which cannot meet the construction requirements such as wellhead elevation, pipeline slope control, and site leveling. 3. Severe obstruction by drilling rigs, joint stations, and large machinery leads to insufficient satellite signal field of view: Oilfield operation areas are generally occupied by large equipment such as drilling rig towers, workover rigs, large vehicles, sand trucks, and fracturing trucks. These tall obstructions block satellite signals, reducing the number of available satellites for the BeiDou system and decreasing RTK calculation efficiency. In summary, existing technologies cannot guarantee stable, reliable, and high-precision BeiDou measurement and positioning in the complex environment of oil fields. Therefore, there is an urgent need to propose a BeiDou-based engineering measurement and positioning method and system that can suppress multipath errors, improve elevation accuracy, and enhance usability in obstructed environments. Summary of the Invention

[0003] Based on the technical problems existing in the background technology, this invention proposes an engineering measurement and positioning method and system based on Beidou, which solves the problems of the dense metal facilities in well sites and stations leading to severe multipath effects, oil fields mostly having undulating terrain, desert and Gobi environments, large Beidou elevation errors, and insufficient satellite signal field of view due to severe obstruction by drilling rigs, joint stations and large machinery.

[0004] The present invention proposes an engineering measurement and positioning method based on BeiDou, comprising the following steps: S1: Establish a three-dimensional model of the oilfield operating environment and perform satellite visibility and multipath prediction; S2: Use a GNSS receiver that supports multi-system reception to receive signals from BeiDou and at least one other satellite navigation system and record raw observations; S3: The received signal is multipath detected and weighted by a dual-antenna array and a multipath discrimination algorithm to obtain the corrected pseudorange and carrier phase observations. S4: Based on the observations output by S3, the BeiDou RTK elevation is model compensated and vertical Kalman filtered by combining the local digital elevation model (DEM) and surface reflection characteristics. S5: When satellite obstruction or insufficient visibility occurs, IMU and visual assistance are activated and tightly coupled with GNSS for filtering and fusion calculation. S6: Performs measurement task control, mode switching, and alarms based on real-time quality indicators; S7: Output engineering layout results, including well site coordinates, pipeline center / side stakes, station equipment layout coordinates, elevation and quality report.

[0005] Preferably, in S1, the established three-dimensional model of the oilfield operating environment includes metal pipelines, tanks, drilling rig towers, station steel structures, roads and terrain elevation data, and generates satellite visibility occlusion maps and multipath reflection candidate sets at each observation point based on the three-dimensional model of the oilfield operating environment. The specific logical steps of S1 are as follows: S101: Collect spatial information of well sites, stations and pipeline areas, including topographic elevation, location and size of metal pipelines, drilling rig towers, station steel structures, three-dimensional dimensions and location information of equipment tank facilities, and location of roads, construction platforms and temporary facilities; S102: The collected data are fused into a unified coordinate system, and a 3D model of the oilfield is constructed using CAD modeling software. The 3D model of the oilfield includes a terrain elevation model, a 3D facility model, and road and construction platform models. S103: Based on the engineering layout requirements, determine the observation points for wellheads, pipeline center stakes, side stakes, and station layout points, and mark each observation point in the 3D model for subsequent satellite visibility and multipath analysis; S104: Perform ray tracing for each observation point, determine whether the satellite signal is obstructed by three-dimensional facilities or terrain, and output a satellite visibility occlusion map; S105: For each observation point, use a 3D model to calculate possible satellite signal reflecting surfaces. Satellite signal reflecting surfaces include, but are not limited to, metal pipelines, tanks and steel structures with high reflectivity, as well as ground or road reflecting surfaces. Calculate the geometric relationship between the satellite and the reflecting surface to generate a multipath reflection candidate set, which serves as the input for subsequent dual-antenna multipath discrimination or weight allocation. S106: Optimize the selection of observation points based on satellite visibility and multipath prediction results, avoid high multipath or severely occluded locations, and mark points that must be measured in severely occluded areas as having IMU / visual assistance enabled.

[0006] Preferably, in step S104, when performing satellite visibility analysis, it is necessary to calculate the satellite-to-observation-point vector: based on the observation point and satellite position Calculate the vector: ,in Let be the vector pointing from the observation point to the i-th satellite; When performing occlusion detection, the surface normal vectors of objects in the scene are used. and observation points Indicates the reflecting / obstructing surface, used to determine whether a satellite is obstructed, and the ray... If the ray intersects the surface of the object, the satellite will be blocked, where t represents the position of a point on the satellite signal path, a scalar that determines the direction vector along the ray originating from the observation point. How far has been traveled? Its determination function is: ; in Indicates a set of visible satellites. , and Let be the coordinates of the three vertices that constitute the surface of the j-th object.

[0007] Preferably, the specific logical steps of S105 are as follows: S1051: Calculate the patch normal and unit normal: for a patch ( ), , ,in For dough pieces The normal vector, It is the unit normal vector; S1052: Use the method of mirrors to find possible specular reflection points. : Put the satellite About the dough Mirroring , For the i-th satellite Regarding the j-th facet The coordinates of the mirror point are obtained by first taking a point on the plane and denoting it as... The vector projection of the satellite onto the plane is ,in For satellite up to the j-th face The length of the vector projection; Mirror point: Next, find the ray. Intersection parameters with the plane , Intersection coordinates: Next, determine whether the intersection point lies within the face. If the intersection point lies within the face, then set... Otherwise, it is assumed that the surface does not produce a specular reflection point; in Indicates from the observation point Point to mirror point The ray equation, For the nth piece The coordinates of the reflection point that produces effective specular reflection; S1053: Verifying the geometry of specular reflection: Defining the incident direction Define the launch direction And calculate the angle difference between the two: ; like If the mirror geometry matches well, then the surface is considered to be a good mirror reflector; otherwise, it can be considered a non-primary mirror reflector. The reference threshold angle is used to determine whether it is a valid specular reflective surface; S1054: Calculate path increment: Calculate the additional path length generated by multipathing. The corresponding time delay is: Where c is the speed of electromagnetic wave propagation, and takes the value of ; Then calculate the angle of incidence. , ; S1055: Estimating Reflection Intensity and Weight , , The larger the value, the better the quality of the surface. For satellite The stronger the multipath impact, the more likely it is to be identified, weighted down, or eliminated. in Let be the reflection coefficient of the j-th surface patch. This is the path length attenuation coefficient; S1056: Satellite All that meet the threshold condition constituting a satellite Multipath reflection candidate set , This information is then used as input for subsequent dual-antenna discrimination and observation weighting / rejection.

[0008] Preferably, in step S2, the original observations include pseudorange observations. Carrier phase observation Doppler frequency shift observation Signal strength observation And channel quality metrics, among which Indicates the satellite number, Indicates frequency point; ,in For geometric distance, , For the clock difference between the receiver and the satellite, For ionospheric delay, For tropospheric delay, For pseudo-range multipath, For the purpose of observing noise, c is the propagation speed of electromagnetic waves; ,in The wavelength at that frequency point, For carrier cumulative phase, The initial ambiguity; Its actual recorded content includes Doppler frequency shift values ​​and carrier lock status; ,in Indicates that satellite i is at frequency point The received carrier power is estimated by the correlator output. represents the thermal noise power spectral density, represents the noise intensity per unit bandwidth, is the noise floor estimated by the receiver front end, and 10 is the scaling factor when converting the linear power ratio to dB. ( () is a base-10 logarithm used to convert the linear power ratio to dB-Hz units.

[0009] Preferably, the specific logical steps of S3 are as follows: S301: Obtain pseudorange observations by simultaneously receiving satellite signals from the same epoch via the main antenna A1 and the secondary antenna A2. Carrier phase observation Doppler frequency shift observation Signal strength observation and channel quality metrics; S302: Calculate the satellite direction unit vector based on the receiver coordinates and satellite coordinates. ,in Let be the unit direction vector pointing from the observation point to the i-th satellite; Furthermore, the baseline vectors of the two antennas are obtained based on the antenna installation parameters: ,in These represent the three-dimensional coordinates of the main antenna A1 and the secondary antenna A2, respectively. S303: Calculate the ideal geometric path difference based on the satellite orientation and baseline vector, and convert it into the theoretical carrier phase difference. ,in Let be the ideal theoretical carrier phase difference corresponding to the i-th satellite; S304: Differentiate the carrier phases of two antennas at the same frequency to obtain the observed phase difference. ; S305: Compare the observed phase difference with the theoretical phase difference to obtain the phase residual reflecting the effects of multipath propagation. ; S306: Obtain pseudorange residuals from epoch filtering or differential calculation. ,in Predict pseudorange for the model; S307: Calculate the drop based on the current signal-to-noise ratio and the average signal-to-noise ratio of the sliding window. The higher this index, the higher the probability that the signal is blocked or interfered with. For the i-th satellite at frequency point The average carrier-to-noise ratio of the sliding window at that location. For the i-th satellite at frequency point The current carrier-to-noise ratio at the location; S308: Combine phase residuals, pseudorange residuals, SNR decrease, and integer jumps according to weights. ,in For weight parameters, To convert the phase residual to meters or cycles, For phase residual, For pseudo-distance residual index, This is a carrier integer cycle transition indicator; 0 indicates no transition, and 1 indicates a transition. For frequency point The wavelength; S309: A joint suppression strategy of SNR and MP is adopted to obtain the weight of each satellite observation. ,in Can be made by C / Convert to linear value As a regulating factor, ,when Hour, When it is close to 1, When it is big, Close to 0, among which For the i-th satellite, the multipath comprehensive index is used. S310: Weighted correction of satellite observation noise When severe multipathing exists, the weight is significantly reduced; if If so, the observations of that satellite are directly discarded. This represents the original observation noise variance. Let $\mathbf{i}$ be the variance of the observation noise of the $i$-th satellite after weight correction. Thresholds are removed from the observations; S311: Input the retained observations and their corresponding noise matrices into the Kalman filter to complete the RTK solution for this epoch.

[0010] Preferably, the specific logical steps of S4 are as follows: S401: Ensure the elevation directly provided by the RTK solution. Elevation values ​​obtained by interpolation at the observation points from the DEM If the reference systems are not the same, a datum transformation is performed. If the DEM resolution is low, bilinear interpolation is performed around the observation point. ; S402: Calculate the initial deviation of DEM-RTK ; S403: Calculate the average signal quality of the satellite signal at the observation point. Where N is the number of satellites participating in the averaging. For the average carrier-to-noise ratio, C / Let C be the carrier-to-noise ratio of the i-th satellite, and C be the signal carrier power. For noise power spectral density, C / This describes the ratio of signal strength to noise for a single satellite. The reflectance factor r is estimated based on the average SNR. , where parameters The parameters of the sigmoid function are set empirically. Constructing a surface reflection correction term by combining with DEM-RTK bias , ,in This is the scaling factor; S404: Structural Integrated Elevation Correction Observation , ,in and To normalize the weighting, a confidence level is set, which is determined by the accuracy of the DEM and the reliability of the reflection term. If the DEM accuracy is high, then... If the SNR is significantly low and the multipath probability is high, increase the value of the multipath probability. ; S405: Establish the state-space model of the vertical Kalman filter (V-KF): ,in Vertical height For the height deviation rate, Let V be the state vector of the vertical Kalman filter V-KF at time k. S406: Determine the observation noise covariance R: ,in For RTK observation uncertainty variance, For the variance of the DEM correction term, ,in This refers to the accuracy variance of the DEM itself. S407: Perform Kalman filtering time update, predict state and covariance. , ,in Let F be the filter state vector, F be the state transition matrix, and Q be the process noise covariance. Let K be the predicted state vector at time k. Let be the predicted state covariance matrix at time k; S408: Smoothed elevation output after filtering , ,in Let k be the vertical height obtained after filtering at time k. And record the height error estimate. ; S409: If If so, mark the DEM as abnormal and reduce it. ,like Then increase If the filtering innovation exceeds the threshold, a retest will be triggered or the system will fall back to RTK-only mode and logs will be recorded. The threshold for DEM anomaly detection.

[0011] Preferably, the specific logical steps of S5 are as follows: S501: When GNSS availability is insufficient, IMU pre-integration and visual odometry are used as short-term alternatives and constraints. Acceleration and angular velocity data of IMU, image frames and feature points of visual sensor are collected respectively, and noise processing, time synchronization and coordinate system unification are performed. Insufficient GNSS availability is manifested as a small number of satellites, high PDOP or poor observation quality. S502: By integrating the acceleration and angular velocity of the IMU, the predicted values ​​of the current position, velocity, and attitude are obtained, forming the prior state for fusion calculation. The formula used for prediction update is: ,in Let be the posterior state at time k. The filtering state at time k-1, ( ) is an IMU prediction model. For the raw IMU measurements, The IMU sampling period; S503: Extracts current position observations from pseudorange, carrier phase, and differential positioning results generated by BeiDou and other available GNSS satellites, and forms a visual observation vector from the relative displacement or attitude change output by the visual odometry. The two together constitute the fused solution of the observation value. ; S504: Calculate the fusion residual based on the difference between the IMU predicted state and the GNSS / visual observations, and obtain the Kalman gain for posterior state update. This yields high-precision position, velocity, and attitude estimates after multi-source information fusion, among which... Let Kalman gain be at time k. For observation models, The observation vectors for fusion solution; S505: The fused positioning result is used as the current accurate state of the system. It is used to output a stable solution from the measurement operation control module and fed back to the next epoch as an initial condition, so that the system can maintain continuous, stable and high-precision positioning capability even when satellite visibility is insufficient.

[0012] Preferably, in step S6, the real-time quality indicators include GNSS satellite quality indicators, multipath and occlusion indicators, IMU stability indicators, and fusion solution reliability indicators. All indicators can be uniformly expressed as: ,in Let k be the real-time quality index vector at time k; The specific logical steps are as follows: S601: Acquires GNSS, IMU, and visual data, and calculates the number of satellites. With position accuracy factor Carrier-to-noise ratio Multi-path intensity comprehensive index residual variance IMU bias drift and fusion covariance ; S602: Determines whether the current RTK can work normally. Its execution conditions are: If all conditions are met, the system enters high-precision RTK mode; if any condition is not met, the system enters weak field-of-view mode. in The threshold for the number of visible satellites. This is the threshold for the position accuracy factor. The carrier-to-noise ratio threshold. For multipath intensity threshold, The threshold for GNSS residual variance; S603: Meets all conditions in S602, triggering high-precision RTK mode. In high-precision RTK mode, GNSS is the primary method, with IMU / vision as a secondary method. The triggering conditions for weak field-of-view mode are: When the weak field-of-view mode is triggered, the IMU weights are automatically increased, and visual constraints are used to maintain horizontal accuracy. The triggering condition for the GNSS interruption mode is as follows: Then it completely switches to IMU + visual navigation, and GNSS enters the fusion filter with zero weight; S604: Execute alarm logic: Audio-visual prompts; The system prompts for mobile base station / equipment adjustment. Warning: Degraded positioning reliability. The system prompts for reinitialization. For RTK, a fixed solution ratio is used; S605: When the quality indicators are normal, the system performs measurement tasks, including controlling the sampling frequency, acquiring control points, and measuring the control trajectory. When the indicators decrease, the system delays sampling to maintain reliable data results.

[0013] This invention also proposes a BeiDou-based engineering measurement and positioning system, including a receiving unit, a dual-antenna multipath detection unit, a DEM and environmental model database, an IMU / visual assistance unit, a differential enhancement unit, a data fusion and solution unit, a measurement operation control unit, and an output unit; The receiving unit is used to receive observations from BeiDou and at least one other GNSS system. The dual-antenna multipath detection unit includes hardware / software modules for SNR acquisition, phase difference calculation, multipath degree calculation, and reflection source matching based on the environment model, and is used for multipath identification and elimination based on the phase and amplitude difference of the dual antennas. The DEM and environmental model database is used to store the three-dimensional structure and topographic elevation of the oilfield and provide elevation compensation data. The DEM and environmental model database supports offline import and rapid on-site modeling. The model can be updated by laser scanning or UAV aerial survey to adapt to construction changes. The IMU / visual aid unit includes a three-axis MEMS or fiber optic gyroscope IMU, a visual camera, and a visual inertial navigation unit. It also performs online estimation and correction of the IMU bias to reduce cumulative errors and provides short-term positioning when GNSS is blocked. The differential enhancement unit is used to access CORS / RTK or PPP enhancement services, and supports VRS / network RTK, base station differential, and real-time differential data processing interface based on CMR / RTCM. The data fusion and solution unit is used to perform weighted fusion of various observations and output engineering measurement results; The measurement operation control unit is used to perform real-time quality assessment and measurement strategy switching; The output unit is used to integrate the multi-source fusion positioning / elevation data and real-time quality indicators, and output the final coordinate results and quality report.

[0014] Compared with existing technologies, the beneficial effects of this invention are: 1. By constructing a three-dimensional environment model of the well site / station and predicting satellite obstruction and reflection paths, combined with the multi-path direction discrimination and weight suppression algorithm of the dual-antenna array, the carrier phase jump caused by reflection from metal facilities can be significantly reduced, thereby improving the fixed solution success rate and solution stability of BeiDou RTK. 2. Based on the local digital elevation model (DEM) and surface reflection characteristics, a BeiDou elevation error compensation model is constructed, and a vertical Kalman filter is used to dynamically constrain the elevation. This can effectively reduce the vertical error of BeiDou under extreme surface conditions such as deserts and Gobi, and improve the elevation accuracy from decimeter level to centimeter level. This enables the BeiDou model to maintain centimeter-level accuracy and high stability in typical oilfield environments with strong reflection, high obstruction, and large topographic undulations. 3. To address the insufficient satellite field of view caused by drilling rig towers, workover rigs, sand-adding trucks, fracturing trucks, etc., the system automatically switches to GNSS / IMU / visual tightly coupled mode when the BeiDou signal is weak or temporarily lost. Continuous positioning is maintained through high-frequency inertial navigation compensation and visual feature constraints, avoiding calculation interruptions and coordinate jumps, thereby improving system availability and anti-blocking capabilities. 4. The positioning quality is evaluated in real time based on indicators such as the number of visible satellites, carrier-to-noise ratio, multipath intensity, fusion covariance, and IMU bias drift. It can automatically perform multi-mode switching such as high-precision RTK, weak field-of-view GNSS / IMU fusion, and GNSS interruption visual inertial navigation, and trigger stability alarms, thereby ensuring the continuity and reliability of measurement tasks under complex working conditions. This invention effectively reduces reflection interference caused by dense metal facilities in well sites and stations by constructing a multi-path geometric model, using dual-antenna array discrimination and weighted suppression algorithms. By fusing local DEM terrain models and surface reflection features and introducing vertical Kalman filtering, it significantly improves the stability and accuracy of BeiDou elevation calculation. In scenarios with insufficient visibility, such as when drilling rigs or large equipment obstruct the view, it maintains continuous and highly reliable positioning capabilities through a tightly coupled fusion strategy of IMU and visual assistance. This enables highly robust and high-precision BeiDou engineering measurement and positioning in complex oilfield environments, significantly improving the problems of traditional RTK such as easy loss of lock, easy drift, unstable elevation, and unavailability in obstructed areas, thus significantly improving engineering measurement efficiency and construction control quality. Attached Figure Description

[0015] Figure 1 This is a flowchart of an engineering measurement and positioning method based on BeiDou proposed in this invention; Figure 2 This is a flowchart of the dual-antenna multipath discrimination and weighted elimination process in a BeiDou-based engineering measurement and positioning method proposed in this invention. Figure 3 This is a flowchart of elevation compensation and V-KF in a BeiDou-based engineering measurement and positioning method proposed in this invention. Figure 4 This is a flowchart of IMU / visual short-term substitution and fusion under occlusion conditions in a BeiDou-based engineering measurement and positioning method proposed in this invention. Figure 5 This is a schematic diagram of a typical application scenario in an oilfield using a BeiDou-based engineering measurement and positioning method proposed in this invention. Figure 6 This is a block diagram of an engineering measurement and positioning system based on BeiDou proposed in this invention. Detailed Implementation

[0016] The present invention will be further explained below with reference to specific embodiments.

[0017] Example 1 Reference Figure 1-5 This embodiment proposes an engineering measurement and positioning method based on BeiDou, including the following steps: S1: Establish a 3D model of the oilfield operating environment and perform satellite visibility and multipath prediction. The established 3D model of the oilfield operating environment includes metal pipelines, tanks, drilling rig towers, station steel structures, roads and terrain elevation data. Based on the 3D model of the oilfield operating environment, a satellite visibility occlusion map and a multipath reflection candidate set are generated at each observation point. The specific logical steps are as follows: S101: Collect spatial information of well sites, stations, and pipeline areas, including terrain elevation, location and size of metal pipelines, drilling rig towers, station steel structures, three-dimensional dimensions and location information of equipment tanks and facilities, and locations of roads, construction platforms, and temporary facilities. The location error of the collected metal pipelines must be ≤0.05m and the height error of the drilling rig towers must be ≤0.1m. If these ranges are exceeded, it will cause the subsequent S104 ray tracing occlusion judgment deviation to exceed 10%, resulting in misjudgment of satellite visibility. S102: The collected data are fused into a unified coordinate system, and a 3D model of the oilfield is constructed using CAD modeling software. The 3D model of the oilfield includes a terrain elevation model, a 3D facility model, and road and construction platform models. S103: Based on the engineering layout requirements, determine the observation points for wellheads, pipeline center stakes, side stakes, and station layout points, and mark each observation point in the 3D model for subsequent satellite visibility and multipath analysis; S104: Perform ray tracing for each observation point, determine whether the satellite signal is obstructed by three-dimensional facilities or terrain, and output a satellite visibility occlusion map; When performing satellite visibility analysis, it is necessary to calculate the satellite-to-observation-point vector: based on the observation point... and satellite position Calculate the vector: ,in Let be the vector pointing from the observation point to the i-th satellite; When performing occlusion detection, the surface normal vectors of objects in the scene are used. and observation points Indicates the reflecting / obstructing surface, used to determine whether a satellite is obstructed, and the ray... If the ray intersects the surface of the object, the satellite will be blocked, where t represents the position of a point on the satellite signal path, a scalar that determines the direction vector along the ray originating from the observation point. How far has been traveled? Its determination function is: ; in Indicates a set of visible satellites. , and Let the coordinates of the three vertices that constitute the surface of the j-th object be given. S105: For each observation point, calculate possible satellite signal reflecting surfaces using a 3D model. Satellite signal reflecting surfaces include, but are not limited to, highly reflective objects such as metal pipelines, tanks, and steel structures, as well as ground or road reflecting surfaces. Calculate the geometric relationship between the satellite and the reflecting surface to generate a multipath reflection candidate set, which serves as input for subsequent dual-antenna multipath discrimination or weight allocation. The specific logical steps are as follows: S1051: Calculate the patch normal and unit normal: for a patch ( ), , ,in For dough pieces The normal vector, It is the unit normal vector; S1052: Use the method of mirrors to find possible specular reflection points. : Put the satellite About the dough Mirroring , For the i-th satellite Regarding the j-th facet The coordinates of the mirror point are obtained by first taking a point on the plane and denoting it as... The vector projection of the satellite onto the plane is ,in For satellite up to the j-th face The length of the vector projection; Mirror point: Next, find the ray. Intersection parameters with the plane , Intersection coordinates: Next, determine whether the intersection point lies within the face. If the intersection point lies within the face, then set... Otherwise, it is assumed that the surface does not produce a specular reflection point; in Indicates from the observation point Point to mirror point The ray equation, For the nth piece The coordinates of the reflection point that produces effective specular reflection; S1053: Verifying the geometry of specular reflection: Defining the incident direction Define the launch direction And calculate the angle difference between the two: ; like If the mirror geometry matches well, then the surface is considered to be a good mirror reflector; otherwise, it can be considered a non-primary mirror reflector. The reference threshold angle is used to determine whether it is a valid specular reflective surface; S1054: Calculate path increment: Calculate the additional path length generated by multipathing. The corresponding time delay is: Where c is the speed of electromagnetic wave propagation, and takes the value of ; Then calculate the angle of incidence. , ; S1055: Estimating Reflection Intensity and Weight , , The larger the value, the better the quality of the surface. For satellite The stronger the multipath impact, the more likely it is to be identified, weighted down, or eliminated. in Let be the reflection coefficient of the j-th surface patch. This is the path length attenuation coefficient; S1056: Satellite All that meet the threshold condition constituting a satellite Multipath reflection candidate set , And it serves as the input for subsequent dual-antenna discrimination and observation weighting / rejection; S106: Based on satellite visibility and multipath prediction results, optimize the selection of observation points, avoid high multipath or severely obstructed locations, and mark points that must be measured in severely obstructed areas as having IMU / visual assistance enabled; S2: Use a GNSS receiver that supports multi-system reception to receive signals from BeiDou and at least one other satellite navigation system and record raw observations; The original observations include pseudorange observations. Carrier phase observation Doppler frequency shift observation Signal strength observation And channel quality metrics, among which Indicates the satellite number, Indicates frequency point; ,in For geometric distance, , For the clock difference between the receiver and the satellite, For ionospheric delay, For tropospheric delay, For pseudo-range multipath, For the purpose of observing noise, c is the propagation speed of electromagnetic waves; ,in The wavelength at that frequency point, For carrier cumulative phase, The initial ambiguity; Its actual recorded content includes Doppler frequency shift values ​​and carrier lock status; ,in Indicates that satellite i is at frequency point The received carrier power is estimated by the correlator output. represents the thermal noise power spectral density, represents the noise intensity per unit bandwidth, is the noise floor estimated by the receiver front end, and 10 is the scaling factor when converting the linear power ratio to dB. ( () is the base-10 logarithm, used to convert the linear power ratio to dB-Hz units; S3: The received signal is multipath detected and weighted by a dual-antenna array and a multipath discrimination algorithm to obtain the corrected pseudorange and carrier phase observations. The specific logical steps are as follows: S301: Obtain pseudorange observations by simultaneously receiving satellite signals from the same epoch via the main antenna A1 and the secondary antenna A2. Carrier phase observation Doppler frequency shift observation Signal strength observation and channel quality metrics; S302: Calculate the satellite direction unit vector based on the receiver coordinates and satellite coordinates. ,in Let be the unit direction vector pointing from the observation point to the i-th satellite; Furthermore, the baseline vectors of the two antennas are obtained based on the antenna installation parameters: ,in These represent the three-dimensional coordinates of the main antenna A1 and the secondary antenna A2, respectively. Additionally, the baseline vector b of the primary and secondary antennas should be kept horizontal (tilt angle ≤ 1°). If the tilt angle reaches 3°, it will lead to a theoretical carrier phase difference. A calculation error exceeding 5%, for example, a baseline length of 0.8m and a tilt of 3°, will result in the corresponding... Deviation 0.1rad; S303: Calculate the ideal geometric path difference based on the satellite orientation and baseline vector, and convert it into the theoretical carrier phase difference. ,in Let be the ideal theoretical carrier phase difference corresponding to the i-th satellite; S304: Differentiate the carrier phases of two antennas at the same frequency to obtain the observed phase difference. ; S305: Compare the observed phase difference with the theoretical phase difference to obtain the phase residual reflecting the effects of multipath propagation. ; S306: Obtain pseudorange residuals from epoch filtering or differential calculation. ,in Predict pseudorange for the model; S307: Calculate the drop based on the current signal-to-noise ratio and the average signal-to-noise ratio of the sliding window. The higher this index, the higher the probability that the signal is blocked or interfered with. For the i-th satellite at frequency point The average carrier-to-noise ratio of the sliding window at that location. For the i-th satellite at frequency point The current carrier-to-noise ratio at the location; S308: Combine phase residuals, pseudorange residuals, SNR decrease, and integer jumps according to weights. ,in For weight parameters, To convert the phase residual to meters or cycles, For phase residual, For pseudo-distance residual index, This is a carrier integer cycle transition indicator; 0 indicates no transition, and 1 indicates a transition. For frequency point The wavelength; S309: A joint suppression strategy of SNR and MP is adopted to obtain the weight of each satellite observation. ,in Can be made by C / Convert to linear value As a regulating factor, ,when Hour, When it is close to 1, When it is big, Close to 0, among which For the i-th satellite, the multipath comprehensive index is used. S310: Weighted correction of satellite observation noise When severe multipathing exists, the weight is significantly reduced; if If so, the observations of that satellite are directly discarded. This represents the original observation noise variance. Let $\mathbf{i}$ be the variance of the observation noise of the $i$-th satellite after weight correction. Thresholds are removed from the observations; S311: Input the retained observations and their corresponding noise matrices into the Kalman filter to complete the RTK solution for this epoch; S4: Based on the observations output by S3, and combined with the local digital elevation model (DEM) and surface reflection characteristics, model compensation and vertical Kalman filtering are performed on the BeiDou RTK elevation. In this implementation, the specific logical steps are as follows: S401: Ensure the elevation directly provided by the RTK solution. Elevation values ​​obtained by interpolation at the observation points from the DEM If the reference systems are not the same, a datum transformation is performed. If the DEM resolution is low, bilinear interpolation is performed around the observation point. ; S402: Calculate the initial deviation of DEM-RTK ; S403: Calculate the average signal quality of the satellite signal at the observation point. Where N is the number of satellites participating in the averaging process. For the average carrier-to-noise ratio, C / Let C be the carrier-to-noise ratio of the i-th satellite, and C be the signal carrier power. For noise power spectral density, C / This describes the ratio of signal strength to noise for a single satellite. The reflectance factor r is estimated based on the average SNR. , where parameters The parameters of the sigmoid function are set empirically. Constructing a surface reflection correction term by combining with DEM-RTK bias , ,in This is the scaling factor; S404: Structural Integrated Elevation Correction Observation , ,in and To normalize the weighting, a confidence level is set, which is determined by the accuracy of the DEM and the reliability of the reflection term. If the DEM accuracy is high, then... If the SNR is significantly low and the multipath probability is high, increase the value of the multipath probability. ; S405: Establish the state-space model of the vertical Kalman filter (V-KF): ,in Vertical height For the height deviation rate, Let V be the state vector of the vertical Kalman filter V-KF at time k. S406: Determine the observation noise covariance R: ,in For RTK observation uncertainty variance, For the variance of the DEM correction term, ,in This refers to the accuracy variance of the DEM itself. S407: Perform Kalman filtering time update, predict state and covariance. , ,in Let F be the filter state vector, F be the state transition matrix, and Q be the process noise covariance. Let K be the predicted state vector at time k. Let be the predicted state covariance matrix at time k; S408: Smoothed elevation output after filtering , ,in Let k be the vertical height obtained after filtering at time k. And record the height error estimate. ; S409: If If so, mark the DEM as abnormal and reduce it. ,like Then increase If the filtering innovation exceeds the threshold, a retest will be triggered or the system will fall back to RTK-only mode and logs will be recorded. The threshold for DEM anomaly detection; S5: When satellite obstruction or insufficient visibility occurs, IMU and visual assistance are activated and tightly coupled with GNSS for filtering and fusion calculation. The specific logical steps are as follows: S501: When GNSS availability is insufficient, IMU pre-integration and visual odometry are used as short-term alternatives and constraints. Acceleration and angular velocity data of IMU, image frames and feature points of visual sensor are collected respectively, and noise processing, time synchronization and coordinate system unification are performed. Insufficient GNSS availability is manifested as a small number of satellites, high PDOP or poor observation quality. S502: By integrating the acceleration and angular velocity of the IMU, the predicted values ​​of the current position, velocity, and attitude are obtained, forming the prior state for fusion calculation. The formula used for prediction update is: ,in Let be the posterior state at time k. The filtering state at time k-1, ( ) is an IMU prediction model. For the raw IMU measurements, The IMU sampling period; S503: Extracts current position observations from pseudorange, carrier phase, and differential positioning results generated by BeiDou and other available GNSS satellites, and forms a visual observation vector from the relative displacement or attitude change output by the visual odometry. The two together constitute the fused solution of the observation value. ; S504: Calculate the fusion residual based on the difference between the IMU predicted state and the GNSS / visual observations, and obtain the Kalman gain for posterior state update. This yields high-precision position, velocity, and attitude estimates after multi-source information fusion, among which... Let Kalman gain be at time k. For observation models, The observation vectors for fusion solution; S505: The fused positioning result is used as the current accurate state of the system, which is used to output a stable solution by the measurement operation control module. At the same time, it is fed back to the next epoch as an initial condition, so that the system can maintain continuous, stable and high-precision positioning capability even when satellite visibility is insufficient. S6: Performs measurement task control, mode switching, and alarms based on real-time quality indicators; Its real-time quality metrics include GNSS satellite quality metrics, multipath and occlusion metrics, IMU stability metrics, and fusion solution reliability metrics. All metrics can be uniformly expressed as: ,in Let k be the real-time quality index vector at time k; The specific logical steps involved in implementing this process are as follows: S601: Acquires GNSS, IMU, and visual data, and calculates the number of satellites. With position accuracy factor Carrier-to-noise ratio Multi-path intensity comprehensive index residual variance IMU bias drift and fusion covariance ; S602: Determines whether the current RTK can work normally. Its execution conditions are: If all conditions are met, the system enters high-precision RTK mode; if any condition is not met, the system enters weak field-of-view mode. in The threshold for the number of visible satellites. This is the threshold for the position accuracy factor. The carrier-to-noise ratio threshold. For multipath intensity threshold, The threshold for GNSS residual variance; S603: Meets all conditions in S602, triggering high-precision RTK mode. In high-precision RTK mode, GNSS is the primary method, with IMU / vision as a secondary method. The triggering conditions for weak field-of-view mode are: When the weak field-of-view mode is triggered, the IMU weights are automatically increased, and visual constraints are used to maintain horizontal accuracy. The triggering condition for the GNSS interruption mode is as follows: Then it completely switches to IMU + visual navigation, and GNSS enters the fusion filter with zero weight; S604: Execute alarm logic: Audio-visual prompts; The system prompts for mobile base station / equipment adjustment. Warning: Degraded positioning reliability. The system prompts for reinitialization. For RTK, a fixed solution ratio is used; S605: When the quality indicators are normal, the system performs measurement tasks, including controlling the sampling frequency, acquiring control points, and measuring the control trajectory. When the indicators decrease, the system delays sampling to maintain reliable data results. S7: Output engineering layout results, including well site coordinates, pipeline center / side stakes, station equipment layout coordinates, elevation and quality report; The output coordinate results are labeled with "planar accuracy", "elevation accuracy" and "corresponding quality indicators" (e.g., planar accuracy ±2.5cm, elevation accuracy ±3.0cm, PDOP=2.8). =0.95), to avoid the construction party misjudging the position due to the lack of precision marking; This embodiment effectively reduces reflection interference caused by dense metal facilities in well sites and stations by constructing a multi-path geometric model, dual-antenna array discrimination, and weighted suppression algorithms. By fusing local DEM terrain models and surface reflection features and introducing vertical Kalman filtering, it significantly improves the stability and accuracy of BeiDou elevation calculation. In scenarios with insufficient visibility, such as drilling rigs and large equipment obstruction, it maintains continuous and highly reliable positioning capabilities through a tightly coupled fusion strategy of IMU and visual assistance. This enables highly robust and high-precision BeiDou engineering measurement and positioning in complex oilfield environments, significantly improving the problems of traditional RTK such as easy loss of lock, easy drift, unstable elevation, and unavailability in obstructed areas, thus significantly improving engineering measurement efficiency and construction control quality.

[0018] Example 2 Based on the same inventive concept as the BeiDou-based engineering measurement and positioning method described in the foregoing embodiments, such as... Figure 6 As shown, this embodiment proposes an engineering measurement and positioning system based on BeiDou, including a receiving unit, a dual-antenna multipath detection unit, a DEM and environmental model database, an IMU / visual assistance unit, a differential enhancement unit, a data fusion and solution unit, a measurement operation control unit, and an output unit; The receiving unit is used to receive observations from BeiDou and at least one other GNSS system; The dual-antenna multipath detection unit includes hardware / software modules for SNR acquisition, phase difference calculation, multipath degree calculation, and environmental model-based reflection source matching, which are used for multipath identification and rejection based on the phase and amplitude difference of the dual antennas. The DEM and environmental model database is used to store the three-dimensional structure and topographic elevation of the oilfield and provide elevation compensation data. The DEM and environmental model database supports offline import and rapid on-site modeling. The model can be updated through laser scanning or UAV aerial surveying to adapt to construction changes. The IMU / visual aid unit includes a three-axis MEMS or fiber optic gyroscope IMU, a visual camera, and a visual inertial navigation unit. It also performs online estimation and correction of IMU bias to reduce cumulative errors and provides short-term positioning when GNSS is blocked. The differential enhancement unit is used to access CORS / RTK or PPP enhancement services, and supports VRS / network RTK, base station differential, and real-time differential data processing interface based on CMR / RTCM; The data fusion and solution unit is used to perform weighted fusion of various observations and output engineering measurement results; The measurement operation control unit is used to perform real-time quality assessment and measurement strategy switching; The output unit is used to integrate multi-source fused positioning / elevation data and real-time quality indicators, and output the final coordinate results and quality report; The BeiDou-based engineering measurement and positioning system provided in this embodiment can execute the BeiDou-based engineering measurement and positioning method provided in any embodiment of the present invention, and has the corresponding functional modules and beneficial effects of the execution method.

[0019] Example 3 Based on the BeiDou-based engineering measurement and positioning method in Embodiment 1 and the BeiDou-based engineering measurement and positioning system in Embodiment 2, the specific implementation process of this method in a complex environment of obstruction (derrick, oil storage tank) and multipath (metal structure reflection) of an offshore oilfield drilling platform includes the following steps: S01: Satellite Visibility and Multipath Candidate Set Screening: Input Platform Observation Point Coordinates: Real-time acquisition of coordinates of 6 BeiDou satellites And calculate the satellite direction vector. Combining the surface normal vectors of the platform derrick and oil storage tank Two satellites were identified. Obscured by the derrick and marked as invisible, calculate the area of ​​the metal oil storage tank from the four visible satellites. mirror point To obtain a multi-path candidate set ,in through Additional path length of reflection Corresponding time delay ; S02: Dual-antenna multipath discrimination and weighted elimination: Calculating satellite Theoretical carrier phase difference Observation phase difference The phase residual is obtained. Calculate pseudorange residuals signal-to-noise ratio decrease The fusion yields a multi-path comprehensive index. ,right Observation weights are adjusted to Remove Observations; S03: Elevation Compensation and Vertical Kalman Filtering: BeiDou RTK Output Elevation DEM interpolation elevation initial deviation And calculate the average signal-to-noise ratio. The surface reflectance factor is obtained through the sigmoid function. Reflection correction term Comprehensive elevation correction: After V-KF filtering, the output is a smoothed elevation. Elevation error ; S04: Short-term substitution and fusion of IMU / visual signals under occlusion: Multi-source fusion is triggered when the platform equipment moves into the derrick shadow area. (1) IMU predicts prior state ; (2) Visual odometry provides observation vectors ; (3) Kalman gain By fusing IMU and visual data, the localization result is output. Once the number of satellites is restored, switch back to the GNSS-dominated fusion mode.

[0020] The implementation process achieved the following results: (1) Positioning accuracy: 2.5cm in the horizontal direction and 3.2cm in the vertical direction; (2) In occluded scenarios, the positioning drift is controlled within 0.1cm; (3) After multipath suppression, the pseudorange error was reduced from 0.5m to 0.1m.

[0021] The above description is only a preferred embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any equivalent substitutions or modifications made by those skilled in the art within the scope of the technology disclosed in the present invention, based on the technical solution and inventive concept of the present invention, should be covered within the scope of protection of the present invention.

Claims

1. A BeiDou-based engineering measurement and positioning method, characterized in that, Includes the following steps: S1: Establish a three-dimensional model of the oilfield operating environment and perform satellite visibility and multipath prediction; S2: Use a GNSS receiver that supports multi-system reception to receive signals from BeiDou and at least one other satellite navigation system and record raw observations; S3: The received signal is multipath detected and weighted by a dual-antenna array and a multipath discrimination algorithm to obtain the corrected pseudorange and carrier phase observations. S4: Based on the observations output by S3, the BeiDou RTK elevation is model compensated and vertical Kalman filtered by combining the local digital elevation model (DEM) and surface reflection characteristics. S5: When satellite obstruction or insufficient visibility occurs, IMU and visual assistance are activated and tightly coupled with GNSS for filtering and fusion calculation. S6: Performs measurement task control, mode switching, and alarms based on real-time quality indicators; S7: Output engineering layout results, including well site coordinates, pipeline center / side stakes, station equipment layout coordinates, elevation and quality report.

2. The engineering measurement and positioning method based on BeiDou according to claim 1, characterized in that, In S1, the established three-dimensional model of the oilfield operation environment includes metal pipelines, tanks, drilling rig towers, station steel structures, roads and terrain elevation data, and generates satellite visibility occlusion maps and multipath reflection candidate sets at each observation point based on the three-dimensional model of the oilfield operation environment. The specific logical steps of S1 are as follows: S101: Collect spatial information of well sites, stations and pipeline areas, including topographic elevation, location and size of metal pipelines, drilling rig towers, station steel structures, three-dimensional dimensions and location information of equipment tank facilities, and location of roads, construction platforms and temporary facilities; S102: The collected data are fused into a unified coordinate system, and a 3D model of the oilfield is constructed using CAD modeling software. The 3D model of the oilfield includes a terrain elevation model, a 3D facility model, and road and construction platform models. S103: Based on the engineering layout requirements, determine the observation points for wellheads, pipeline center stakes, side stakes, and station layout points, and mark each observation point in the 3D model for subsequent satellite visibility and multipath analysis; S104: Perform ray tracing for each observation point, determine whether the satellite signal is obstructed by three-dimensional facilities or terrain, and output a satellite visibility occlusion map; S105: For each observation point, use a 3D model to calculate possible satellite signal reflecting surfaces. Satellite signal reflecting surfaces include, but are not limited to, metal pipelines, tanks and steel structures with high reflectivity, as well as ground or road reflecting surfaces. Calculate the geometric relationship between the satellite and the reflecting surface to generate a multipath reflection candidate set, which serves as the input for subsequent dual-antenna multipath discrimination or weight allocation. S106: Optimize the selection of observation points based on satellite visibility and multipath prediction results, avoid high multipath or severely occluded locations, and mark points that must be measured in severely occluded areas as having IMU / visual assistance enabled.

3. The engineering measurement and positioning method based on BeiDou according to claim 2, characterized in that, In step S104, satellite visibility analysis requires calculating the satellite-to-observation-point vector: based on the observation point... and satellite position Calculate the vector: ,in Let be the vector pointing from the observation point to the i-th satellite; When performing occlusion detection, the surface normal vectors of objects in the scene are used. and observation points Indicates the reflecting / obstructing surface, used to determine whether a satellite is obstructed, and the ray... If the ray intersects the surface of the object, the satellite will be blocked, where t represents the position of a point on the satellite signal path, a scalar that determines the direction vector along the ray originating from the observation point. How far has been traveled? Its determination function is: ； in Indicates a set of visible satellites. , and Let be the coordinates of the three vertices that constitute the surface of the j-th object.

4. The engineering measurement and positioning method based on BeiDou according to claim 2, characterized in that, The specific logical steps of S105 are as follows: S1051: Calculate the patch normal and unit normal: for a patch ( ), , ,in For dough pieces The normal vector, It is the unit normal vector; S1052: Use the method of mirrors to find possible points of reflection on a mirror surface. : Put the satellite About the dough Mirroring , For the i-th satellite Regarding the j-th facet The coordinates of the mirror point are obtained by first taking a point on the plane and denoting it as... The vector projection of the satellite onto the plane is ,in For satellite up to the j-th face The length of the vector projection; Mirror point: Next, find the ray. Intersection parameters with the plane , Intersection coordinates: Next, determine whether the intersection point lies within the face. If the intersection point lies within the face, then set... Otherwise, it is assumed that the surface does not produce a specular reflection point; in Indicates from the observation point Point to mirror point The ray equation, For the nth piece The coordinates of the reflection point that produces effective specular reflection; S1053: Verifying the geometry of specular reflection: Defining the incident direction Define the launch direction And calculate the angle difference between the two: ； like If the mirror geometry matches well, then the surface is considered to be a good mirror reflector; otherwise, it can be considered a non-primary mirror reflector. The reference threshold angle is used to determine whether it is a valid specular reflective surface; S1054: Calculate path increment: Calculate the additional path length generated by multipathing. The corresponding time delay is: Where c is the speed of electromagnetic wave propagation, and takes the value of ; Then calculate the angle of incidence. , ; S1055: Estimating Reflection Intensity and Weight , , The larger the value, the better the quality of the surface. For satellite The stronger the multipath impact, the more likely it is to be identified, weighted down, or eliminated. in Let be the reflection coefficient of the j-th surface patch. This is the path length attenuation coefficient; S1056: Satellite All that meet the threshold condition constituting a satellite Multipath reflection candidate set , This information is then used as input for subsequent dual-antenna discrimination and observation weighting / rejection.

5. The engineering measurement and positioning method based on BeiDou according to claim 1, characterized in that, In S2, the original observations include pseudorange observations. Carrier phase observation Doppler frequency shift observation Signal strength observation And channel quality metrics, among which Indicates the satellite number, Indicates frequency point; ,in For geometric distance, , For the clock difference between the receiver and the satellite, For ionospheric delay, For tropospheric delay, For pseudo-range multipath, For the purpose of observing noise, c is the propagation speed of electromagnetic waves; ,in The wavelength at that frequency point, For carrier cumulative phase, The initial ambiguity; Its actual recorded content includes Doppler frequency shift values ​​and carrier lock status; ,in Indicates that satellite i is at frequency point The received carrier power is estimated by the correlator output. represents the thermal noise power spectral density, represents the noise intensity per unit bandwidth, is the noise floor estimated by the receiver front end, and 10 is the scaling factor when converting the linear power ratio to dB. ( () is a base-10 logarithm used to convert the linear power ratio to dB-Hz units.

6. The engineering measurement and positioning method based on BeiDou according to claim 1, characterized in that, The specific logical steps of S3 are as follows: S301: Obtain pseudorange observations by simultaneously receiving satellite signals from the same epoch via the main antenna A1 and the secondary antenna A2. Carrier phase observation Doppler frequency shift observation Signal strength observation and channel quality metrics; S302: Calculate the satellite direction unit vector based on the receiver coordinates and satellite coordinates. ,in Let be the unit direction vector pointing from the observation point to the i-th satellite; Furthermore, the baseline vectors of the two antennas are obtained based on the antenna installation parameters: ,in These represent the three-dimensional coordinates of the main antenna A1 and the secondary antenna A2, respectively. S303: Calculate the ideal geometric path difference based on the satellite orientation and baseline vector, and convert it into the theoretical carrier phase difference. ,in Let be the ideal theoretical carrier phase difference corresponding to the i-th satellite; S304: Differentiate the carrier phases of two antennas at the same frequency to obtain the observed phase difference. ; S305: Compare the observed phase difference with the theoretical phase difference to obtain the phase residual reflecting the effects of multipath propagation. ; S306: Obtain pseudorange residuals from epoch filtering or differential calculation. ,in Predict pseudorange for the model; S307: Calculate the drop based on the current signal-to-noise ratio and the average signal-to-noise ratio of the sliding window. The higher this index, the higher the probability that the signal is blocked or interfered with. For the i-th satellite at frequency point The average carrier-to-noise ratio of the sliding window at that location. For the i-th satellite at frequency point The current carrier-to-noise ratio at the location; S308: Combine phase residuals, pseudorange residuals, SNR decrease, and integer jumps according to weights. ,in For weight parameters, To convert the phase residual to meters or cycles, For phase residual, For pseudo-distance residual index, This is a carrier integer cycle transition indicator; 0 indicates no transition, and 1 indicates a transition. For frequency point The wavelength; S309: A joint suppression strategy of SNR and MP is adopted to obtain the weight of each satellite observation. ,in Can be made by C / Convert to linear value As a regulating factor, ,when Hour, When it is close to 1, When it is big, Close to 0, among which For the i-th satellite, the multipath comprehensive index is used. S310: Weighted correction of satellite observation noise When severe multipathing exists, the weight is significantly reduced; if If so, the observations of that satellite are directly discarded. This represents the original observation noise variance. Let $\mathbf{i}$ be the variance of the observation noise of the $i$-th satellite after weight correction. Thresholds are removed from the observations; S311: Input the retained observations and their corresponding noise matrices into the Kalman filter to complete the RTK solution for this epoch.

7. The engineering measurement and positioning method based on BeiDou according to claim 6, characterized in that, The specific logical steps of S4 are as follows: S401: Ensure the elevation directly provided by the RTK solution. Elevation values ​​obtained by interpolation at the observation points from the DEM If the reference systems are not the same, a datum transformation is performed. If the DEM resolution is low, bilinear interpolation is performed around the observation point. ; S402: Calculate the initial deviation of DEM-RTK ; S403: Calculate the average signal quality of the satellite signal at the observation point. Where N is the number of satellites participating in the averaging. For the average carrier-to-noise ratio, C / Let C be the carrier-to-noise ratio of the i-th satellite, and C be the signal carrier power. For noise power spectral density, C / This describes the ratio of signal strength to noise for a single satellite. The reflectance factor r is estimated based on the average SNR. , where parameters The parameters of the sigmoid function are set empirically. Constructing a surface reflection correction term by combining with DEM-RTK bias , ,in This is the scaling factor; S404: Structural Integrated Elevation Correction Observation , ,in and To normalize the weighting, a confidence level is set, which is determined by the accuracy of the DEM and the reliability of the reflection term. If the DEM accuracy is high, then... If the SNR is significantly low and the multipath probability is high, increase the value of the multipath probability. ; S405: Establish the state-space model of the vertical Kalman filter (V-KF): ,in Vertical height For the height deviation rate, Let V be the state vector of the vertical Kalman filter V-KF at time k. S406: Determine the observation noise covariance R: ,in For RTK observation uncertainty variance, The variance of the uncertainty of the DEM correction term. ,in This refers to the accuracy variance of the DEM itself. S407: Perform Kalman filtering time update, predict state and covariance. , ,in Let F be the filter state vector, F be the state transition matrix, and Q be the process noise covariance. Let K be the predicted state vector at time k. Let be the predicted state covariance matrix at time k; S408: Smoothed elevation output after filtering , ,in Let k be the vertical height obtained after filtering at time k. And record the height error estimate. ; S409: If If so, mark the DEM as abnormal and reduce it. ,like Then increase If the filtering innovation exceeds the threshold, a retest will be triggered or the system will fall back to RTK-only mode and logs will be recorded. The threshold for DEM anomaly detection.

8. The engineering measurement and positioning method based on BeiDou according to claim 1, characterized in that, The specific logical steps of S5 are as follows: S501: When GNSS availability is insufficient, IMU pre-integration and visual odometry are used as short-term alternatives and constraints. Acceleration and angular velocity data of IMU, image frames and feature points of visual sensor are collected respectively, and noise processing, time synchronization and coordinate system unification are performed. Insufficient GNSS availability is manifested as a small number of satellites, high PDOP or poor observation quality. S502: By integrating the acceleration and angular velocity of the IMU, the predicted values ​​of the current position, velocity, and attitude are obtained, forming the prior state for fusion calculation. The formula used for prediction update is: ,in Let be the posterior state at time k. Let k-1 be the filtering state. ( ) is an IMU prediction model. For the raw measurements of the IMU, The IMU sampling period; S503: Extracts current position observations from pseudorange, carrier phase, and differential positioning results generated by BeiDou and other available GNSS satellites, and forms a visual observation vector from the relative displacement or attitude change output by the visual odometry. The two together constitute the fused solution of the observation value. ; S504: Calculate the fusion residual based on the difference between the IMU predicted state and the GNSS / visual observations, and obtain the Kalman gain for posterior state update. This yields high-precision position, velocity, and attitude estimates after multi-source information fusion, among which... Let Kalman gain be at time k. For observation models, The observation vectors for fusion solution; S505: The fused positioning result is used as the current accurate state of the system. It is used to output a stable solution from the measurement operation control module and fed back to the next epoch as an initial condition, so that the system can maintain continuous, stable and high-precision positioning capability even when satellite visibility is insufficient.

9. The engineering measurement and positioning method based on BeiDou according to claim 1, characterized in that, In S6, the real-time quality indicators include GNSS satellite quality indicators, multipath and occlusion indicators, IMU stability indicators, and fusion solution reliability indicators. All indicators can be uniformly expressed as: ,in Let k be the real-time quality index vector at time k; The specific logical steps are as follows: S601: Acquires GNSS, IMU, and visual data, and calculates the number of satellites. With position accuracy factor Carrier-to-noise ratio Multi-path intensity comprehensive index residual variance IMU bias drift and fusion covariance ; S602: Determines whether the current RTK can work normally. Its execution conditions are: If all conditions are met, the system enters high-precision RTK mode; if any condition is not met, the system enters weak field-of-view mode. in The threshold for the number of visible satellites. This is the threshold for the position accuracy factor. The carrier-to-noise ratio threshold. For multipath intensity threshold, The threshold for GNSS residual variance; S603: Meets all conditions in S602, triggering high-precision RTK mode. In high-precision RTK mode, GNSS is the primary method, with IMU / vision as a secondary method. The triggering conditions for weak field-of-view mode are: When the weak field-of-view mode is triggered, the IMU weights are automatically increased, and visual constraints are used to maintain horizontal accuracy. The triggering condition for the GNSS interruption mode is as follows: Then it completely switches to IMU + visual navigation, and GNSS enters the fusion filter with zero weight; S604: Execute alarm logic: Audio-visual prompts; The system prompts for mobile base station / equipment adjustment. Warning: Degraded positioning reliability. The system prompts for reinitialization, where... For RTK, a fixed solution ratio is used; S605: When the quality indicators are normal, the system performs measurement tasks, including controlling the sampling frequency, acquiring control points, and measuring the control trajectory. When the indicators decrease, the system delays sampling to maintain reliable data results.

10. A BeiDou-based engineering measurement and positioning system, used to implement the method described in any one of claims 1-9, characterized in that, It includes a receiving unit, a dual-antenna multipath detection unit, a DEM and environmental model database, an IMU / visual assistance unit, a differential enhancement unit, a data fusion and solution unit, a measurement operation control unit, and an output unit; The receiving unit is used to receive observations from BeiDou and at least one other GNSS system. The dual-antenna multipath detection unit includes hardware / software modules for SNR acquisition, phase difference calculation, multipath degree calculation, and reflection source matching based on the environment model, and is used for multipath identification and elimination based on the phase and amplitude difference of the dual antennas. The DEM and environmental model database is used to store the three-dimensional structure and topographic elevation of the oilfield and provide elevation compensation data. The DEM and environmental model database supports offline import and rapid on-site modeling. The model can be updated by laser scanning or UAV aerial survey to adapt to construction changes. The IMU / visual aid unit includes a three-axis MEMS or fiber optic gyroscope IMU, a visual camera, and a visual inertial navigation unit. It also performs online estimation and correction of the IMU bias to reduce cumulative errors and provides short-term positioning when GNSS is blocked. The differential enhancement unit is used to access CORS / RTK or PPP enhancement services, and supports VRS / network RTK, base station differential, and real-time differential data processing interface based on CMR / RTCM. The data fusion and solution unit is used to perform weighted fusion of various observations and output engineering measurement results; The measurement operation control unit is used to perform real-time quality assessment and measurement strategy switching; The output unit is used to integrate the multi-source fusion positioning / elevation data and real-time quality indicators, and output the final coordinate results and quality report.

Images