A method for terrain data acquisition based on GNSS and total station

Abstract

This invention relates to the field of engineering surveying and positioning technology, specifically a terrain data acquisition method based on the combined use of GNSS and a total station. The method includes: acquiring satellite observation data streams and optical polar coordinate observation data streams of target points; calculating a quality factor matrix based on continuous epoch satellite observation data streams to automatically switch acquisition modes; combining a pre-calibrated spatial eccentricity vector with GNSS dynamically calculated coordinates, total station optical polar coordinates, and inertial measurement unit attitude data for synchronous latching and rigid body geometric translation calculations to obtain underlying homogeneous coordinates; using these coordinates as the starting point, calculating the theoretical horizontal distance and theoretical positive azimuth angle backwards to a physically invariant reference point, and comparing this with synchronously measured values ​​from the total station to determine and eliminate pseudo-fixed solutions affected by multipath effects or non-line-of-sight signals; and simultaneously constructing a dynamic health assessment model to update the temporary reference point pool and eliminate failed points.

Description

Technical Field

[0001] This invention relates to the field of engineering surveying and positioning technology, specifically to a method for acquiring terrain data based on the combination of GNSS and total station. Background Technology

[0002] In the fields of topographic surveying and construction layout, GNSS positioning and optical measurement equipment are usually used in combination to obtain topographic data. However, in complex construction sites (such as scenarios where tall buildings or structures block the view, metal structures reflect the data, or temporary facilities are densely distributed), GNSS signals are affected by multipath effects and interference from non-line-of-sight signals, which can lead to fixed solution artifacts in the solution results. That is, although the receiver outputs a fixed solution status indicator, the actual coordinates have a large deviation. Such gross errors are often difficult to identify through conventional quality indicators such as the receiver's signal-to-noise ratio or accuracy attenuation factor, and are therefore mistakenly used as high-precision valid data. Furthermore, existing joint measurement methods suffer from deficiencies in time synchronization and spatial benchmark unification of coordinate data acquired by different types of sensors. Systematic biases remaining after multi-source data fusion are difficult to detect and correct. Simultaneously, benchmarks relied upon for field measurements often drift or even fail over time due to factors such as ground settlement, mechanical vibration, and human disturbance. Existing methods generally lack means for continuous and quantitative assessment of benchmark reliability, resulting in failed benchmarks not being identified and eliminated in a timely manner and still being used as the starting point for subsequent measurements. This leads to errors continuously propagating and accumulating in the measurement chain, ultimately affecting the accuracy and reliability of the topographic data acquisition results for the entire survey area.

[0003] To address this, a topographic data acquisition method based on the combination of GNSS and total station is proposed. Summary of the Invention

[0004] This invention aims to provide a terrain data acquisition method based on the combination of GNSS and total station. By reverse-calculating the coordinates calculated by GNSS to the physical invariant reference point and comparing and cross-validating them with the measured data of the total station, the method can achieve real-time identification and elimination of pseudo-fixed solutions.

[0005] To achieve the above objectives, the present invention provides the following technical solution: A method for topographic data acquisition based on GNSS and total station integration includes: Acquire satellite observation data stream and optical polar coordinate observation data stream of the target point; calculate GNSS quality factor matrix based on satellite observation data stream of continuous epochs; and automatically switch acquisition modes based on the quality factor matrix. Using a pre-imported site model, a dynamic accuracy threshold matching the current site is extracted based on the current spatial location. Within the same period, the dynamic solution coordinates of GNSS, the optical polar coordinates of the total station, and the attitude data output by the inertial measurement unit are synchronously latched. Rigid body geometric translation calculations are performed by combining the pre-calibrated antenna and the spatial eccentric vector from the omnidirectional prism to the bottom tip, and the underlying homogeneous coordinates are obtained respectively. Using the GNSS coordinates in the underlying homogeneous coordinates as the starting point, the theoretical horizontal distance and theoretical positive azimuth from the starting point to the pre-determined physical invariant reference point are calculated in reverse. The theoretical positive azimuth and theoretical horizontal distance are compared with the azimuth and horizontal distance measured synchronously by the total station to determine whether the GNSS coordinates are pseudo-fixed solutions. A dynamic health assessment model is constructed to calculate the health score of each temporary reference point in the test area. The temporary reference point pool is updated based on the health score, and failed points are removed.

[0006] Preferably, the step of automatically switching the acquisition mode includes: parsing the satellite observation data stream of the continuous epochs, extracting the carrier signal-to-noise ratio, spatial position accuracy factor, and ionospheric delay rate of change of dual-frequency carrier phase observations of each visible satellite; detecting cycle slip markers in the satellite observation data stream using a geometric distance-free combined observation model, and quantifying multipath effect error by combining the ionospheric delay rate of change and the carrier signal-to-noise ratio; constructing a GNSS quality factor matrix using the spatial position accuracy factor, the cycle slip marker, and the multipath effect error as dynamic parameters, and calculating the comprehensive signal health index of the current epoch; and automatically switching the acquisition mode according to the comprehensive signal health index.

[0007] Preferably, the step of extracting the dynamic accuracy threshold matching the current site includes: obtaining the relative local design coordinates in the BIM design file; transforming the relative local design coordinates to the engineering construction coordinate system based on the known geoid elevation anomaly parameters and projection deformation modification parameters of the survey area to obtain a three-dimensional survey area model; dividing the topological boundary of the three-dimensional survey area model into non-overlapping three-dimensional polygon bounding boxes, and establishing a spatial index using an octree structure; using the tolerance values ​​of different accuracy requirement areas of the survey area as attribute fields and binding them to the topological nodes of the corresponding bounding boxes; inputting the current spatial coordinates as probe points into the spatial index, using a point-polygon inclusion test algorithm to lock the target three-dimensional polygon bounding box to which the probe point belongs, and parsing the attribute fields of the bounding box to obtain the dynamic accuracy threshold.

[0008] Preferably, the step of obtaining the underlying homogeneous coordinates includes: using a second pulse signal as a time reference, synchronously acquiring the dynamic calculated coordinates of the GNSS, the optical polar coordinates of the total station, and the attitude data output by the inertial measurement unit, wherein the attitude data includes roll angle, pitch angle, and yaw angle; constructing rotation matrices around the three coordinate axes of the spatial rectangular coordinate system based on the roll angle, pitch angle, and yaw angle, and multiplying the three rotation matrices to obtain the rigid body pose isomorphic mapping matrix; using the rigid body pose isomorphic mapping matrix to perform rotation transformation on the pre-calibrated spatial eccentricity vector from the antenna and the omnidirectional prism to the bottom tip, and calculating the corresponding spatial eccentricity compensation amount; subtracting the spatial eccentricity compensation amount from the GNSS dynamic calculated coordinates and the total station polar coordinates converted to spatial rectangular coordinates to form the underlying homogeneous coordinates.

[0009] Preferably, the reverse calculation of the theoretical horizontal distance and theoretical positive azimuth from the starting point to the pre-determined physical invariant reference point includes: projecting the GNSS coordinates in the underlying source coordinates to the engineering construction plane coordinate system to obtain the two-dimensional plane starting coordinates; extracting the known reference plane coordinates of the pre-determined physical invariant reference point in the engineering construction plane coordinate system; calculating the geometric distance between the two-dimensional plane starting coordinates and each of the known reference plane coordinates as the theoretical horizontal distance; calculating the coordinate azimuth based on the coordinate increment of the two-dimensional plane starting coordinates and each of the known reference plane coordinates, and correcting it by combining the meridian convergence angle and direction modification parameters of the current position of the survey area to obtain the theoretical positive azimuth.

[0010] Preferably, the step of determining whether the GNSS coordinates are pseudo-fixed solutions includes: obtaining the measured horizontal distance and measured azimuth angle of the total station synchronously measured to each of the physical invariant reference points at the same timestamp, comparing them with the corresponding theoretical horizontal distance and theoretical positive azimuth angle, and calculating the spatial position residuals of each reference point; when the maximum value of each of the spatial position residuals is greater than the dynamic accuracy threshold, the GNSS coordinates are determined to be pseudo-fixed solutions and discarded; when the optical mapping data of the total station is valid, the coordinates are output according to the optical polar coordinates of the total station; when the total station loses synchronization, the attitude data output by the inertial measurement unit is fused with the historical epochs to perform short-time coordinate extrapolation.

[0011] Preferably, the step of constructing a dynamic health assessment model and calculating the health score of each temporary benchmark point within the survey area includes: obtaining a material disturbance resistance factor characterizing the physical displacement resistance of the benchmark point based on the geological and engineering structural properties of the carrier to which each temporary benchmark point is attached; obtaining the measured coordinates of each temporary benchmark point periodically re-measured by a total station, and comparing them with the initial coordinates to obtain spatial mapping residuals; applying coordinate residual penalties to the spatial mapping residuals exceeding the allowable mapping tolerance based on the robustness estimation theory, and obtaining a mapping reliability weight negatively correlated with the spatial mapping residuals; constructing an exponential time decay factor characterizing the increase in the probability of natural settlement failure of each temporary benchmark point over time based on the mechanical vibration and soil rheological characteristics of the construction site; multiplying the material disturbance resistance factor, the mapping reliability weight, and the exponential time decay factor to obtain the health score of each temporary benchmark point; updating the temporary benchmark point pool based on the health score and removing failed points.

[0012] Compared with the prior art, the beneficial effects of the present invention are as follows: 1. This invention calculates the GNSS quality factor matrix from satellite observation data streams of consecutive epochs and automatically switches the acquisition mode accordingly. It can dynamically adjust the data acquisition strategy according to the real-time changes in satellite signal quality, avoiding the problem of fixed acquisition modes in traditional methods that cannot adapt to fluctuations in on-site signal quality. This allows for timely response in the early stages of signal quality degradation, providing a more reliable original observation basis for subsequent data processing and improving the adaptability and robustness of the overall acquisition process.

[0013] 2. This invention uses GNSS coordinates as the starting point to calculate the theoretical horizontal distance and theoretical positive azimuth from the physical invariant reference point, and compares them with the azimuth and horizontal distance measured synchronously by the total station. This can directly and quantitatively identify whether there are pseudo-fixed solutions in GNSS coordinates due to multipath effects or non-line-of-sight signals. Compared with traditional discrimination methods that rely on receiver status indicators or accuracy attenuation factors, this method has higher discrimination sensitivity and reliability, and can effectively avoid situations where gross errors are mistakenly output as valid data.

[0014] 3. This invention integrates the synchronously latched multi-source coordinate rigid body geometric translation calculation, the pseudo-fixed solution determination mechanism based on physically invariant reference points, and the temporary reference point dynamic health assessment model. The precise unification of the underlying homogeneous coordinates provides an accurate comparison basis for pseudo-fixed solution determination. The pseudo-fixed solution determination results provide a reliable basis for the health assessment model to eliminate failure points. In turn, the updated reference point pool after health assessment improves the accuracy of subsequent pseudo-fixed solution determination. The three support each other and iterate in a closed loop, significantly improving the overall accuracy, reliability, and field adaptability of the entire terrain data acquisition process, and avoiding the error accumulation problem caused by the isolated operation of each link. Attached Figure Description

[0015] Figure 1 This is a flowchart illustrating the steps of a terrain data acquisition method based on the combination of GNSS and total station according to the present invention. Figure 2 This is a flowchart illustrating the adaptive switching of multi-source acquisition modes according to the present invention. Figure 3 This is a schematic diagram illustrating the process of dynamic health assessment and updating of the baseline point in this invention. Detailed Implementation

[0016] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0017] Please see Figures 1 to 3 This invention provides a topographic data acquisition method based on the combination of GNSS and total station, referring to... Figure 1 Step-by-step flowchart Figure 2 A flowchart illustrating the adaptive switching process of multi-source acquisition modes, and Figure 3 A flowchart illustrating the benchmark dynamic health assessment and update process; the technical solution is as follows: Acquire satellite observation data stream and optical polar coordinate observation data stream of the target point; calculate GNSS quality factor matrix based on satellite observation data stream of continuous epochs; and automatically switch acquisition modes based on the quality factor matrix. Using a pre-imported site model, a dynamic accuracy threshold matching the current site is extracted based on the current spatial location. Within the same period, the dynamic solution coordinates of GNSS, the optical polar coordinates of the total station, and the attitude data output by the inertial measurement unit are synchronously latched. Rigid body geometric translation calculations are performed by combining the pre-calibrated antenna and the spatial eccentric vector from the omnidirectional prism to the bottom tip, and the underlying homogeneous coordinates are obtained respectively. Using the GNSS coordinates in the underlying homogeneous coordinates as the starting point, the theoretical horizontal distance and theoretical positive azimuth from the starting point to the pre-determined physical invariant reference point are calculated in reverse. The theoretical positive azimuth and theoretical horizontal distance are compared with the azimuth and horizontal distance measured synchronously by the total station to determine whether the GNSS coordinates are pseudo-fixed solutions. A dynamic health assessment model is constructed to calculate the health score of each temporary reference point in the test area. The temporary reference point pool is updated based on the health score, and failed points are removed.

[0018] Example 1: This invention can be applied to scenarios such as topographic surveying, deformation monitoring, and construction layout at engineering construction sites. The various links work together in a closed-loop iteration to jointly ensure the reliability and consistency of topographic data collection results under complex construction site conditions.

[0019] First, the automatic switching acquisition mode steps include: parsing the satellite observation data stream of the continuous epochs, extracting the carrier signal-to-noise ratio, spatial position accuracy factor, and ionospheric delay rate of change of dual-frequency carrier phase observations for each visible satellite; detecting cycle slip markers in the satellite observation data stream using a geometric distance-free combined observation model, and quantifying multipath effect error by combining the ionospheric delay rate of change and the carrier signal-to-noise ratio; constructing a GNSS quality factor matrix using the spatial position accuracy factor, the cycle slip marker, and the multipath effect error as dynamic parameters, and calculating the comprehensive signal health index for the current epoch; and automatically switching the acquisition mode based on the comprehensive signal health index.

[0020] Specifically, the GNSS signal processing module built into the acquisition terminal continuously parses the raw satellite observation data stream at fixed epoch intervals (e.g., 1 second). For each complete epoch, it extracts three types of physical observations in parallel from the raw observation data frame as the raw input for subsequent calculations: First, the carrier signal-to-noise ratio of each visible satellite at the current epoch. This value, in decibels and hertz, is stored in the corresponding field of the receiver's raw observation file, reflecting the ratio of satellite signal strength to background noise. Second, the spatial position accuracy factor calculated from the spatial geometric distribution of all visible satellites at the current epoch. This value is a dimensionless scalar, output in real time by the receiver's navigation calculation module and written to the data buffer. Third, the ionospheric delay rate of change obtained by subtracting the dual-frequency carrier phase observations of adjacent consecutive epochs. This rate of change, in centimeters per epoch, reflects the degree of influence of ionospheric activity on the propagation paths of different frequency carrier phases. Its absolute value is usually less than 0.2 centimeters per epoch during periods of ionospheric calm, but can rise to more than 1 centimeter per epoch when solar activity is active or the observation environment is disturbed.

[0021] After the three types of raw observations are extracted, the signal processing module calls the geometric distance-free combined observation model to perform cycle slip detection on the carrier phase observation sequence. Geometric distance-free combined observation eliminates the geometric distance term from the satellite to the receiver by linearly combining the dual-frequency carrier phase observations, ensuring that the time series of the combined values ​​contains only the slowly varying ionospheric delay residual when cycle slip-free. When a sudden change occurs in the combined value at a certain epoch and the change exceeds a preset detection threshold (e.g., 0.1 wavelengths), the cycle slip flag of the corresponding satellite at that epoch is set to valid. The processing module records this flag in binary form in the state array of each satellite. During cycle slip detection... Upon completion of the measurement, the processing module further performs joint analysis of the ionospheric delay change rate of that epoch and the corresponding satellite's carrier signal-to-noise ratio (SNR). Based on the satellite elevation angle, it calls the segmented carrier SNR threshold. When the carrier SNR variance or the duration of decline of a satellite exceeds the corresponding elevation angle threshold for N consecutive epochs (in this embodiment, 5 consecutive epochs are selected, but this can be adjusted according to the actual situation), it is determined that the satellite is currently significantly affected by the multipath effect by combining the geometrically insensitive distance combination abrupt change and the receiver's calculated covariance. The multipath effect error estimate is then quantified and calculated. This estimate is stored in the form of pseudorange residuals in meters.

[0022] After obtaining the spatial position accuracy factor, the cycle slip flag array for each satellite, and the multipath effect error estimate, the signal processing module assembles the three types of parameters into a GNSS quality factor matrix according to a preset weight structure. The row dimension of this matrix corresponds to the sequence number of each visible satellite in the current epoch, and the column dimensions are, in order, the cycle slip flag, the multipath effect error estimate, and its contribution to the overall spatial position accuracy factor for that satellite. All elements in the matrix undergo dimension normalization, so that the numerical range of different physical quantities is uniformly mapped to the interval between 0 and 1. Subsequently, the processing module assigns a set of initial weights to each column of the quality factor matrix according to preset weighting coefficients (where the spatial position accuracy factor has a weight of 0.4, the cycle slip comprehensive influence has a weight of 0.35, and the multipath error has a weight of 0.25). Specifically, it first assigns a set of initial weights to each column of the quality factor matrix. Initially, all observations are assigned equal weights. Then, adjustment calculations are performed using a large amount of measured data to calculate the posterior unit weight variance for each observation. If the posterior variance of a certain type of observation is too large, it indicates that its initial weight allocation is unreasonable. The system will readjust the weights using the variance component calculation formula and proceed to the next iteration. The iteration stops when the posterior unit weight variances of all types of observations are equal. The ratio of the inverses of the variances obtained at this time (which is the weight set above) is then weighted and summed to obtain a comprehensive signal health index reflecting the overall satellite signal availability at the current epoch. This comprehensive index is a dimensionless real number between 0 and 1. The higher the value, the better the current satellite observation environment. The calculation results are written to the processor's status register in real time and read by the mode switching decision module.

[0023] Furthermore, the system constructs a state transition mechanism based on the joint perception of spatial visibility and measurement environment errors, enabling automatic switching between three modes: satellite navigation and positioning system-dominated mode, factor graph optimization joint mode, and total station polar coordinate-dominated mode. The state transition triggering conditions simultaneously incorporate two types of inputs: dynamic accuracy threshold and on-site visibility detection results. This ensures that mode switching decisions no longer solely rely on satellite signal quality indicators, but rather incorporate site accuracy requirements and optical visibility conditions into the joint decision-making process, thereby dynamically matching the acquisition mode with the actual measurement environment of the current work area.

[0024] Specifically, the state transition decision module reads two input data simultaneously at each acquisition epoch: one is the current epoch dynamic precision threshold output by the dynamic precision threshold module, which is input in the form of a single-precision floating-point number in millimeters, reflecting the tolerance requirements of the spatial area where the probe point is currently located; the other is the physical line-of-sight status flag between the total station's optical line of sight and the omnidirectional prism output by the line-of-sight detection module. This flag is stored in binary form, indicating that line-of-sight is valid and obstruction is invalid. The line-of-sight detection module obtains this flag by continuously monitoring the lock status register of the total station's automatic tracking servo. The lock status register is refreshed by the total station's measurement controller at a fixed frequency. Both input data are read and sent to the joint condition judgment unit of the state transition decision module within the same processing cycle.

[0025] When the dynamic accuracy threshold read by the joint condition judgment unit is less than or equal to 5 mm, and the visibility marker is also in a valid state, the decision module determines that the current working point is in a high-precision area with strict coordinate accuracy requirements, such as the core tube shear wall or embedded part positioning, and the optical link between the total station and the prism is unobstructed. Under this condition, even if the satellite signal is in a fixed solution state, its residual error may exceed the tolerance range of this area. Therefore, the decision module issues a weight reduction command to the GNSS solution weight control unit to reduce the weight of the GNSS observation data in the subsequent fusion solution to a minimum value close to zero. The system is forced to enter the total station polar coordinate dominant mode. The subsequent coordinate solution depends entirely on the optical angle measurement and laser ranging data of the total station. After polar coordinate reduction and eccentricity compensation, the underlying homogeneous coordinates are directly output. The GNSS solution channel continues to run, but its output does not participate in the generation of coordinate results for the current epoch. The total station measurement controller synchronously switches the servo tracking mode to active precision tracking to maintain stable locking of the prism.

[0026] When the dynamic accuracy threshold read by the joint condition judgment unit is within the intermediate tolerance range of approximately 10 mm, and the comprehensive evaluation result of the quality factor matrix shows that the current number of visible satellites is at an edge level or the carrier-to-noise ratio fluctuates drastically (e.g., the standard deviation of the carrier-to-noise ratio exceeds the measurement threshold within 3 consecutive epochs), and the visibility markers intermittently fail within the current retest window (e.g., the cumulative number of failures exceeds 3 within the last 10 epochs), the decision module determines that the current operating environment is a complex and difficult scenario, such as a deep foundation pit or high-density steel structure, which intermittently interferes with both GNSS signals and optical visibility. In such cases, a single sensor-dominated mode cannot stably support coordinate output at this accuracy level. The decision module then sends an activation command to the fusion calculation control unit, and the system switches to factor-based mode. The graph optimization joint mode integrates the satellite observation data and total station optical polar coordinate data synchronously latched at the current epoch into independent physical observation nodes. Variable nodes represent the three-dimensional coordinate state variables of the centering rod base and tip to be determined, and factor nodes represent the constraint relationships between the sensor observations and the state variables. The factor graph topology structure of the current epoch is constructed, and a nonlinear least squares solver is called to perform joint adjustment on all observation nodes. During the iterative solution process, the residuals of each sensor observation are synchronously incorporated into the topology inversion, so that the final output underlying homogeneous coordinates can still achieve better solution robustness than the single sensor-dominated mode even under the condition that there is intermittent quality degradation in both types of sensor data. The solution results are written into the coordinate output control module in the form of three-dimensional coordinate triples.

[0027] When the dynamic accuracy threshold read by the joint condition judgment unit is greater than or equal to 30 mm, and the comprehensive index of the quality factor matrix shows that the current satellite signal is in a good state under the open sky field of view (e.g., the spatial position accuracy factor is less than 2 and there is no cycle slip flag trigger), the decision module determines that the current operation point is in an area with relatively relaxed coordinate accuracy requirements, such as flat site or earthwork backfilling, and the satellite signal quality is sufficient to independently support the RTK fixed solution output of this accuracy level. The decision module sends a switching command to the acquisition controller, and the system enters the satellite navigation and positioning system dominant mode. The GNSS solution channel outputs the dynamically calculated coordinates with standard weights. At the same time, the decision module sends a power management command to the total station servo control unit. The total station servo motor switches from active precision tracking to passive following standby mode. The servo drive current drops to the minimum level required to maintain the basic position to reduce unnecessary energy consumption. The total station measurement channel remains in a ready state and can quickly resume active tracking when the decision results of subsequent epochs change.

[0028] By introducing dynamic accuracy threshold and optical visibility status as joint triggering conditions into the mode switching decision, the system can select the acquisition mode that best matches the current work area based on the combination of the two conditions in complex construction environments where site accuracy requirements and sensor availability change simultaneously. This reduces the risk of mode selection mismatch due to differences in site accuracy requirements when relying solely on satellite signal quality indicators for mode switching.

[0029] By performing multidimensional quantitative analysis and three-level threshold judgment on continuous epoch satellite observation data, the system can adjust the data acquisition strategy in a timely manner when the satellite signal quality changes. It can trigger mode switching response in the early stage of signal quality degradation, effectively reducing the probability of pseudo-fixed solutions entering the subsequent processing flow and improving the overall acquisition process's adaptability to complex observation environments.

[0030] Further, the step of extracting the dynamic accuracy threshold matching the current site includes: obtaining the relative local design coordinates in the BIM design file; transforming the relative local design coordinates to the engineering construction coordinate system based on the known geoid elevation anomaly parameters and projection deformation modification parameters of the survey area to obtain a three-dimensional survey area model; dividing the topological boundary of the three-dimensional survey area model into non-overlapping three-dimensional polygon bounding boxes, and establishing a spatial index using an octree structure; using the tolerance values ​​of different accuracy requirement areas of the survey area as attribute fields and binding them to the topological nodes of the corresponding bounding boxes; inputting the current spatial coordinates as probe points into the spatial index, using a point-polygon inclusion test algorithm to lock the target three-dimensional polygon bounding box to which the probe point belongs, and parsing the attribute fields of the bounding box to obtain the dynamic accuracy threshold.

[0031] Specifically, before the data acquisition operation officially begins, the data preprocessing module reads the BIM design file of the current survey area from the engineering digital management platform. This file is stored in the standard IFC format and contains the three-dimensional geometric vertex data of each component of the building structure in the local coordinate system of the project. The coordinate unit is millimeters, and the coordinate origin is fixed at the project design reference point. The vertex coordinates of all components are organized in the form of a triplet list. Each record in the list corresponds to the local design coordinate value of a component vertex. This list is loaded into the processor's memory buffer as the raw input data for this step.

[0032] After reading the local design coordinate list, the preprocessing module performs two-stage coordinate transformations to unify them to the engineering construction coordinate system. The first-stage transformation introduces the known geoid elevation anomaly parameters of the survey area. These parameters are obtained by fitting the difference between the measured orthometric height and the corresponding ellipsoidal height of no less than five known leveling points within the survey area, and are stored in the configuration file in meters. This is used to correct the ellipsoidal height component of each vertex to the corresponding normal height, eliminating the elevation system deviation caused by inconsistencies in the reference surface. The second-stage transformation introduces the projection deformation modification parameters of the survey area. These parameters are pre-calculated based on the longitude difference between the central meridian of the survey area and the standard meridian of the projection zone, as well as the average altitude of the survey area, and are fixed in the configuration file. This is used to apply length normalization correction to the horizontal coordinate components, correcting the plane coordinates under the geodetic coordinate system to the practical plane coordinates under the engineering construction coordinate system, ensuring the consistency between the coordinate results and the actual distances measured on site. After the two-stage transformation is completed, the vertex coordinates of all components are rewritten into memory in the form of three-dimensional coordinate triplets under the engineering construction coordinate system, with the unit unified in millimeters, forming the input dataset for subsequent modeling steps, i.e., the vertex coordinate set of the three-dimensional survey area model.

[0033] After the set of vertex coordinates of the 3D survey area model is generated, the spatial index construction module reads the set, extracts closed polygonal faces one by one according to the topological boundaries of each component, and aligns the calculation axes of the entire set of polygonal faces with the outer rectangle to determine the 3D spatial range of the entire survey area in the engineering construction coordinate system. Then, the octree recursive splitting process is started with this spatial range as the root node. The splitting rule is: when the number of polygonal faces contained in a node exceeds the preset capacity limit (e.g., 8 faces), the node is divided into 8 child nodes along the three coordinate axes, and each polygonal face is assigned to the corresponding child node according to its spatial position. The recursion continues until the number of faces in all leaf nodes does not exceed the capacity limit or the node side length is reduced to the preset minimum resolution (e.g., 0.1 meters). Finally, an octree index structure covering the 3D space of the engineering construction coordinate system is formed. This structure is persistently stored in the local storage module as a tree data object for subsequent query and retrieval.

[0034] After the octree index is constructed, the attribute binding module reads the survey area accuracy partitioning configuration table from the engineering digital management platform. This configuration table is stored in structured text format, and each record contains three fields: area name identifier, corresponding geometric range description, and tolerance value. The tolerance value is in millimeters. For example, the tolerance value for the core tube shear wall area is 3 millimeters, the tolerance value for the main structure beam and column area is 5 millimeters, and the tolerance value for the site leveling area is 30 millimeters. Based on the geometric range description of each record, the attribute binding module locates all leaf nodes that intersect with the range in the octree and writes the corresponding tolerance value into the attribute field of each leaf node in the form of key-value pairs, thus completing the binding of the tolerance value with the spatial area. The binding result is synchronously updated to the octree persistent file in the local storage module.

[0035] During the data acquisition process, the current spatial coordinates are output in real time by the GNSS calculation module or the total station coordinate reduction module and transmitted to the spatial query module in the form of a three-dimensional coordinate triplet under the engineering construction coordinate system. This coordinate triplet is the input data of the probe point. The spatial query module first performs an inclusion judgment on the spatial range of the probe point and the root node of the octree, and then recursively traverses along the direction from the tree to the leaf node. At each level, it only enters the child node that intersects with the spatial position of the probe point until it reaches the leaf node. Then, it performs a closed body classification judgment within the leaf node: it prioritizes the parity judgment of the ray-triangle intersection point for candidate three-dimensional closed bodies; if the construction area is defined by floor or work surface, it first performs preliminary screening by elevation interval, and then projects it onto the two-dimensional plane to perform the judgment of the two-dimensional point within the polygon. After the test is passed, the leaf node to which the target polygon belongs is locked, and the tolerance value stored in the attribute field of the leaf node is read. This value is used as the dynamic precision threshold of the current epoch and output to the subsequent pseudo-fixed solution determination module. The output format is a single-precision floating-point number in millimeters. This output value is automatically updated in each acquisition epoch as the spatial position of the probe point changes.

[0036] By precisely aligning the BIM design coordinates with the engineering construction coordinate system and constructing an octree spatial index, the system can automatically match the tolerance requirements of the corresponding area based on the real-time spatial location of the acquisition device. This avoids the lag and misconfiguration risks caused by manually configuring the accuracy threshold for partitioning, and ensures that the accuracy standard for pseudo-fixed solution determination remains dynamically consistent with the actual construction requirements on site.

[0037] Further, the step of obtaining the underlying homogeneous coordinates includes: using a second pulse signal as a time reference, synchronously acquiring the dynamic calculated coordinates of the GNSS, the optical polar coordinates of the total station, and the attitude data output by the inertial measurement unit, wherein the attitude data includes roll angle, pitch angle, and yaw angle; constructing rotation matrices around the three coordinate axes of the spatial rectangular coordinate system based on the roll angle, pitch angle, and yaw angle, and multiplying the three rotation matrices to obtain the rigid body pose isomorphic mapping matrix; using the rigid body pose isomorphic mapping matrix to perform rotation transformation on the pre-calibrated spatial eccentricity vector from the antenna and the omnidirectional prism to the bottom tip, and calculating the corresponding spatial eccentricity compensation amount; subtracting the spatial eccentricity compensation amount from the GNSS dynamic calculated coordinates and the total station polar coordinates converted to spatial rectangular coordinates respectively to form the underlying homogeneous coordinates.

[0038] Specifically, the integrated centering rod physically integrates three types of sensors: a GNSS antenna, an omnidirectional prism, and an inertial measurement unit. These three sensors are installed at different heights on the rod, and their respective geometric centers have a fixed spatial offset from the bottom tip of the centering rod. Before operation, the spatial eccentricity vectors from the phase center of the GNSS antenna to the bottom tip and from the reflection center of the omnidirectional prism to the bottom tip must be pre-determined in a static calibration site using precise coordinate measurement. Both sets of vectors are stored in the configuration file of the acquisition terminal in the form of three-dimensional components in the coordinate system of the centering rod body, with the unit being millimeters. The calibration results remain fixed throughout the entire operation cycle and serve as the original input parameters for subsequent rotation transformations.

[0039] In this embodiment, the total station is set up at a fixed station. The station is set up at one of the precisely measured physical invariant reference points within the survey area, and its three-dimensional coordinates are pre-included in the reference point database. After the total station's measurement controller completes the initial station orientation, its station coordinates and orientation azimuth remain unchanged throughout the entire operation cycle. At a normal epoch, the total station servo system drives the telescope to track the omnidirectional prism of the centering rod, acquiring the homogeneous coordinate components of the rod tip's bottom layer. Within the verification time window triggered by the residual comparison module, the total station servo system automatically rotates to the remaining physical invariant reference points within the survey area according to the preset observation sequence for automatic repetition observation. At this time, the measured horizontal distance and measured azimuth are both based on the fixed station coordinates. Correspondingly, the theoretical value reverse calculation module does not directly compare the measured values ​​with the GNSS pole tip coordinates as the starting point. Instead, it first performs a geometric conversion between the GNSS bottom-level source coordinates and the known coordinates of the total station fixed station in the same coordinate system to obtain the theoretical azimuth and theoretical horizontal distance independently calculated from the total station fixed station through the GNSS observation link. Then, it compares the theoretical azimuth and theoretical horizontal distance with the measured azimuth and measured horizontal distance of the total station at the reference point with the same reference, thereby ensuring that the observation starting points of the two sides are consistent.

[0040] The time synchronization module uses the second pulse signal output by the GNSS receiver as a unified time reference. When each pulse rises, it synchronously triggers three data latching commands to latch the dynamic calculated coordinates output by the GNSS calculation module at the current epoch, the original horizontal and vertical angles and slant range observations of the total station measurement controller before polar coordinate transformation, and the roll, pitch, and heading angles output by the attitude calculation module of the inertial measurement unit. The three data streams are written into the synchronization data buffer after sharing the same timestamp, ensuring that the data from each sensor used in subsequent fusion calculations strictly come from the same physical time and eliminating the time inconsistency bias introduced by asynchronous acquisition.

[0041] The attitude calculation module reads the roll, pitch, and yaw angles of the current epoch from the synchronization data buffer and inputs them into the rotation matrix construction unit in single-precision floating-point format. This unit constructs rotation matrices based on the roll angle about the vertical axis of the centering rod body coordinate system, rotation matrices based on the pitch angle about the horizontal axis, and rotation matrices based on the yaw angle about the vertical axis. The three 3×3 rotation matrices are stored in the processor register in standard orthogonal matrix form. Then, they are multiplied on the right in a fixed order of yaw angle rotation matrix, pitch angle rotation matrix, and roll angle rotation matrix to obtain a rigid body pose isomorphism mapping matrix that can map any vector in the centering rod body coordinate system to the engineering construction coordinate system. This matrix is ​​also written into the processor register in 3×3 orthogonal matrix form for subsequent eccentric vector rotation transformation.

[0042] The eccentricity compensation calculation module reads the rigid body pose isomorphic mapping matrix from the processor register and the pre-calibrated spatial eccentricity vectors from the GNSS antenna phase center to the bottom tip and from the omnidirectional prism reflection center to the bottom tip from the configuration file. It performs matrix and column vector multiplication operations on the two eccentricity vectors and the rigid body pose isomorphic mapping matrix respectively, and rotates and transforms the eccentricity vectors in the body coordinate system to the engineering construction coordinate system. It outputs two sets of spatial eccentricity compensation quantities, each of which is a three-dimensional component column vector in the engineering construction coordinate system, with the unit being millimeters. The rotation transformation results are written to the eccentricity compensation data area for use in subtraction operations.

[0043] The total station coordinate reduction module reads the original observation values ​​of horizontal angle, vertical angle, and slant distance latched by the total station from the synchronous data buffer. Based on the known coordinates and orientation azimuth of the station, it reduces the polar coordinate observation values ​​to three-dimensional rectangular coordinates in the engineering construction coordinate system. The reduction result is written into the input data area of ​​the coordinate fusion module in the form of three-dimensional coordinate triplets. The coordinate fusion module simultaneously reads the GNSS dynamically calculated coordinates from the synchronous data buffer and converts them from WGS-84 ellipsoidal coordinates to three-dimensional coordinate triplets in the engineering construction coordinate system. Then, it performs component-wise subtraction on the GNSS three-dimensional coordinate triplets and the corresponding GNSS antenna eccentricity compensation, and performs component-wise subtraction on the total station reduced three-dimensional coordinate triplets and the corresponding omnidirectional prism eccentricity compensation. The output results of the two sets of subtraction operations are three-dimensional coordinate triplets of the centering rod tip in the engineering construction coordinate system. The two together constitute the underlying homogeneous coordinates and are written into the input buffer of the subsequent pseudo-fixed solution determination module in the form of structured data records. The records contain three fields: timestamp, underlying homogeneous coordinates from GNSS source, and underlying homogeneous coordinates from total station source.

[0044] By driving the synchronous latching of three-channel sensor data with a second pulse signal and combining the real-time attitude of the inertial measurement unit to dynamically rotate and compensate the installation eccentricity vectors of the two types of sensors, the systematic deviations introduced by the asynchronous acquisition of sensors and installation eccentricity are effectively reduced. This ensures that the coordinates of GNSS and total station are consistent in terms of time and space references, providing a reliable comparison basis for subsequent pseudo-fixed solution determination.

[0045] Further, the reverse calculation of the theoretical horizontal distance and theoretical positive azimuth from the starting point to the pre-determined physical invariant reference point includes: projecting the GNSS coordinates in the underlying homogeneous coordinate system to the engineering construction plane coordinate system to obtain the two-dimensional plane starting coordinates; extracting the known reference plane coordinates of the pre-determined physical invariant reference point in the engineering construction plane coordinate system; calculating the geometric distance between the two-dimensional plane starting coordinates and each of the known reference plane coordinates as the theoretical horizontal distance; calculating the coordinate azimuth based on the coordinate increment of the two-dimensional plane starting coordinates and each of the known reference plane coordinates, and correcting it by combining the meridian convergence angle and direction modification parameters of the current position of the survey area to obtain the theoretical positive azimuth.

[0046] Specifically, the theoretical value reverse calculation module reads the three-dimensional coordinate triplet from the GNSS source to the bottom tip of the centering rod at the current epoch from the underlying homogeneous coordinate buffer. This triplet is stored in the form of geodetic latitude, geodetic longitude, and ellipsoidal height components in the WGS-84 geodetic coordinate system, with units of degrees and meters, respectively. This serves as the original input data for the projection transformation in this step. The projection transformation submodule converts the three-dimensional coordinate triplet in the geodetic coordinate system into two-dimensional plane coordinate pairs in the engineering construction plane coordinate system based on the pre-configured Gauss-Kruger projection parameters (including central meridian longitude, projection scale factor, and plane coordinate addition constant) of the survey area. In this step, the elevation component is separately extracted and stored in the elevation buffer for later use. The two-dimensional plane starting coordinates output by the transformation are written into the input register of the theoretical value calculation module in the form of double-precision floating-point pairs, with the unit being millimeters.

[0047] After receiving the two-dimensional plane starting coordinates, the theoretical value calculation module reads the known reference plane coordinates of all pre-determined physical invariant reference points in the survey area in batches from the reference point database in the local storage module. The database is stored in structured text format, and each record contains three fields: reference point number, north coordinate component, and east coordinate component, in millimeters. The database is entered and encrypted by the survey control engineer before the operation and is only available for read-only access during the operation. The theoretical value calculation module loads all reference point coordinate records into a memory array in sequence. The array length corresponds to the total number of reference points, and each array element is a reference point coordinate record.

[0048] The theoretical horizontal distance calculation submodule takes the two-dimensional plane starting coordinates of the current epoch as the starting point, calculates the north coordinate increment and east coordinate increment between the starting point and the reference point for each reference point coordinate record in the memory array, squares the two increments respectively, sums them and then takes the square root to obtain the planar geometric distance between the two points, which is used as the theoretical horizontal distance corresponding to the reference point. The theoretical horizontal distance is stored in the corresponding subscript position of the theoretical horizontal distance array in the form of a single-precision floating-point number, with the unit being millimeters. The array length is the same as the total number of reference points. After the theoretical horizontal distance calculation of all reference points is completed, the array is transferred to the residual comparison module for later use.

[0049] The theoretical positive azimuth calculation submodule starts with the two-dimensional plane coordinates of the current epoch and ends with the physical constant reference point. It extracts the north and east coordinate increments (i.e., reference point coordinates minus the starting point coordinates) for each reference point coordinate record in the memory array, and calls the arctangent function to calculate the coordinate azimuth from the starting point to the reference point, normalizing the value range to 0 to 360 degrees. Subsequently, the azimuth correction submodule applies meridian convergence angle correction and direction modification correction to the coordinate azimuth sequentially. Since the calculation direction is consistent with the measured direction of the total station from the station to the reference point, no 180-degree reverse modification is applied. The final theoretical positive azimuth is stored as a single-precision floating-point number at the corresponding index position in the theoretical positive azimuth array, in degrees. The array length is the same as the total number of reference points. After the theoretical positive azimuth calculation for all reference points is completed, this array and the theoretical horizontal distance array are passed to the residual comparison module for point-by-point comparison with the total station's synchronous measured values.

[0050] By reverse-calculating the theoretical geometric relationship of each physical invariant reference point after rigorous projection transformation and orientation modification of the GNSS underlying source coordinates, a calculation benchmark consistent with the field coordinate system is provided for the subsequent quantitative residual comparison between measured values ​​and theoretical values. This enables the determination of pseudo-fixed solutions to have clear geometric traceability and reduces the interference of systematic errors introduced by coordinate system inconsistency on the determination results.

[0051] Further, the step of determining whether the GNSS coordinates are pseudo-fixed solutions includes: obtaining the measured horizontal distance and measured azimuth angle of the total station synchronously measured to each of the physical invariant reference points at the same timestamp, comparing them with the corresponding theoretical horizontal distance and theoretical positive azimuth angle, and calculating the spatial position residuals of each reference point; when the maximum value of each spatial position residual is greater than the dynamic accuracy threshold, the GNSS coordinates are determined to be pseudo-fixed solutions and discarded; when the total station optical mapping data is valid, the coordinates are output according to the optical polar coordinates of the total station; when the total station loses synchronization, the attitude data output by the inertial measurement unit and the historical epochs are fused to perform short-time coordinate extrapolation, including: using the three-axis acceleration and three-axis angular velocity output in real time by the inertial measurement unit, and the initial velocity and initial position of the continuous historical epochs before the GNSS coordinates are determined to be pseudo-fixed solutions, and combining the rigid body pose isomorphic mapping matrix to perform inertial navigation integration calculation, to obtain the continuously recursively extrapolated short-time extrapolated coordinates during the period of total station loss of synchronization.

[0052] Specifically, within the same verification time window, the residual comparison module drives the total station to automatically measure at least two physically invariant reference points in the survey area according to a preset observation sequence. This data is organized in a structured record format, with each record containing four fields: reference point number, measured horizontal angle, measured vertical angle, and measured slope distance. The total station measurement controller converts the above-mentioned original observation values ​​into measured horizontal distance and measured azimuth in real time based on the current station coordinates and orientation azimuth. The observation values ​​of each reference point are then time-aligned with the underlying source coordinates of the GNSS source according to the corresponding acquisition timestamp through an interpolation model or time offset compensation. Subsequently, they are written into the measured horizontal distance array and the measured azimuth array, respectively. The array subscripts correspond one-to-one with the reference point numbers, and the array length is equal to the total number of physically invariant reference points currently participating in the comparison. The above two arrays, as the measured input data of the residual comparison module, are aligned with the theoretical horizontal distance array and theoretical positive azimuth array transmitted by the theoretical value calculation module within the same processing cycle.

[0053] The residual calculation submodule subtracts element-by-element from the measured horizontal distance array and the theoretical horizontal distance array to obtain the horizontal distance residual for each benchmark point. Simultaneously, it subtracts element-by-element from the measured azimuth array and the theoretical positive azimuth array, normalizing the difference to the range of -180 degrees to +180 degrees to obtain the azimuth residual for each benchmark point. Then, it calculates the horizontal distance residual and azimuth residual for each benchmark point according to preset distance-angle joint weights (e.g., horizontal distance residual weight 0.6, azimuth residual weight converted to arc length weight 0.4), extracts historical remeasurement data from a large number of benchmark points using the total station, and calculates the "distance component" and "angle-converted arc length component" respectively. The variance of the measured unit weight is calculated. Since the variance of distance measurement is usually slightly smaller than the equivalent variance of angle measurement, after iterative calculation, the optimal relative weight naturally converges to a level of about 6:4. This indicates that in most conventional construction scenarios, the optimal spatial residual merging can be achieved with extremely low computing power. The spatial position residual of the benchmark point is obtained by weighted merging and stored in the spatial position residual array in the form of single-precision floating-point numbers. After the spatial position residuals of all benchmark points are calculated, the judgment module extracts the maximum value from the array and compares it with the dynamic precision threshold passed by the dynamic precision threshold module in the current epoch.

[0054] When the maximum value in the spatial position residual array exceeds the dynamic accuracy threshold, the determination module marks the current epoch GNSS source underlying homologous coordinates as a pseudo-fixed solution, sends a discard command to the coordinate output control module, and simultaneously queries the tracking lock flag in the total station status register. When the total station tracking lock flag is in a valid state, the coordinate output control module switches the output source, reads the three-dimensional coordinate triplet of the current epoch calculated from the total station optical polar coordinates to the bottom tip of the centering rod from the total station source underlying homologous coordinate buffer, and uses this coordinate triplet to replace the GNSS source coordinates as the valid output coordinates of the current epoch. It writes the coordinates into the acquisition results database and adds a data source identifier field, setting the identifier value to "total station independent mode" for post-processing verification.

[0055] When the total station tracking lock flag is simultaneously in a unlocked state, the coordinate output control module activates the inertial navigation short-time calculation submodule. This submodule reads the real-time three-axis acceleration and three-axis angular velocity output from the inertial measurement unit data buffer. The three-axis acceleration is stored as a triplet in meters per square second, and the three-axis angular velocity is stored as a triplet in radians per second. At the same time, this submodule reads the effective initial velocity sequence and initial position sequence from the historical epoch database for the most recent consecutive epochs (e.g., 5 epochs) before the GNSS coordinates were determined to be pseudo-fixed solutions. The initial velocity is stored as a three-dimensional velocity triplet in meters per second, and the initial position is stored as a three-dimensional coordinate triplet in the engineering construction coordinate system in millimeters. The above data is used as the initial state input for the inertial navigation integration calculation.

[0056] The inertial navigation integration submodule reads the rigid body pose isomorphic mapping matrix calculated for the current epoch from the attitude calculation module's register. It then transforms the three-axis acceleration output by the inertial measurement unit from the body coordinate system to the engineering construction coordinate system. After compensating for the gravitational acceleration component, it obtains the motion acceleration for the current epoch. Subsequently, it performs a single numerical integration on the motion acceleration starting from the initial velocity of the historical epoch to obtain the estimated velocity triplet for the current epoch. Then, it performs a second numerical integration on the estimated velocity starting from the initial position of the historical epoch to obtain the short-time estimated coordinates for the current epoch. These coordinates are written to the coordinate output control module in the form of a three-dimensional coordinate triplet in the engineering construction coordinate system. The data source identifier field is set to "IMU short-time estimation mode". The accuracy of the above estimated coordinates gradually decreases due to the characteristic of inertial navigation integration error accumulating over time. When the number of consecutive lost-lock epochs exceeds the preset limit (e.g., 10 consecutive epochs), the submodule pushes a failure warning to the acquisition terminal display interface to remind the operators to verify the status of the on-site equipment. The short-time estimated coordinates are no longer written to the acquisition results database.

[0057] By constructing a cross-validation mechanism based on an externally physical, unchanging benchmark, and using inertial navigation integration for short-term coordinate recursion as a fallback when both GNSS and total station fail, a multi-level degradation output strategy is formed. This enables the acquisition system to maintain the continuity of coordinate output even in complex interference environments, effectively reducing the risk of acquisition interruption or gross error output due to the failure of a single sensor.

[0058] Furthermore, the step of constructing a dynamic health assessment model and calculating the health score of each temporary benchmark point within the survey area includes: obtaining a material disturbance resistance factor characterizing the physical displacement resistance of each benchmark point based on the geological and engineering structural properties of the carrier to which it is attached; obtaining the measured coordinates of each temporary benchmark point periodically re-measured by a total station, and comparing them with the initial coordinates to obtain spatial mapping residuals; applying coordinate residual penalties to the spatial mapping residuals that exceed the allowable tolerance of mapping based on robust estimation theory, and obtaining a mapping reliability weight that is negatively correlated with the spatial mapping residuals; and constructing an exponential time decay factor characterizing the increase in the probability of natural settlement failure of each temporary benchmark point over time based on the mechanical vibration and soil rheological characteristics of the construction site. The health score of each temporary reference point is obtained by multiplying the material anti-disturbance factor, the mapping reliability weight, and the exponential time decay factor. The temporary reference point pool is updated according to the health score, and failed points are removed. This includes: in each dynamic re-measurement cycle of the survey area control network, the temporary reference points are optimized and sorted according to the health score, and a set of stable physical measurement points with qualified health scores is extracted; the three-dimensional coordinates of the set of stable physical measurement points are included in the stable temporary measurement point level, updated in the station setting and starting parameters of the total station, and used as a weighted verification point with a lower weight than the permanent physical control point in the subsequent pseudo-fixed solution determination; temporary reference points with unqualified health scores are discarded from the topology network of the current spatial intersection measurement.

[0059] Specifically, the first input source for the health assessment model is the material disturbance resistance factor. The acquisition process of this factor is completed by the benchmark point attribute entry module during the work preparation stage. After the workers complete the initial installation of each temporary benchmark point on site, they select the category of the carrier to which the benchmark point is attached through the attribute entry interface of the acquisition terminal. The carrier category is stored in the configuration file in the form of enumerated values, including reinforced concrete embedded parts, steel pipe piles, wooden piles, masonry attachment points, and temporary observation piers, etc. Each category corresponds to a preset initial value of the material disturbance resistance factor. For example, the initial value of reinforced concrete embedded parts is 0.95, the initial value of steel pipe piles is 0.85, and the initial value of wooden piles is 0.6. The above initial values ​​and subsequent health thresholds are all obtained by fitting the residual samples of historical retests. The initial value reflects the relative ability of this type of carrier to maintain spatial position stability under construction vibration and foundation disturbance conditions. After the entry is completed, it is written into the material disturbance resistance factor field of the benchmark point database in the form of a single-precision floating-point number. This field remains unchanged in each subsequent retest cycle and is only manually updated by the workers when the carrier category changes physically.

[0060] The second input source for the health assessment model is the surveying reliability weight. This weight is calculated by the total station re-survey submodule in each dynamic re-survey cycle. The total station re-survey submodule automatically retrieves the observation sequence of all temporary benchmarks in the survey area by the total station according to the preset re-survey time interval (e.g., once every 4 hours). For each temporary benchmark, the horizontal angle, vertical angle and slope distance are measured in sequence, and the measured polar coordinates are converted into three-dimensional re-survey coordinates under the engineering construction coordinate system and written into the re-survey coordinate database in the form of three-dimensional coordinate triplets. The residual calculation submodule reads the re-survey coordinates of the current re-survey from the re-survey coordinate database, and at the same time reads the initial coordinate triplets of the benchmark after the first installation by the benchmark database. The two sets of triplets are subtracted component by component and the modulus is taken to obtain the spatial surveying residual of the current re-survey cycle, in millimeters. This residual is written into the residual record array as a single-precision floating point number.

[0061] The robustness penalty submodule reads the spatial mapping residuals of each temporary benchmark point from the residual record array and reads the corresponding mapping tolerance limit value (e.g., 5 mm for the core area and 15 mm for the general area) from the survey area configuration file. For benchmark points whose spatial mapping residuals do not exceed the tolerance limit, the mapping reliability weight is linearly mapped to the range of 0.7 to 1.0 within a preset range according to the residual size. The smaller the residual, the higher the corresponding weight. For benchmark points whose spatial mapping residuals exceed the tolerance limit, the robustness penalty submodule applies a nonlinear coordinate residual penalty to the excess part according to the robustness estimation theory. The penalty coefficient decreases rapidly as the excess increases, so that the mapping reliability weight of the corresponding benchmark point is reduced to below 0.7 until it approaches 0. The penalized mapping reliability weight is written to the corresponding field of the benchmark point database as a single-precision floating-point number for use in the current health score calculation.

[0062] The third input source for the health assessment model is the exponential time decay factor. This factor is calculated by the time decay modeling submodule at the beginning of each retest cycle. This submodule reads the initial burial timestamp of each temporary benchmark from the benchmark database, calculates the time difference between the current retest time and the initial burial time in hours, and reads the mechanical vibration intensity level and soil rheological coefficient of the current test area from the test area construction environment parameter configuration file. The mechanical vibration intensity level is represented by integers from 1 to 5 (1 for low vibration environment and 5 for strong vibration environment), and the soil rheological coefficient is stored in dimensionless decimal form. The two together determine the rate parameter of exponential decay. The time decay modeling submodule substitutes the time difference and rate parameter into the negative exponential decay model and outputs the exponential time decay factor of the benchmark at the current time. The value range is a single-precision floating-point number between 0 and 1. The larger the time difference, the higher the vibration intensity, or the larger the soil rheological coefficient, the smaller the decay factor. This factor is written into the corresponding field of the benchmark database.

[0063] After the calculation of the three types of factors in each retesting cycle is completed, the material anti-disturbance factor, mapping reliability weight, and exponential time decay factor of each temporary benchmark point are read one by one from the benchmark point database. The three are then multiplied element by element to obtain the health score of the benchmark point in the current retesting cycle. The value range is a single-precision floating-point number between 0 and 1. Here, a lower limit protection rule is introduced: After the health score synthesis submodule completes the calculation of the three factors, it reads the spatial mapping residual records of the most recent K consecutive retests of the temporary benchmark point (3 are selected in this embodiment, and can be adjusted according to the actual situation). If the K residuals are all less than the mapping tolerance limit of the corresponding area, the health score obtained by this multiplication is compared with the preset protection threshold (e.g., 0.75), and the larger value is taken as the final health score of the benchmark point written into the database. That is, the score decrease caused by the exponential time decay factor is limited to above the protection threshold in this case. If there are records in the K retest residuals that exceed the tolerance limit, the protection rule is not triggered, and the health score is directly written according to the calculation result of the three factors multiplication. A higher health score indicates that the current spatial location of the benchmark is more reliable. The score calculation result is written into the health score field of the benchmark database and a timestamp of the current retest is attached to form a time series record of the health score of each benchmark.

[0064] The benchmark pool update module starts immediately after the health score is written. It reads the current health score of all temporary benchmarks from the benchmark database, sorts them in descending order of score, and includes temporary benchmarks with health scores greater than or equal to a preset threshold (e.g., 0.75, the value selected in this embodiment, which can be adjusted according to the actual situation) into the stable physical measurement point set. This method establishes a three-level benchmark control network: permanent physical control points, stable temporary measurement points, and candidate / abandoned temporary points. The pseudo-fixed solution determination prioritizes cross-verification using permanent physical control points. The coordinates of the aforementioned re-measured and confirmed stable physical measurement points are written into the station setup starting parameter data area of ​​the total station controller as stable temporary measurement points, and used as weighted verification points when participating in the pseudo-fixed solution determination. After receiving the update instruction, the total station controller will execute subsequent station orientation and polar coordinate reduction with the new starting parameters from the next acquisition epoch. For temporary benchmark points with health scores below the preset threshold, the abandonment processing submodule removes their numbers from the node list of the current spatial intersection measurement topology network. The benchmark point array of the subsequent residual comparison module will no longer contain the record corresponding to the number, thereby avoiding the observation data of failed benchmark points from entering the comparison calculation of the pseudo-fixed solution determination. The abandonment operation result is synchronously written to the status field of the benchmark point database, the status value is set to "abandoned", and an abandonment timestamp is attached for post-processing review.

[0065] By integrating three quantitative factors—material physical properties, periodic remeasurement residual penalty, and time decay—into a unified health score, and driving the dynamic screening of the benchmark pool and the synchronous update of the total station's starting parameters in each remeasurement cycle, the survey area control network can continuously eliminate benchmarks with declining reliability and supplement them with verified stable measuring points. This effectively suppresses the risk of errors caused by failed benchmarks being transmitted step by step in the measurement chain, and improves the long-term stability of the topographic data acquisition results of the entire survey area.

[0066] In summary, this embodiment constructs a complete technical chain for terrain data acquisition in complex construction sites through four mutually supportive functional links: signal quality adaptive mode switching, multi-source data synchronous latching and rigid body eccentricity compensation, cross-verification judgment based on physically invariant benchmarks, and dynamic health management of temporary benchmarks. The output data of each link serves as the input basis for the next link, forming a closed-loop data flow from raw satellite observation data to reliable coordinate results. This effectively reduces the impact of various error sources such as multipath effects, asynchronous sensor acquisition, installation eccentricity, and benchmark failure on the accuracy of the final results, and realizes the reliability of multi-sensor joint acquisition in terrain mapping and construction layout scenarios.

[0067] Example 2: During the steel structure hoisting operation in the core tube shear wall area of ​​a super high-rise building, the large-scale obstruction caused by the frequent rotation of the tower crane caused the number of visible GNSS satellites in the survey area to fluctuate drastically in a short period of time. At the same time, the hoisting vibration and foundation disturbance caused the health scores of multiple temporary reference points in the survey area to decline continuously over several consecutive retesting cycles. The number of available stable measurement points in the existing reference point pool approached the minimum lower limit required for spatial intersection calculation. The system faced a complex and difficult scenario with the simultaneous superposition of insufficient reference point support and intermittent degradation of GNSS signals. This embodiment provides a detailed explanation of the dynamic weight adaptive adjustment mechanism of the quality factor matrix and the emergency expansion strategy of the reference point pool for the above scenario.

[0068] Specifically, under normal operating conditions, the weighting coefficients of the three columns of spatial position accuracy factor, cycle slip comprehensive influence, and multipath error in the quality factor matrix operate with a fixed configuration. However, in the occlusion scenario caused by tower crane hoisting operations, the occlusion direction changes continuously with the tower crane rotation angle, causing the visible satellite set to frequently change between adjacent epochs. The contribution of the spatial position accuracy factor component to the overall signal health comprehensive index under the fixed weighting coefficients undergoes a non-physical abrupt change at the moment the occlusion direction changes, causing the mode switching decision module to generate frequent mode jitter. The acquisition controller repeatedly switches between RTK fixed solution acquisition mode and tightly coupled solution acquisition mode, resulting in a decrease in the continuity of the coordinate output sequence.

[0069] To address the aforementioned issues, a dynamic weight adjustment submodule for the quality factor matrix is ​​activated. This submodule continuously monitors the flip frequency of mode switching decisions over several consecutive epochs (e.g., 10 epochs). When the flip frequency exceeds a preset stability threshold (e.g., 4 or more mode switches within 10 epochs), it is determined that the system is currently in a dynamic occlusion disturbance state. The submodule reads the column weight configuration array of the current quality factor matrix from the processor memory, lowers the weight of the spatial position accuracy factor column from 0.4 to 0.25, and raises the weight of the multipath error column from 0.25 to 0.40. The weight of the cycle slip comprehensive impact column remains unchanged at 0.35. The adjusted weight configuration array is written back to the processor memory and takes effect immediately, making the calculation of the signal health comprehensive index more sensitive to changes in multipath error and relatively narrowing the response amplitude to transient changes in satellite geometry, thereby suppressing the jitter frequency of mode switching decisions. The adjusted weight configuration is written to the local storage module in log form with the current epoch timestamp appended, for post-processing stage verification of weight change history.

[0070] Furthermore, in scenarios where the number of available stable measurement points in the benchmark pool approaches the lower limit, the benchmark pool emergency expansion submodule is triggered synchronously. This submodule reads all temporary benchmark records with the current status field set to "abandoned" from the benchmark database, extracts their health score and abandonment timestamp from the most recent retest, and filters out candidate revival points whose abandonment time is no more than a preset backtracking window (e.g., 2 retest cycles) and whose health score is below a preset threshold but not below the emergency activation lower limit (e.g., health score between 0.55 and 0.75). The numbers of the above candidate revival points are written into the emergency candidate array.

[0071] After the emergency candidate array is generated, the total station resurvey submodule immediately performs a temporary encrypted resurvey on each candidate resurrection point in the array. After the resurvey is completed, the residual calculation submodule recalculates the current spatial mapping residual of each candidate point and updates its mapping reliability weight. The health score synthesis submodule recalculates the health score of each candidate point with the updated mapping reliability weight. For candidate resurrection points whose health scores meet the temporary activation threshold (e.g., not less than 0.65) after recalculation, the benchmark pool update module temporarily changes its status field from "abandoned" to "emergency activation" and includes its three-dimensional coordinates in the node list of the current spatial intersection measurement topology network. At the same time, the coordinate record of this point is added to the station setting and starting parameter data area of ​​the total station measurement controller, and the participation weight is marked as reduced participation (e.g., in the weighted calculation of the subsequent residual comparison, the weight of this point is taken as 0.6 times that of the normal stable point), so that the topology network can still maintain the necessary geometric strength when the number of stable measurement points is insufficient, ensuring the normal progress of spatial intersection solution.

[0072] The reference point in the emergency activation state continues to be periodically retested and updated with health score by the total station in each subsequent regular retest cycle. If the health score rises above the normal stable threshold for several consecutive retest cycles (e.g., 3 consecutive cycles), the reference point pool update module upgrades its status field from "emergency activation" to normal stable state, and the weight is restored to the standard value. If the health score continues to decline after emergency activation and falls below the emergency activation lower limit, the submodule resets its status field to "abandoned" and removes it from the topology node list. All the above status change records are written to the local storage module in the form of timestamp logs.

[0073] This embodiment effectively suppresses frequent jitter during mode switching in dynamic occlusion scenarios by introducing a dynamic adaptive adjustment mechanism for the column weights of the quality factor matrix. At the same time, by implementing a conditional emergency revival and deweighting participation strategy for abandoned reference points, it provides a controllable geometric support supplement to the spatial intersection measurement topology when the number of stable measurement points is insufficient. This enables the system to maintain continuous operation of the acquisition process even in the complex and difficult scenario of simultaneous superposition of GNSS signal degradation and limited available reference points, thus improving the overall solution's adaptability to extreme construction interference environments.

[0074] Example 3: In large-scale construction sites, multiple surveying teams simultaneously conduct topographic data acquisition. Each team uses independent integrated centering rod equipment. Various types of heavy machinery, such as tower cranes, pile drivers, and concrete pump trucks, operate concurrently within the survey area. The dynamic operating areas of these machines frequently overlap with the work areas of the surveying teams, leading to risks such as equipment being obstructed by robotic arms, being moved by construction workers, or the total station's line of sight being cut off by temporary structures. Simultaneously, temporary reference points within the survey area are shared by multiple teams. If each team independently maintains its own reference point status, the health assessment results of the same reference point from different teams will conflict due to asynchronous retesting times. A failed reference point may be reused by another team even after being abandoned by one team. This embodiment provides a detailed explanation of the spatiotemporal occupancy management mechanism for the shared temporary reference point pool and the real-time alert function for collectable sectors under dynamic mechanical sensing, addressing the aforementioned scenario of multiple teams working together.

[0075] Specifically, the system establishes a unified shared temporary reference point pool database at the survey area-level management server. This database runs on an edge computing server deployed on-site in the survey area. The acquisition terminals of each group maintain real-time data synchronization with the server via a wireless LAN. The table structure of the shared reference point pool database adds three fields to the independent operation scenario of each group: the current occupying group number field, the occupying timestamp field, and the operation area spatial range field. The operation area spatial range field is stored in the form of a three-dimensional polygon vertex list under the engineering construction coordinate system, describing the spatial range currently actually covered by the group's operation. When each group's acquisition terminal enters the operation state, it sends an occupation registration request to the server. The request data includes the group number, the expected spatial range polygon of the operation area, and the expected operation duration. After receiving the request, the server performs spatial overlay detection in the shared reference point pool database to determine whether the operation area declared by the group overlaps with the operation areas of other groups that have already registered and occupied. If an overlap is detected, a conflict warning is pushed to the requesting group's terminal, marking the spatial range of the conflict area and the current occupying group number. The operators coordinate and adjust the registration and resubmit the registration. After the conflict detection passes, the server writes the group's occupation record into the database and broadcasts an update notification to all online terminals.

[0076] The occupancy lock management module continues to run after each group completes its occupancy registration. Its core function is to prevent data conflicts caused by concurrent retesting and status modification of the same temporary benchmark point by multiple groups. When a group's total station retesting submodule initiates a retest write request for a specific temporary benchmark point to the server, the occupancy lock management module first checks whether the benchmark point is currently locked by other groups' retesting tasks. If it is locked, the request is added to a waiting queue. After the lock is released, the write is executed sequentially according to the queue order, ensuring that only one group's retest result is written to the database for the same benchmark point at any given time. The health score synthesis submodule of each group calculates the score after... After completion, the shared benchmark pool is updated through the unified write interface on the server. When the server receives the health scores of the same benchmark from different groups, it merges the evaluation results of multiple groups in a weighted average manner. The weight is determined based on the number of consecutive epochs locked by the total station when each group submits the retest data. The more consecutive epochs locked, the higher the merge weight of the results submitted by that group. The merged health score is written back to the shared benchmark pool database and broadcast to all terminals. After receiving the broadcast, each group's terminal synchronously refreshes its local benchmark status cache to ensure that the benchmark reliability evaluation results called by each group always come from the latest comprehensive evaluation across the entire survey area.

[0077] Furthermore, the mechanical dynamic perception module obtains real-time operating status data of various heavy machinery from the tower crane operation management system and pile driver positioning system at the construction site. The tower crane operation management system pushes the current slewing angle, boom length, and hook height of each tower crane to the edge computing server at fixed time intervals (e.g., every 2 seconds). The pile driver positioning system pushes the current planar coordinates and operating influence radius of each pile driver. After receiving the above data, the mechanical dynamic perception module on the server calculates the spatial sector currently covered by the boom of each tower crane in real time under the engineering construction coordinate system based on the tower crane's slewing angle and boom length. It stores the sector in the mechanical occupancy space cache in the form of a list of vertex coordinates of a three-dimensional sector and synchronously writes the operating influence range of the pile driver into the cache in the form of a cylindrical space bounding box. The above mechanical occupancy space data is refreshed with each push.

[0078] The real-time calculation submodule for the collectable sector reads the current space occupancy data of all machines from the mechanical space occupancy buffer at each acquisition epoch. Combined with the current spatial coordinates reported by each group's acquisition terminal, it calculates the line-of-sight vectors pointing from each terminal location to each temporary reference point and the physical invariant reference point. It then checks whether each line-of-sight vector crosses any mechanical space. Crossing is determined by the intersection test between the ray and the 3D geometry. If there is an intersection, the line-of-sight in that direction is currently obstructed by machinery, and the corresponding reference point direction is marked as uncollectible. If there is no intersection, it is marked as collectable. All directions are considered collectable. The data collection results are organized into a list of azimuth intervals as collectable sector description data. This data is written to the local status cache of each terminal using the current epoch timestamp as an index, and is simultaneously pushed to the display interface of the corresponding group's data collection terminal. The display interface renders a 360-degree azimuth ring diagram centered on the current terminal position. Collectable sectors are highlighted in green, non-collectable sectors obstructed by machinery are highlighted in red, and areas affected by locks occupied by other groups are highlighted in yellow. Based on this, operators can determine the effective measurement direction range of the current station in real time and reasonably select the measurement time or adjust the station position.

[0079] When the real-time calculation submodule of the acquireable sector determines that the line-of-sight directions pointing to all physically invariant reference points within the current epoch are in an unacquireable state, causing the pseudo-fixed solution determination module to be unable to obtain effective measured comparison data, the module sends a pause determination command to the acquisition controller. The GNSS coordinates of the current epoch will not be subject to pseudo-fixed solution determination for the time being. The coordinate output control module marks the output of the current epoch as [pending verification and temporary storage] and writes it into the temporary storage queue. The comparison and verification will be performed after the acquireable sector is restored. If the data in the temporary storage queue is not verified within the preset waiting time limit (e.g., 30 seconds), a timeout warning will be pushed to the operators, indicating that the current station is in a state of mechanical obstruction for a long time, and it is recommended to change the station position or wait for the mechanical operation interval before continuing to acquire data.

[0080] This embodiment establishes a spatiotemporal occupancy registration and occupancy lock management mechanism for a shared temporary reference point pool at the survey area-level edge computing server, and couples the real-time operating status of heavy machinery at the construction site with the calculation of collectable sectors. This enables multiple surveying teams to share a unified reference point reliability assessment result under cross-operation conditions, avoiding data conflicts caused by each team independently maintaining the reference point status. At the same time, it provides real-time visualization prompts on the availability of measurement directions for each team's operators.

[0081] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A method for acquiring topographic data based on a combination of GNSS and total station, characterized in that, include: Acquire satellite observation data stream and optical polar coordinate observation data stream for the target point; The GNSS quality factor matrix is ​​calculated based on the satellite observation data stream of consecutive epochs, and the acquisition mode is automatically switched based on the quality factor matrix. Using a pre-imported site model, a dynamic accuracy threshold matching the current site is extracted based on the current spatial location. Within the same period, the dynamic solution coordinates of GNSS, the optical polar coordinates of the total station, and the attitude data output by the inertial measurement unit are synchronously latched. Rigid body geometric translation calculations are performed by combining the pre-calibrated antenna and the spatial eccentric vector from the omnidirectional prism to the bottom tip, and the underlying homogeneous coordinates are obtained respectively. Using the GNSS coordinates in the underlying homogeneous coordinates as the starting point, the theoretical horizontal distance and theoretical positive azimuth from the starting point to the pre-determined physical invariant reference point are calculated in reverse. The theoretical positive azimuth and theoretical horizontal distance are compared with the azimuth and horizontal distance measured synchronously by the total station to determine whether the GNSS coordinates are pseudo-fixed solutions. Calculate the health score of each temporary reference point in the survey area, update the temporary reference point pool based on the health score, and remove failed points.

2. The terrain data acquisition method based on GNSS and total station combined according to claim 1, characterized in that, The steps for automatically switching acquisition modes include: The satellite observation data stream of the continuous epochs is analyzed to extract the carrier signal-to-noise ratio, spatial position accuracy factor, and ionospheric delay rate of change of dual-frequency carrier phase observations for each visible satellite. A geometrically distance-free combined observation model is used to detect cycle slip indicators in the satellite observation data stream, and the multipath effect error is quantified by combining the ionospheric delay rate of change and the carrier signal-to-noise ratio. The spatial position accuracy factor, cycle slip indicator, and multipath effect error are used as dynamic parameters to construct a GNSS quality factor matrix, and the comprehensive signal health index for the current epoch is calculated. The acquisition mode is automatically switched based on the comprehensive signal health index.

3. The topographic data acquisition method based on GNSS and total station combined according to claim 1, characterized in that, The step of extracting the dynamic accuracy threshold that matches the current site includes: Obtain the relative local design coordinates from the BIM design file. Based on the known geoid elevation anomaly parameters and projection deformation modification parameters of the survey area, transform the relative local design coordinates to the engineering construction coordinate system to obtain a three-dimensional survey area model. According to the topological boundary of the three-dimensional survey area model, divide it into non-overlapping three-dimensional polygon bounding boxes and establish a spatial index using an octree structure. The tolerance values ​​of different accuracy requirements in the test area are used as attribute fields and bound to the topology nodes of the corresponding bounding boxes. The current spatial coordinates are used as probe points and input into the spatial index. The point-polygon inclusion test algorithm is used to lock the target three-dimensional polygon bounding box to which the probe point belongs, and the attribute fields of the bounding box are parsed to obtain the dynamic accuracy threshold.

4. The terrain data acquisition method based on GNSS and total station combined according to claim 1, characterized in that, The steps for obtaining the underlying homogeneous coordinates include: Using a second pulse signal as a time reference, the dynamic calculated coordinates of the GNSS, the optical polar coordinates of the total station, and the attitude data output by the inertial measurement unit are synchronously acquired. The attitude data includes roll angle, pitch angle, and yaw angle. Based on the roll angle, pitch angle, and yaw angle, rotation matrices are constructed around the three coordinate axes of the spatial rectangular coordinate system, and the three rotation matrices are multiplied to obtain the rigid body attitude isomorphic mapping matrix. The rigid body attitude isomorphic mapping matrix is ​​used to perform rotation transformation on the pre-calibrated spatial eccentricity vector from the antenna and the omnidirectional prism to the bottom tip, and the corresponding spatial eccentricity compensation is calculated. The spatial eccentricity compensation is subtracted from the GNSS dynamic calculated coordinates and the total station polar coordinates converted to spatial rectangular coordinates to form the underlying homogeneous coordinates.

5. A method for acquiring terrain data based on a combination of GNSS and total station as described in claim 1, characterized in that, The theoretical horizontal distance and theoretical positive azimuth from the reverse calculation starting point to the pre-determined physical invariant reference point are included: Project the GNSS coordinates in the underlying homogeneous coordinate system onto the engineering construction plane coordinate system to obtain the two-dimensional plane starting coordinates; extract the known reference plane coordinates of the pre-determined physical invariant reference points in the engineering construction plane coordinate system; calculate the geometric distance between the two-dimensional plane starting coordinates and each of the known reference plane coordinates as the theoretical horizontal distance; calculate the coordinate azimuth angle based on the coordinate increment of the two-dimensional plane starting coordinates and each of the known reference plane coordinates, and correct it by combining the meridian convergence angle and direction modification parameters of the current position of the survey area to obtain the theoretical positive azimuth angle.

6. The terrain data acquisition method based on GNSS and total station combined according to claim 1, characterized in that, The steps for determining whether the GNSS coordinates are pseudo-fixed solutions include: The measured horizontal distance and measured azimuth of the total station synchronously measured to each of the physical invariant reference points at the same timestamp are obtained and compared with the corresponding theoretical horizontal distance and theoretical positive azimuth to calculate the spatial position residual of each reference point. When the maximum value of each spatial position residual is greater than the dynamic accuracy threshold, the GNSS coordinates are determined to be pseudo-fixed solutions and discarded. When the optical mapping data of the total station is valid, the coordinates are output according to the optical polar coordinates of the total station. When the total station loses synchronization, the attitude data output by the inertial measurement unit and the historical epochs are fused to perform short-time coordinate extrapolation.

7. The terrain data acquisition method based on GNSS and total station combined according to claim 1, characterized in that, The steps for calculating the health score of each temporary reference point within the survey area include: Based on the geological and engineering structural properties of the carriers to which each temporary benchmark point is attached, a material disturbance resistance factor characterizing the benchmark point's resistance to physical displacement is obtained; the measured coordinates of each temporary benchmark point are obtained by periodic remeasurement by a total station, and compared with the initial coordinates to obtain spatial mapping residuals; according to robust estimation theory, coordinate residual penalties are applied to the spatial mapping residuals that exceed the allowable limits of mapping, and a mapping reliability weight that is negatively correlated with the spatial mapping residual is obtained. Based on the mechanical vibration and soil rheological characteristics at the construction site, an exponential time decay factor is constructed to characterize the increase in the probability of natural settlement failure of each temporary benchmark point over time. The material disturbance resistance factor, the surveying reliability weight, and the exponential time decay factor are multiplied to obtain the health score of each temporary benchmark point. The temporary benchmark point pool is updated according to the health score, and failed points are removed.

Images