Unmanned aerial vehicle field surveying method without ground control points by fusing real-time differential positioning and visual odometry

Abstract

The application discloses a kind of unmanned aerial vehicle field free control point surveying method of fusing real-time differential positioning and visual odometer.The method is with uniform time base synchronous acquisition differential positioning epoch, exposure feedback pulse and visual frame sequence;From visual odometer, interframe pose increment is calculated and visual quality factor is calculated;With the satellite positioning baseline of fixed solution epoch interval, visual scale factor and its relative mean error are calibrated online;With the visual pose increment after scale reduction, exposure time is motion compensated and rod arm is corrected;Then with solution, PDOP, lock loss time, visual quality factor and the scale relative mean error, variance model is constructed to weight, wherein the first term of lock loss time corresponds to random walk, quadratic term corresponds to scale drift;Further, four kinds of observation values are included in unified objective function to make free control point beam method adjustment.The application can inhibit regional network deformation and height system error caused by lock loss and time synchronization residual without laying control point.

Description

Technical Field

[0001] This invention relates to the field of aerial photogrammetry, and in particular to a UAV field mapping method that integrates real-time differential positioning and visual odometry without image control. Background Technology

[0002] Low-altitude photogrammetry using unmanned aerial vehicles (UAVs) is now widely used in large-scale topographic mapping (1:500 to 1:2000), 3D modeling, and engineering surveying and design. Traditionally, to integrate imagery into the national geodetic coordinate system and control regional network deformation, ground control points (hereinafter referred to as image control points) must be established and connected at a certain density within the survey area. The selection, installation, marking, and satellite positioning of image control points typically account for more than half of the total fieldwork for a surveying project. In steep cliffs, tidal flats, waterways, mining subsidence areas, high-voltage power transmission corridors, and densely vegetated forest areas, image control points are often impossible to establish. Even if they are established, personnel face considerable safety risks. Therefore, image control-free mapping has become an urgent production need.

[0003] Existing image control-free technology primarily uses airborne GNSS real-time dynamic differential positioning or post-processing dynamic differential positioning to obtain the camera projection center coordinates of each exposure point. These coordinates are then used as weighted absolute observations and directly incorporated into bundle adjustment, replacing ground-based image control points to provide an absolute reference and geometric constraints for the regional network—a process known as POS-assisted bundle adjustment. In terms of nominal specifications, real-time dynamic differential positioning can achieve a planar accuracy of ±8mm + 1ppm and an elevation accuracy of ±15mm + 1ppm, theoretically sufficient to support large-scale mapping.

[0004] Several improvement methods have been developed in this field. In image-controlled aerial triangulation, the camera principal distance and photographic distance are highly correlated. To weaken this correlation, existing research has employed multi-level or intersecting flight paths at different altitudes to supplement the geometric conditions for self-calibration adjustment. Related literature also reveals that the elevation systematic error in image-controlled results is approximately linearly related to the degree to which the principal distance estimate deviates from the true value. Regarding the uncertainty of exposure time, existing technologies send exposure test signals to the camera and receive its feedback synchronization signals to determine the exposure delay and accurate exposure time. Other studies have incorporated exposure delay as an additional systematic error parameter into GPS-assisted bundle adjustment models, taking into account the differences in delay caused by varying instantaneous velocities at different exposure points. In terms of observation weighting, classifying exposure point location observations into several categories—fixed solutions, floating-point solutions, and single-point solutions—based on solution status identifiers and assigning different empirical weights, as well as using weighted iterative methods to detect gross errors, are all standard practices in this field.

[0005] However, in real-world production scenarios where differential signals frequently lose lock, the aforementioned existing technologies still have the following unresolved issues.

[0006] Existing weighting methods do not align with the actual evolution of lock-up errors. In mountainous areas, canyons, dense forests, urban high-rise clusters, and near high-voltage power transmission corridors, satellite signals are easily blocked and reflected, leading to integer ambiguity lock-up. After lock-up, the solution state degenerates from a fixed solution to a floating-point solution or even a single-point solution, with positioning accuracy deteriorating from centimeter-level to decimeter-level or even meter-level. Reinitialization typically requires several seconds to tens of seconds. Furthermore, the positioning error caused by lock-up does not jump to a certain value immediately upon lock-up, but rather accumulates continuously as lock-up persists. Existing weighting methods based on solution state tiers, regardless of the number of tiers, result in a step change in weight at the moment of solution state switching, while remaining constant during the duration of the same solution state. This does not correspond to the actual situation where errors accumulate continuously over time. In the early stages of lock-up, observations are excessively downweighted, wasting usable positioning information; in the later stages, the weights of observations are excessively high, with severely degraded observations still participating in adjustment with significant weights, creating local distortions and elevation steps in the regional network. Gross errors can even spread to adjacent flight zones.

[0007] Extrapolating the exposure time lacks observations that can describe the actual maneuvering process. Airborne GNSS receivers typically output positioning epochs at frequencies between 1 Hz and 10 Hz, while the camera's exposure time falls randomly between two adjacent epochs. Existing methods mostly perform linear interpolation between adjacent epochs; however, when UAVs turn, encounter gusts of wind, or make attitude adjustments, the maneuvering acceleration can reach 2 m / s². 2 As mentioned above, linear interpolation cannot describe the acceleration term within the envelope interval, and its residual error can reach several centimeters to tens of centimeters at low epoch frequencies. The aforementioned practice of incorporating exposure delay as a systematic error parameter into the adjustment still relies on the velocity vector at the instant of exposure, which is essentially a first-order approximation and does not introduce independent observations that can truly characterize the maneuvering process within the envelope interval. Furthermore, this approach leaves the error to be absorbed by parameters during the adjustment stage, rather than eliminating it during the observation generation stage. Once this systematic error is related to the exterior azimuth element and the principal point offset, it will be transformed into an imperceptible systematic error in the results.

[0008] Purely visual methods cannot be introduced alone. Monocular visual odometry can recover the pose increment between adjacent images with high relative accuracy, which can make up for the two shortcomings mentioned above; however, its scale is unobservable, and it is prone to tracking failure and pose drift in areas with weak or repetitive textures, such as water, sand, and snow. Therefore, for visual pose increments to participate in adjustment as quantitative observations, two issues must be addressed: how to calibrate the scale, and how to quantify the accuracy of the calibrated scale itself and how it affects the observations derived from it. Existing technologies that combine visual or inertial methods with satellite positioning mostly focus on real-time navigation and positioning in satellite-denied environments, aiming to obtain continuous position solutions, without feeding back the calibration accuracy of the visual scale to the stochastic model of photogrammetric adjustment; even when visual methods are introduced, the weighting basis is still the solution state label itself, and the visual methods and the stochastic model are isolated from each other.

[0009] Current technology lacks a UAV mapping method that uses the short-term relative pose of visual odometry for extrapolation of the exposure time and quantitatively feeds back the calibration accuracy of the visual method itself to the adjustment stochastic model, thus ensuring the accuracy of the results even under conditions of frequent differential signal loss and no image control points.

[0010] This invention aims to solve the technical problems existing in the prior art. To this end, it proposes a UAV field mapping method that integrates real-time differential positioning and visual odometry. Summary of the Invention

[0011] The purpose of this invention is to provide a UAV field mapping method that integrates real-time differential positioning and visual odometry to solve the technical problems existing in the prior art.

[0012] To achieve the above objectives, the present invention provides the following technical solution:

[0013] A UAV field mapping method integrating real-time differential positioning and visual odometry includes the following steps:

[0014] S1. Airborne Multi-Source Data Synchronous Acquisition: A unified time base is established using satellite second pulses. Under this unified time base, four observation data streams are synchronously acquired. The four observation data streams include: differential positioning epoch data output by the airborne GNSS receiver, which contains positioning coordinates, solution status identifiers, and PDOP values; body attitude angles output by the airborne attitude measurement unit and synchronized with the differential positioning epoch; exposure pulse times obtained by latching the exposure feedback pulses from the mapping camera; and a continuous visual frame sequence with a frame rate higher than the output frequency of the differential positioning epoch.

[0015] S2. Visual odometry front-end calculation and quality quantization: Feature points are extracted and tracked frame by frame in the continuous visual frame sequence. The relative rotation matrix and the relative translation vector with undetermined scale between adjacent frames are calculated. The relative rotation matrix and the relative translation vector are then converted to the body coordinate system based on the pre-calibrated fixed connection between the visual sensor and the body. At the same time, the visual quality factor is calculated from the number of tracked feature points, the average disparity, and the average reprojection residual.

[0016] S3. Online calibration of visual scale factor: Select a continuous epoch interval where the solution state is identified as a fixed solution and the visual quality factor is not less than a set quality threshold as the scale calibration interval. Use the GNSS baseline length between adjacent epochs in the scale calibration interval as the length reference, perform weighted least squares estimation on the magnitude of the relative translation vector to obtain the visual scale factor, and simultaneously calculate the relative standard error of the visual scale factor from the residual of the weighted least squares.

[0017] S4. Exposure Time Projection Center and Attitude Calculation: The exposure time obtained after applying shutter delay correction to the exposure pulse time is taken as the target time. The differential positioning epochs that enclose the target time are retrieved. The rotation matrix formed by the positioning coordinates of the epoch and the body attitude angle is taken as the starting value. The inter-frame pose increments in the interval from the epoch to the target time are accumulated and multiplied. After being restored to the physical scale by the visual scale factor, motion compensation is performed on the positioning coordinates of the epoch. At the same time, lever correction related to the exposure time attitude is applied to obtain the camera projection center coordinates and attitude of each exposure point.

[0018] S5. Adaptive weighting considering the evolution of lockout: Using the lockout duration calculated from the most recent fixed solution epoch as the independent variable, a variance model is constructed for the observations of each exposure point position. This variance model consists of three superimposed terms: the first term is a differential positioning term that increases with the solution state coefficient and PDOP value; the second term is a random walk term that is directly proportional to the lockout duration and amplified as the visual quality factor decreases; and the third term is a scale drift term that is proportional to the square of the product of the relative mean square error of the visual scale factor, the flight speed, and the lockout duration. The ratio of the baseline variance under the fixed solution state to this variance is used as the weight of the observation of the exposure point position.

[0019] S6. Image Control-Free Joint Bundle Adjustment: Without using ground image control points, the image point coordinates are the main observations, the camera projection center coordinates of the exposure points obtained in S4 are the absolute observations and participate in the adjustment according to the weights obtained in S5, the relative position increments and relative attitude increments derived from the inter-frame pose increments between adjacent exposure points are the relative observations, and the prior values ​​of the camera self-calibration additional parameters are the virtual observations. A unified weighted least squares objective function is constructed and iteratively solved to obtain the exterior orientation elements of each exposure point and the ground coordinates of the densification points.

[0020] S7. Accuracy Assessment and Output Production: Statistically measure the internal accuracy of densification points and the mean error of check points, and generate digital surface models, digital elevation models, digital orthophoto maps, and digital line maps based on the internal and external orientation elements obtained from the adjustment.

[0021] By adopting the above technical solution, the present invention has the following beneficial effects:

[0022] The core idea of ​​this invention is to use visual odometry as an observation method with quantitatively assessable accuracy, rather than merely as a backup when satellite lock-on is lost. The visual scale factor is calibrated online using a fixed baseline from differential positioning, simultaneously providing the relative mean square error of this scale factor. This relative mean square error, along with the lock-on duration and flight velocity, constitutes the scale drift term in the variance of the exposure point location observations. Thus, the relationship between satellite positioning and visual odometry is not one of data substitution or weighted averaging, but rather one where satellite positioning calibrates the visual scale, and the accuracy of the visual scale calibration, in turn, determines the weighting of the observations. This is what distinguishes this invention from existing technologies.

[0023] As a further aspect of the present invention: the visual quality factor in S2 and the visual scale factor in S3 are calculated according to formula (1).

[0024]

[0025] In the formula, η k N represents the visual quality factor corresponding to the k-th epoch; k d k r k The values ​​are, in order, the number of tracking feature points at that epoch, the average disparity in pixels, and the average reprojection residual in pixels; N0, d0, and r0 are, in order, the baseline values ​​for the number of tracking feature points, the average disparity, and the average reprojection residual; ŝ is the estimated visual scale factor; Ω is the scale calibration interval; ω k Δt is the weighting coefficient for the k-th epoch, which is the ratio of the visual quality factor to the PDOP value for that epoch; k,k+1 b X is a relative translation vector with an undetermined scale between the k-th epoch and the (k+1)-th epoch; G This provides the GNSS positioning coordinates for this epoch. The visual quality factor incorporates the number of features, baseline geometric strength, and matching residuals. The number of tracked feature points and average disparity are saturated, with no further gain beyond the baseline value to avoid inflated evaluations in highly textured regions. The average reprojection residual is negatively exponential, causing the quality factor to decay rapidly with even slight increases in residuals, making it highly sensitive to mismatches and motion blur. The scale factor is estimated by weighting the visual quality factor to the PDOP value, with epochs exhibiting high geometric strength and reliable visual resolution dominating the calibration process.

[0026] As a further aspect of the present invention: the scale calibration interval in S3 is taken from a straight flight segment where the change in body attitude angle between adjacent epochs is less than a set attitude stability threshold, so that the difference between the displacement of the GNSS antenna phase center and the displacement of the camera projection center within this interval is negligible; the number of epochs within the scale calibration interval is not less than 20, and the corresponding cumulative track length is not less than 50 meters; after removing outliers from the scale residuals of each epoch within this interval according to the median absolute deviation criterion, the scale is iteratively resolved until the relative change of the visual scale factor is less than 1×10. −4 The relative mean square error of the visual scale factor is taken as the ratio of the mean square error of the visual scale factor to its estimated value, wherein the mean square error of the visual scale factor is obtained by the unit weight mean square error of the weighted least squares method and the weight coefficients corresponding to the visual scale factor according to the error propagation law. Specifying the minimum epoch number and minimum track length of the calibration interval is to ensure that the length reference itself has sufficient strength: when the track is too short, the relative error of the GNSS baseline length is amplified, and the scale factor used for calibration will carry unacceptable systematic bias.

[0027] As a further aspect of the present invention: the camera projection center coordinates and attitude of the exposure point mentioned in S4 are calculated according to equation (2).

[0028]

[0029] In the formula, t e For the exposure moment; t p δ is the exposure pulse time; t is the shutter lag correction; k To enclose the differential positioning epoch time of this exposure time; X S X represents the coordinates of the camera projection center in the mapping coordinate system. G R represents the coordinates of the GNSS antenna phase center in the mapping coordinate system. b w The rotation matrix is ​​Δt from the body coordinate system to the mapping coordinate system. k,e b With ΔR k,e b The relative translation vector and relative rotation matrix accumulated within the interval are, in order; bThis is the representation of the lever vector pointing from the camera projection center to the GNSS antenna phase center in the body coordinate system. The formula first extrapolates the coordinates of the differential positioning epoch to the exposure time using the pose increment accumulated by the visual odometry within the envelope interval. Then, it uses the attitude at the exposure time to normalize the lever vector from the GNSS antenna phase center to the camera projection center. Since the frame rate of the visual odometry is much higher than the output frequency of the differential positioning epoch, this extrapolation is based on the actual pose change accumulated frame by frame within the envelope interval, rather than a first-order extrapolation from the velocity vector at a single moment. Therefore, it can realistically depict the maneuvering process within the interval. Compared to the approach of treating the exposure delay as a system error parameter to be absorbed during adjustment, this scheme eliminates this error during the observation generation stage, avoiding its correlation with exterior azimuth elements and principal point offset. Furthermore, since the attitude is also synchronously derived from the visual pose increment, the lever correction and motion compensation are strictly self-consistent at the same time.

[0030] As a further aspect of the present invention: In S4, the camera projection center coordinates obtained by solving according to equation (2) are subtracted from the camera projection center coordinates obtained by linear interpolation of adjacent epochs. When the difference is greater than the set limit, bidirectional interpolation is performed on both sides of the envelope epoch with the visual quality factor as the weight, and the weighted mean is taken for recalculation. When the difference between the two calculation results is still greater than the set limit, the observation value of the exposure point position is reduced to a weakly constrained observation value. Further, after the first adjustment convergence in S6, the exterior azimuth angle elements of each exposure point obtained by adjustment are used to replace the body attitude angles collected in S1. The exposure point coordinates are recalculated in S4 and then adjusted again. The iteration continues until the difference in the error of the check point between two adjacent adjustments is less than the set convergence limit. The significance of this back-substitution iteration lies in the following: the motion compensation amount of S4 needs to be rotated to the mapping coordinate system based on the attitude of the starting epoch. The angular error of the airborne attitude measurement unit will be proportionally transmitted to the exposure point coordinates according to the length of the compensation amount. Since the compensation interval is shorter than one positioning epoch and the compensation amount usually does not exceed the product of the flight speed and the epoch interval, the error is initially at the centimeter level, and can be reduced to the millimeter level after one back-substitution. Thus, the accuracy requirement of this method for the airborne attitude measurement unit is reduced to the flight control level. Therefore, a self-calibration mechanism is established within the solution process, which can automatically identify periods of intense maneuvering and visual degradation without relying on external checkpoints.

[0031] As a further aspect of the present invention: the variance model described in S5 is Equation (3),

[0032]

[0033] In the formula, σ k κ represents the standard error of the observation value of the k-th exposure point position; k τ is the solution state coefficient; σ0 is the reference standard error under the fixed solution state; PDOP0 is the reference PDOP value; τ kη is the duration of lockout calculated from the most recent fixed solution epoch; k ε is the visual quality factor. s denoted as the relative mean square error of the visual scale factor; v is the flight speed; α is the visual random walk coefficient; β is the scale drift coefficient; p k This represents the weight of the observed value at the exposure point. This variance model is the core of this invention. The first term is the differential positioning term, which characterizes the accuracy level of differential positioning itself; the second term is the random walk term, which is directly proportional to the duration of the loss of lock. This is because the pose recursion of the visual odometry during the loss of lock is a frame-by-frame cumulative random process, and its cumulative variance is proportional to the number of recursion steps, i.e., proportional to the duration. Moreover, this term is amplified as the visual quality factor decreases, thus incorporating the visual working state into the random model; the third term is the scale drift term, which is proportional to the square of the product of the scale relative error, flight speed, and duration of the loss of lock. This is because the relative error of the scale factor exhibits a systematic deviation that is linearly amplified along the flight path when recursively applied along the flight path. Since the flight path is equal to the product of the flight speed and the duration of the loss of lock, this term is quadratic with respect to the duration of the loss of lock. After the three terms are superimposed, the weight of the exposure point location observation decreases continuously, monotonically, and physically based as the loss of lock occurs and continues, overcoming the inherent weight step of weighting according to the solution state. At the same time, since the random error term and the systematic error term have different orders with respect to the duration of the loss of lock, they dominate in different stages of the loss of lock, so that the variance model can give appropriate weights in both short-term and long-term loss of lock cases.

[0034] As a further aspect of the present invention: the solution state coefficients are taken as 1, 3, 10, and 30 respectively when the solution state is identified as a fixed solution, a floating-point solution, a pseudo-range difference decomposition, and a single-point solution; the visual random walk coefficients are taken as 1×10. −4 Up to 1×10 −3 The scale drift coefficient is between 0.5 and 2.0, measured in square meters per second. The values ​​of these coefficients can be determined by least-squares fitting of the calibration field flight data, or by taking empirical values ​​based on equipment specifications.

[0035] As a further aspect of this invention: In S5, the flight path is divided into strongly constrained segments and weakly constrained segments based on the continuity of the solution state and the visual quality factor. Each end of any weakly constrained segment contains at least three consecutive fixed solution exposure points, providing an absolute reference for the strongly constrained segments at both ends of the weakly constrained segment to bridge the gap. When a weakly constrained segment does not meet this condition, the variance of the observed values ​​at the exposure point positions of that segment is amplified according to the maximum unlock duration of that weakly constrained segment, and a supplementary flight prompt for that segment is output. When the visual quality factor is less than a set visual degradation threshold, the time period is determined to be a visual degradation time period. No relative observation values ​​are formed within this time period, and instead, a cross-segment constraint is established using at least five strongly constrained exposure points before and after this time period. The function of the double-ended bridging is equivalent to transforming a branch line with only a starting-end constraint into a bounded line with both ends attached, which allows the drift amount of visual recursion to be effectively absorbed and assigned during adjustment.

[0036]

[0037] As a further aspect of the present invention: the unified weighted least squares objective function and the error equation of the relative observations mentioned in S6 are given by equation (4).

[0038] In the formula, Φ is the objective function; V p V G V R V c The residual vectors are, in order, the image point coordinate observations, the absolute exposure point position observations, the visual relative observations, and the prior virtual observations of the additional parameters; P p P G P R P c The weight matrices are, in order, the weight matrices corresponding to the four types of observations mentioned above; V R,ij X is the residual vector of the visual relative observations between the i-th exposure point and the j-th exposure point, with the upper half being the relative position residual and the lower half being the relative pose residual; S,i X S,j The coordinates of the camera projection centers for the two exposure points are, respectively; Δt ij b With ΔR ij b These are the cumulative relative translation vector and relative rotation matrix between the two exposure points, respectively; log(·) ∨ This represents the logarithmic mapping operation that maps a rotation matrix to a rotation vector. See also Figure 7After the above four types of observations are combined into a normal equation, the exterior orientation elements of each exposure point, the ground coordinates of the densification points, the interior orientation elements, and the self-calibration additional parameters are uniformly calculated. When the scale calibration interval is short or spans multiple flights, the correction value dŝ of the visual scale factor can also be listed as an unknown, and then substituted back into the relative observation value after the overall adjustment is solved. The upper part of the error equation is the relative position residual, which compares the difference between the projection center coordinates of adjacent exposure points with the visual relative translation vector after scale restoration; the lower part is the relative attitude residual, which maps the difference of the rotation matrix into a three-dimensional rotation vector through logarithmic mapping, thereby ensuring the rigor of the attitude residual on the manifold and the correctness of the linearization.

[0039] As a further aspect of the present invention, a UAV field mapping system integrating real-time differential positioning and visual odometry is provided, comprising an onboard acquisition unit and a ground data processing unit. The onboard acquisition unit includes: a GNSS differential positioning module for receiving differential correction data and outputting differential positioning epoch data containing positioning coordinates, solution status identifiers, and PDOP values; an attitude measurement module for outputting the aircraft's attitude angle synchronized with the differential positioning epoch; a synchronization time base module for establishing a unified time base using satellite second pulses, triggering camera exposure, and latching the exposure feedback pulse time; a mapping camera with a feedback interface to send exposure pulses back to the synchronization time base module; a visual sensor with a frame rate higher than the output frequency of the differential positioning epoch and rigidly connected to the mapping camera; and an onboard storage module for recording four observation data streams under the same unified time base. The ground data processing unit includes: a visual odometry calculation module, used to calculate the inter-frame relative rotation matrix and the scale-determined relative translation vector, reduce them to the body coordinate system, and calculate the visual quality factor; a scale factor calibration module, used to estimate the visual scale factor based on the GNSS baseline length of a fixed solution epoch interval and simultaneously calculate its relative mean square error; an exposure point calculation module, used to perform visual motion compensation and attitude-related lever correction for exposure times to obtain the camera projection center coordinates and attitude of each exposure point; an adaptive weighting and flight strip segmentation module, used to construct a variance model and weight it by superimposing three terms: the lockout duration as the independent variable, the differential positioning term, the random walk term, and the scale drift term, and to divide it into strongly constrained segments and weakly constrained segments; an image control-free joint adjustment module, used to incorporate four types of observations—image point coordinate observations, absolute observations of exposure point positions, visual relative observations, and prior virtual observations of additional parameters—into a unified objective function for iterative calculation; and a results production module, used to generate digital surface models, digital elevation models, digital orthophoto maps, and digital line maps. Among the aforementioned modules, the scale relative error output by the scale factor calibration module is fed into the adaptive weighting and flight strip segmentation module for variance calculation. The input to the scale factor calibration module is taken from the fixed solution epoch of the GNSS differential positioning module. These three modules thus constitute the aforementioned accuracy closed loop. This system can be implemented on existing UAV aerial photography platforms by adding a visual sensor and upgrading the software, without requiring the replacement of the high-precision inertial measurement unit or the addition of a lidar payload.

[0040] Compared with the prior art, the beneficial effects of the present invention are:

[0041] The weighting model aligns with the actual evolution of the loss-of-lock error. The first term of the loss-of-lock duration corresponds to the random walk of visual recursion, while the second term corresponds to scale drift. Consequently, the weights of the exposure point location observations decrease continuously and monotonically with the continuation of the loss-of-lock, eliminating the step-like fluctuations seen in tiered weighting. In actual measurements, the mean square error of the exposure point location within a continuous loss-of-lock segment on the order of 30 seconds decreased from 0.4m to 0.6m in free recursion to within 0.08m, essentially eliminating local distortions and elevation steps in the regional network.

[0042] The position extrapolation error at the exposure time is eliminated during the observation generation stage. The compensation is based on the cumulative true pose change within the envelope interval, which can characterize maneuvers that linear interpolation and first-order velocity extrapolation cannot. Since this error is eliminated before adjustment, it will not be correlated with exterior azimuth elements or principal point offset, thus avoiding becoming a systematic error in the results. The measured time synchronization residual has been reduced from 3cm to tens of centimeters in existing schemes to less than 1cm, and it is not worsened by the intensity of motion.

[0043] The two observation methods create a mutual constraint on accuracy. When satellite positioning is good, scale calibration is accurate, and visual observations can maintain a high weight even during periods of loss of lock. When satellite positioning is poor, scale calibration accuracy decreases accordingly, and the weight of the observations decays rapidly, preventing the introduction of new systematic errors by extrapolating from an unreliable scale.

[0044] The implementation cost is low. Tasks such as the establishment and connection of control points can be completely eliminated, which usually accounts for more than half of the fieldwork. There is no need to install high-precision inertial measurement units or lidar on the airborne end. It is only necessary to use the existing exposure feedback interface of the surveying camera and add a conventional vision sensor. Existing models can be implemented through software upgrades. Attached Figure Description

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

[0046] Figure 1 This is a flowchart illustrating the overall process of the UAV field mapping method that integrates real-time differential positioning and visual odometry as described in this invention.

[0047] Figure 2 This is a schematic diagram of the relationship between the variable altitude flight path layout and ground coverage in this invention, where (a) is a planar layout diagram and (b) is an elevation profile and photographic beam diagram.

[0048] Figure 3 This is a schematic diagram illustrating the geometric relationship between the projection center motion compensation and the lever arm correction at the exposure time in this invention.

[0049] Figure 4 This is a timing diagram of the synchronous acquisition of airborne multi-source data in this invention.

[0050] Figure 5 This is a schematic diagram illustrating the division of the strong constraint segment and the weak constraint segment of the flight strip in this invention, as well as the double-end bridging constraint.

[0051] Figure 6 This is a curve showing the change of the weight of the exposure point location observation value as a function of the unlock duration, as described in this invention, and it is compared with the existing grading and weighting methods.

[0052] Figure 7 This is a diagram showing the relationship between observed values ​​and unknown constraints in the image-controlled joint adjustment of this invention.

[0053] Figure 8 This is a block diagram of the composition of the UAV field image-free control and mapping system that integrates real-time differential positioning and visual odometry as described in this invention.

[0054] In the diagram: 100, UAV airborne data acquisition unit; 110, GNSS differential positioning module; 120, synchronous time base module; 130, mapping camera; 140, visual sensor; 150, airborne storage module; 160, attitude measurement module; 200, ground data processing unit; 210, visual odometry calculation module; 220, scale factor calibration module; 230, exposure point calculation module; 240, adaptive weighting and flight strip segmentation module; 250, image control-free joint adjustment module; 260, results production module; 300, differential reference station or network RTK service. Detailed Implementation

[0055] To make the technical problems to be solved, the technical solutions, and the beneficial effects of the present invention clearer, the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the present invention and are not intended to limit the present invention.

[0056] Example 1

[0057] See Figure 1 This embodiment uses a 1:500 topographic mapping task in a suburban area of ​​a county town as an example to illustrate the complete implementation process of this method. The survey area is approximately 1.2 km². 2 The terrain is a mixture of flat land and hills, with a maximum elevation difference of 86m in the survey area. A 220kV high-voltage transmission corridor runs diagonally from southwest to northeast through the middle of the survey area. There are also two large groves of trees on both sides of the corridor. The differential correction data of the real-time dynamic differential positioning of the network is intermittently interrupted in the airspace above the corridor and the forest area.

[0058] The operational equipment used constitutes the UAV-borne data acquisition unit 100, the composition of which is described in [reference needed]. Figure 8 Specifically: one hexacopter UAV aerial photography platform with a maximum takeoff weight of 6.2 kg and an endurance of 42 min; a full-frame digital camera 130 with 45 megapixels, a pixel size of 4.4 μm, and a nominal lens focal length of 35 mm; the hot shoe interface on the fuselage constitutes the exposure feedback interface; a GNSS differential positioning module 110 is a dual-frequency multi-system receiver supporting BeiDou, GPS, and GLONASS, with a positioning epoch output frequency of 5 Hz, simultaneously outputting solution status indicators and PDOP values; its differential correction data is broadcast by the network RTK service 300; and a visual sensor 140... A global shutter industrial camera with a resolution of 1280×800 and a frame rate of 30Hz is rigidly connected to the surveying camera 130 via a rigid bracket. The relative pose between the two cameras is pre-calibrated in the indoor calibration field. The synchronization time base module 120 uses GNSS second pulses for time synchronization, and its time latching resolution is better than 0.1ms. The above-mentioned observation data streams, along with the aircraft attitude angles output by the attitude measurement module 160, total four streams, are all recorded by the onboard storage module 150. The attitude measurement module 160 is a microelectromechanical system attitude reference system integrated into the flight control system, with an attitude angle mean square error of 0.5°, and outputs data synchronously with the differential positioning epoch at 5Hz. Data processing is completed by the ground data processing unit 200.

[0059] S1. Simultaneous acquisition of airborne multi-source data. See also Figure 4 The synchronization time base module 120 establishes a unified time base using GNSS second pulses and triggers the mapping camera 130 to expose according to the flight path; the hot shoe feedback pulse of the mapping camera 130 is latched by the synchronization time base module 120 to obtain the exposure pulse time t. p The GNSS differential positioning module 110 outputs positioning epoch data at 5Hz, with each epoch containing positioning coordinates, solution status identifier, and PDOP value; the visual sensor 140 outputs visual frame sequences at 30Hz, along with frame timestamps under the same time base. All four observation data streams are written to the onboard storage module 150.

[0060] The flight path is laid out alternately at two high-altitude levels, see [reference]. Figure 2The flight altitude for odd-numbered flight strips was 120m, and for even-numbered flight strips it was 135m, with an interlayer height difference of 15m. The ground resolution corresponding to a flight altitude of 120m was 1.51cm, and that of a flight altitude of 135m was 1.70cm, both better than the ground resolution requirements for 1:500 mapping. The forward overlap was set to 80%, and the lateral overlap was set to 70%, with a corresponding forward photographic baseline of 16.5m and a flight strip spacing of 37.1m. It should be noted that using multi-layer flight paths at different flight altitudes to improve self-calibration geometry is a known technique in the art. This embodiment uses this flight path configuration only to provide a complete and feasible operational plan; the technical contribution of this invention does not lie in the flight path configuration itself. This aerial photography was completed in two sorties, with a total effective operation time of approximately 75 minutes, acquiring 1965 mapping images, approximately 22,500 positioning epochs, and approximately 135,000 visual frames. The fixed solution rate for the entire measurement area was 84.7%, with a total of 11 integer ambiguity loss events. The longest continuous loss of lock duration was 18.4 seconds, which occurred above the high-voltage transmission corridor.

[0061] S2. Visual odometry front-end calculation and quality quantification. The visual odometry calculation module 210 extracts feature points frame by frame from the continuous visual frame sequence and tracks them between adjacent frames. The relative rotation matrix and the relative translation vector of undetermined scale are calculated using the sliding window local bundle adjustment method. The magnitude of the relative translation vector is expressed in the internal units of the visual odometry. At the same time, the number of tracked feature points, the average disparity and the average reprojection residual are counted frame by frame, and the visual quality factor is calculated according to the first equation of (1).

[0062] In this embodiment, the baseline value N0 for the number of tracked feature points is 200, the baseline value d0 for the average disparity is 8 pixels, the baseline value r0 for the average reprojection residual is 1.0 pixels, and the visual quality threshold η0 is 0.60. In areas with good texture, such as built-up areas, roads, and farmland, the number of tracked feature points is approximately 420, the average disparity is approximately 11.3 pixels, and the average reprojection residual is approximately 0.28 pixels. Substituting these values ​​into equation (1) yields a visual quality factor of approximately 0.76. During periods of occlusion and motion blur in forest areas, the number of tracked feature points decreases to approximately 90, the average disparity decreases to approximately 3.2 pixels, and the average reprojection residual increases to approximately 0.85 pixels. Substituting these values ​​into equation (1) yields a visual quality factor of only approximately 0.08. The difference between the two conditions is nearly an order of magnitude, and this quality factor is sufficient to distinguish the quality of the visual solution.

[0063] S3. Online calibration of visual scale factor. The scale factor calibration module 220 traverses the entire measurement area to locate epochs, selecting continuous epoch intervals with a solution state identifier of fixed solution and a visual quality factor of not less than 0.60, and requiring that the attitude angle change between adjacent epochs be less than 0.5°. This constitutes the scale calibration interval. The requirement of attitude stability is to ensure that the difference between the displacement of the antenna phase center and the displacement of the camera projection center within the interval is negligible. In this embodiment, the boom arm length is 0.149m, and the difference in displacement caused by a 0.5° attitude change between adjacent epochs is only 1.3mm, which is negligible relative to the approximately 1.6m displacement between epochs within the interval. A total of 27 calibration intervals meeting the conditions were selected in this aerial survey, of which 19 qualified intervals had no less than 20 epochs and a cumulative track length of not less than 50m, totaling 4380 epochs and a cumulative track length of 7008m. Using the GNSS baseline length between adjacent epochs within each segment as the length benchmark, the magnitude of the relative translation vector is estimated using weighted least squares according to the second equation of (1). The weighting coefficient is the ratio of the visual quality factor to the PDOP value. The initial solution yields a visual scale factor of 2.4691. Subsequently, 58 epochs with scale residuals exceeding three times the median absolute deviation are removed according to the median absolute deviation criterion. After three iterations and resolving, convergence occurs, and the final visual scale factor ŝ is 2.4713. According to the error propagation law, the standard error of this scale factor is calculated as 4.4 × 10⁻⁶ using the unit weight standard error and the corresponding weight coefficient. −3 The ratio of this to the scale factor estimate is the relative standard error ε. s , equals 1.8 × 10 −3 .

[0064] S4. Calculation of the projection center and attitude at the exposure time. See also... Figure 3 The exposure point calculation module 230 performs the following processing on each mapping image. First, it reads the exposure pulse time t of the image. p A shutter retardation correction δ, obtained from the calibration field, is applied. In this embodiment, δ is 8.6 ms, with a calibration residual jitter of ±0.4 ms, corresponding to a positional uncertainty of approximately 3 mm at a speed of 8 m / s. Based on this, the exposure time t is obtained. e Then, the envelope t was retrieved. e Differential positioning epoch t k The positioning coordinates of this epoch are read, and the starting rotation matrix R is constructed using the attitude angles output by the attitude measurement module 160 at this epoch. b w (t k ). Then t k To t eThe pose increments of each frame within the interval are multiplied sequentially to obtain the cumulative relative rotation matrix and the cumulative relative translation vector with undetermined scale for that interval. The latter is then restored to the physical scale using ŝ obtained from S3. Finally, motion compensation is performed on the positioning coordinates of that epoch according to equation (2), and the lever arm vector l is adjusted based on the pose at the exposure time. b The phase center of the GNSS antenna is reduced to the projection center of the camera. In this embodiment, the three components of the boom vector are 0.021m, -0.008m, and -0.147m, respectively.

[0065] To perform self-calibration, the camera projection center coordinates obtained by solving according to equation (2) were subtracted from the camera projection center coordinates obtained by linear interpolation of adjacent epochs, and the difference value Δd was obtained, with a limit of 0.02m. Statistics show that Δd generally does not exceed 0.006m in straight flight segments; however, in the flight path turning segment and during gust disturbance periods, Δd reaches a maximum of 0.037m, with 63 images exceeding the limit, accounting for 3.2% of the total number of images. For these 63 images, bidirectional interpolation was performed on both sides of the envelope epoch using the visual quality factor as weight, and the weighted average was taken for recalculation. Among them, the difference between the two calculation results of 59 images was reduced to within the limit, and the remaining 4 images were downweighted to weakly constrained observations.

[0066] S5. Adaptive weighting that takes into account the evolution of lock loss. See also Figure 6 The variance of the observation values ​​at each exposure point position is calculated and weighted by the adaptive weighting and flight strip segmentation module 240 according to equation (3). In this embodiment, the reference mean square error σ0 under the fixed solution state is 0.02m, the reference PDOP value PDOP0 is 2.0, and the visual random walk coefficient α is 4×10 −4 m 2 / s, scale drift coefficient β is set to 1.0, and scale relative error ε s Take the 1.8 × 10 from S3. −3 Take the flight speed v as 8 m / s; solve for the state coefficient κ. k The values ​​are 1, 3, 10, and 30 for the fixed solution, floating-point solution, pseudo-distance difference decomposition, and single-point solution, respectively.

[0067] For example, for an exposure point with a fixed solution state, a PDOP value of 1.6, and a lockout duration of 0, the mean square error of its position observation calculated according to equation (3) is approximately 0.026m, and the weight is approximately 0.61; while for an exposure point with a floating-point solution state, a PDOP value of 3.4, a lockout duration of 12s, and a visual quality factor of 0.71, the difference positioning term in its variance is 1.40 × 10⁻⁶. −2 m 2 The random walk term is 6.76 × 10⁻⁶. −3 m 2 The scale drift term is 2.99 × 10⁻⁶. −2 m 2The mean square error after superimposing the three terms is approximately 0.225m, and the weight is approximately 0.0079. See also... Figure 6 The weights obtained by this invention continuously and monotonically decrease with the duration of the lockout. In the early stage of the lockout, the random walk term dominates, and in the later stage, the scale drift term dominates. The existing weight curves that are classified according to the solution state are stepped, with a step at the solution state switching point and remain constant within the same solution state. The difference between the two is particularly significant after the lockout duration exceeds 10 seconds.

[0068] See Figure 5 Subsequently, the flight path was divided into strongly constrained and weakly constrained segments based on the continuity of the solution state and visual quality factor. This aerial survey resulted in 12 strongly constrained segments and 11 weakly constrained segments. For each weakly constrained segment, it was checked whether each end contained at least three consecutive fixed solution exposure points: 10 segments met this condition, and double-ended bridging constraints were established for them; the remaining segment was located at the edge of the survey area and only one end met the condition. The variance of the observed exposure point position for this segment was amplified by a factor of four based on its maximum unlock duration, and a re-flight prompt for this segment was output.

[0069] S6. Image-controlled joint bundle adjustment of the zone network. See also: Figure 7 The unified weighted least squares objective function shown in equation (4) is constructed by the image control-free joint adjustment module 250 without using any ground control points. Among them, the image point coordinate observation values ​​are obtained by automatic matching of corresponding points between images. In this embodiment, a total of 4.12 million connection points are matched and linearized according to the collinearity condition equation to form the error equation; the absolute observation value of the exposure point position is taken from the camera projection center coordinates obtained in S4, and participates in the adjustment according to the weight obtained in S5; the visual relative observation value is composed according to the second equation of equation (4), and its weight is determined according to the average value of the segment length and the visual quality factor of the segment; the a priori virtual observation value of the additional parameter is taken from the camera factory calibration value, which is used to suppress parameter divergence in the early stage of iteration. The self-calibration additional parameters are initially set as principal distance, two components of image principal point coordinates, three terms of radial distortion coefficient and two terms of eccentric distortion coefficient, a total of 8. The normal equation is stored as a sparse matrix and solved by the preconditional conjugate gradient method.

[0070] During the adjustment process, significance and correlation tests were performed on the self-calibration additional parameters according to conventional practices in the field. Additional parameters with significance statistics less than the critical value or absolute values ​​of correlation coefficients with exterior orientation elements not less than the correlation threshold were removed. Gross errors were detected using a weighted iterative method. In this embodiment, the significance level was set to 0.05, the corresponding critical value was set to 1.96, and the correlation threshold was set to 0.85. After the first adjustment in this embodiment, the third-order radial distortion coefficient and the second eccentric distortion coefficient were removed. After reconstructing the normal equation, the second adjustment converged. The unit weight mean square error after convergence was 0.42 pixels. A total of more than 17,000 image point gross errors and 9 exposure point position gross errors were detected and processed.

[0071] After the initial adjustment convergence, the exterior azimuth elements of each exposure point obtained from the adjustment are used to replace the aircraft attitude angles output by the attitude measurement module 160. The process then returns to S4 to recalculate the exposure point coordinates and adjust again. In this embodiment, the attitude angle mean square error of the attitude measurement module 160 is 0.5°, while the compensation amount of S4 does not exceed the product of the flight speed and the epoch interval, i.e., 1.6m. Therefore, the initial exposure point coordinate error introduced by the initial attitude error is approximately 1.4cm. After one back-substitution, this error is reduced to 1.4mm, the mean square error of the checkpoint elevation is improved from 7.1cm to 6.7cm, and the difference between the mean square errors of the checkpoints in the second and third adjustments is less than 0.2cm, thus convergence is determined.

[0072] S7. Accuracy Assessment and Results Production. Twenty-four high-precision checkpoints were evenly distributed within the survey area. These checkpoints were used solely for accuracy assessment and did not participate in the adjustment. Their horizontal and vertical accuracy was better than ±0.01m. Statistical results showed that the horizontal mean square error of the densified points was ±3.8cm, the vertical mean square error was ±6.7cm, and the maximum vertical residual was 0.14m. In contrast, for the same batch of image data, with all other processing steps identical except for replacing the continuous variance model in S5 with the existing solution-based weighting method, the horizontal mean square error was ±6.2cm, and the vertical mean square error was ±23.4cm. Furthermore, the results exhibited a significant bowl-shaped deformation, with the central part of the survey area sinking approximately 0.21m relative to the surrounding areas. Since the difference between the two adjustments lies solely in the weighting method, this comparison demonstrates that the approximately 71% reduction in the vertical mean square error is indeed contributed by the continuous variance model of this invention.

[0073] The aforementioned accuracy meets the accuracy requirements for densification points in flat and hilly areas as specified in the 1:500 topographic map aerial photogrammetry documentation. Finally, the results production module 260 generates, in sequence, a digital surface model, a digital elevation model, a digital orthophoto map with a ground resolution of 0.02m, and a 1:500 digital line map based on the internal and external orientation elements obtained from the adjustment.

[0074] Example 2

[0075] This embodiment focuses on a scenario with severe differential signal obstruction in a densely built-up urban area, illustrating the role of adaptive weighting and double-ended bridging constraints in the proposed method. The test area is a central business district of a certain city, covering an area of ​​approximately 0.86 km². 2 The area is densely populated with high-rise buildings, the tallest of which is about 160 meters high, forming a typical urban canyon between the buildings; the real-time dynamic differential positioning of the network frequently loses lock on the flight path between the buildings. The equipment used is the same as in Example 1.

[0076] In S1, considering the height of the tallest building and safety clearance, the two flight altitudes were set at 180m and 202m respectively. The forward overlap was set at 85%, and the lateral overlap at 75%, with a corresponding forward photogrammetric baseline of 18.5m. A total of 1001 images were acquired during this aerial survey, with 27 instances of image loss. The fixed resolution was only 61.3%, with the longest continuous image loss lasting 43 seconds, occurring between three supertall buildings on the southeast side of the survey area. The flight speed for this sortie was set at 9 m / s.

[0077] The implementation methods for S2 to S4 are the same as in Example 1. The built-up area has rich texture and clear geometric structure; the average visual quality factor across the entire measurement area reaches 0.81, indicating good working condition of the visual odometry; the relative mean error ε of the visual scale factor obtained in S3... s 1.5×10 −3 .

[0078] In S5, the longest 43s weakly constrained segment contains 21 exposure points. The solution states sequentially go through floating-point solutions, single-point solutions, and then back to floating-point solutions. Taking an exposure point in the middle of this segment with a floating-point solution state, a PDOP value of 3.6, a lockout duration of 21.5s, and a visual quality factor of 0.81 as an example, the variance of this point is calculated according to equation (3) to have a difference positioning term of 1.53 × 10⁻⁶. −2 m 2 The random walk term is 1.06 × 10⁻⁶. −2 m 2 The scale drift term is 8.42 × 10⁻⁶. −2 m 2 The mean square error after superimposing the three terms is 0.33m, and the weight is 3.6×10. −3 The weight of the exposure point relative to the fixed solution is reduced to about 1 / 170; if the existing grading and weighting scheme is used, the weight of the floating point solution is reduced to only 1 / 9 of that of the fixed solution, which will directly introduce the decimeter-level positioning error into the regional network.

[0079] This embodiment further illustrates the role of the aforementioned accuracy closed-loop. If, during the same period of lockout, the relative mean square error of the visual scale factor deteriorates to 5.0 × 10⁻⁶ due to insufficient calibration interval... −3 Then, according to equation (3), the scale drift term at that exposure point will be 8.42 × 10 −2 m 2 Increased to 9.36 × 10 −1 m 2 The mean square error increased from 0.33m to 0.99m, and the weights correspondingly increased from 3.6×10 −3 Reduced to 4.0×10 −4Once the accuracy of the scale calibration deteriorates, the weight of the observations decreases accordingly without human intervention, thus avoiding the introduction of new systematic errors into the regional network through unreliable scale extrapolation; this automatic weight reduction effect is not available in existing technologies that only weight according to solution state.

[0080] See Figure 5 The weakly constrained segment has 5 and 7 consecutive fixed exposure points at its two ends, respectively, satisfying the condition of at least 3. Therefore, a double-ended bridging constraint is established for it. After verification, the visual odometry drift at the end of the free recursion of this segment is 0.58m. After incorporating it into the joint adjustment according to formula (4) and using the strong constrained segments at both ends to provide the absolute reference, the mean error of the camera projection center coordinates of each exposure point in this segment is 7.6cm, which is about 87% better than the free recursion.

[0081] The implementation methods for S6 and S7 are the same as in Example 1. The accuracy was ultimately determined using 18 checkpoints, with a plane mean square error of ±5.1 cm and an elevation mean square error of ±8.9 cm for the densified points. In contrast, when the same batch of data was adjusted using an existing grading and weighting scheme, significant regional network distortion was observed in the area containing the 43-second weakly constrained segment, with a maximum elevation deviation of 0.78 m at the checkpoints, failing to meet mapping requirements. Even in scenarios with degraded differential signals, such as densely built-up urban areas, this method can still provide image control-free results that meet mapping accuracy requirements of 1:1000 to 1:500 scales.

[0082] Example 3

[0083] This embodiment addresses the scenario of a strip-shaped survey area where weak geometry and texture coexist, illustrating the suppression effect of visual relative observations on flight strip distortion and the handling of visual degradation periods in this method. The survey area is a 110kV transmission line corridor with a length of 18.6km and a survey bandwidth of 300m, passing through a reservoir, mudflats, and a large bamboo forest in sequence. The equipment used is the same as in Embodiment 1.

[0084] The geometric conditions of the strip-shaped survey area are inherently weak: limited by the corridor width, only one main flight strip and two secondary flight strips can be deployed, with a lateral overlap of only about 40%. The regional network lacks lateral constraints in the longitudinal direction, making it prone to flight strip bending and twisting, with errors accumulating along the length of the flight strips. In S1, a stepped elevation flight path design is adopted along the flight direction, that is, the flight altitude is switched between 150m and 170m every about 2km. The ground resolution corresponding to a flight altitude of 150m is 1.89cm, and the ground sheet of a single image is 154m×103m; the flight direction overlap is taken as 80%, the corresponding flight direction photographic baseline is 20.6m, the spacing between the three flight strips is 92.7m, and the total coverage width is 339m, which can completely cover the 300m mapping bandwidth; based on the two baselines at each end of the flight path extending beyond the survey area, a total of 2722 images are acquired along the entire line.

[0085] In S2, in the reservoir water surface section (approximately 1.2 km long, with approximately 58 images per navigation strip) and some mudflat sections, the images exhibit typical weak and repetitive texture features. The number of tracked feature points once dropped below 30, and the average reprojection residual rose above 1.6 pixels. The visual quality factor calculated according to formula (1) was below 0.05. In this embodiment, the visual degradation threshold is set to 0.15, and the above-mentioned time period is determined to be the visual degradation period.

[0086] In S5, visual relative observations are not included for periods deemed to be visually degraded to avoid introducing drift errors into the adjustment; instead, cross-segment constraints are established using at least 5 strongly constrained exposure points before and after the period, and the variance of the exposure point position observations in that segment is amplified according to the maximum unlock duration of that period. Since the water surface segment itself does not contain any land features that need to be mapped, this processing method does not affect the completeness of the final results. In the bamboo forest segment, although the visual quality factor drops to between 0.31 and 0.44, it is still higher than the visual degradation threshold, so the random walk term in equation (3) is calculated according to η. k The reciprocal of η is automatically amplified after variance is added and then used in the adjustment; k Taking a value of 0.35 and a lockout duration of 8 seconds as an example, the variance of this term is 9.14 × 10⁻⁶. −3 m 2 , approximately η k Taking 0.81 is 2.3 times that under the same conditions, which reflects the adaptive ability of the variance model to the visual working state.

[0087] In S6, visual relative observations play a crucial role in stiffening the strip network: the relative positions and attitudes between adjacent exposure points are constrained with centimeter-level and angle-level precision, which is equivalent to laying a continuous high-precision guide along the flight strip direction, thereby suppressing the cumulative bending and twisting of the flight strip in the longitudinal direction.

[0088] In S7, 32 checkpoints were evenly distributed along the route for accuracy assessment. The results were: horizontal mean square error ±4.6cm, and vertical mean square error ±7.8cm; the vertical deviation between the two ends and the middle of the flight strip was reduced from 0.42m in the existing grading and weighting scheme to 0.06m. This embodiment shows that this method is also applicable to various linear engineering survey areas such as highways, railways, power transmission lines, and oil and gas pipelines, and does not require the establishment of any image control points along the route, which can reduce the fieldwork costs and operational risks of linear engineering surveying.

[0089] Example 4

[0090] This embodiment illustrates the composition and implementation of the image-free control mapping system of the present invention. (See also...) Figure 8 The system includes an onboard data acquisition unit 100 and a ground data processing unit 200, and is provided with differential correction data by a differential reference station or network RTK service 300.

[0091] The UAV-borne data acquisition unit 100 is mounted on the UAV aerial photography platform. The GNSS differential positioning module 110 receives differential correction data broadcast by the differential reference station or network RTK service 300, calculates and outputs differential positioning epoch data containing positioning coordinates, solution status identifiers and PDOP values; the synchronization time base module 120 synchronously outputs the body attitude angle with the attitude measurement module 160 and the GNSS differential positioning module 110, providing the starting rotation matrix for motion compensation in S4; the synchronization time base module 120 establishes a unified time base based on the second pulse output by the GNSS differential positioning module 110, one path sends an exposure trigger command to the mapping camera 130 according to the flight path plan, and the other path latches the exposure pulse time fed back by the mapping camera 130 through the exposure feedback interface; the hot shoe interface of the mapping camera 130 constitutes the exposure feedback interface; the visual sensor 140 outputs a visual frame sequence at a frame rate higher than the differential positioning epoch output frequency, and is rigidly connected to the mapping camera 130 via a rigid bracket; the airborne storage module 150 records the above four observation data streams under the same unified time base and adds their respective time stamps.

[0092] The ground data processing unit 200 is a group of software modules deployed on the graphics workstation. Among them, the visual odometry calculation module 210 executes S2, outputting the inter-frame pose increment and visual quality factor; the scale factor calibration module 220 executes S3, outputting the visual scale factor and its relative mean square error; the exposure point calculation module 230 executes S4, outputting the camera projection center coordinates and attitude of each exposure point; the adaptive weighting and flight strip segmentation module 240 executes S5, outputting the weights of the observation values ​​of each exposure point position and the division results of strong and weak constraint segments; the image control-free joint adjustment module 250 executes S6, constructing the normal equations and iteratively solving them; and the results production module 260 executes S7, outputting the densification point results and various surveying and mapping products.

[0093] It should be noted that the scale relative error output by the scale factor calibration module 220 is not only used for quality reporting, but is also fed into the adaptive weighting and flight strip segmentation module 240 as a coefficient of the scale drift term in equation (3) to participate in the variance calculation; and the input of the scale factor calibration module 220 comes from the fixed solution epoch of the GNSS differential positioning module 110. These three components thus form the aforementioned accuracy closed loop, which is the essence that distinguishes this system from existing systems.

[0094] The airborne part of this system only requires adding a visual sensor 140 to the existing UAV aerial photography platform and exporting the existing exposure feedback signal of the mapping camera 130 to the synchronization time base module 120. The ground part is implemented in software, without the need to replace the high-precision inertial measurement unit or add a lidar payload.

[0095] It should also be noted that the relative flight altitude, inter-layer elevation difference, various thresholds and coefficient values ​​in the above embodiments are all exemplary values ​​given for ease of explanation. Those skilled in the art can make adaptive adjustments to the above values ​​according to the map scale, terrain type, equipment indicators and operating environment, based on equations (1) to (4), and such adjustments all fall within the protection scope of this invention.

[0096] It should be noted that, in this document, the terms "comprising," "including," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, article, or apparatus. Unless otherwise specified, an element defined by the phrase "comprising one..." does not exclude the presence of other identical elements in the process, method, article, or apparatus that includes that element.

[0097] The above are merely preferred embodiments of the present invention and do not limit the scope of the patent. Any equivalent structural or procedural transformations made based on the description and drawings of the present invention, or direct or indirect applications in other related technical fields, are similarly included within the scope of patent protection of the present invention.

Claims

1. A method for unmanned aerial vehicle (UAV) field mapping without image control, integrating real-time differential positioning and visual odometry, characterized in that, Includes the following steps: S1. Airborne Multi-Source Data Synchronous Acquisition: A unified time base is established using satellite second pulses. Under this unified time base, four observation data streams are synchronously acquired. The four observation data streams include: differential positioning epoch data output by the airborne GNSS receiver, which contains positioning coordinates, solution status identifiers, and PDOP values; body attitude angles output by the airborne attitude measurement unit and synchronized with the differential positioning epoch; exposure pulse times obtained by latching the exposure feedback pulses from the mapping camera; and a continuous visual frame sequence with a frame rate higher than the output frequency of the differential positioning epoch. S2. Visual odometry front-end calculation and quality quantization: Feature points are extracted and tracked frame by frame in the continuous visual frame sequence. The relative rotation matrix and the relative translation vector with undetermined scale between adjacent frames are calculated. The relative rotation matrix and the relative translation vector are then converted to the body coordinate system based on the pre-calibrated fixed connection between the visual sensor and the body. At the same time, the visual quality factor is calculated from the number of tracked feature points, the average disparity, and the average reprojection residual. S3. Online calibration of visual scale factor: Select a continuous epoch interval where the solution state is identified as a fixed solution and the visual quality factor is not less than a set quality threshold as the scale calibration interval. Use the GNSS baseline length between adjacent epochs in the scale calibration interval as the length reference, perform weighted least squares estimation on the magnitude of the relative translation vector to obtain the visual scale factor, and simultaneously calculate the relative standard error of the visual scale factor from the residual of the weighted least squares. S4. Exposure Time Projection Center and Attitude Calculation: The exposure time obtained after applying shutter delay correction to the exposure pulse time is taken as the target time. The differential positioning epochs that enclose the target time are retrieved. The rotation matrix formed by the positioning coordinates of the epoch and the body attitude angle is taken as the starting value. The inter-frame pose increments in the interval from the epoch to the target time are accumulated and multiplied. After being restored to the physical scale by the visual scale factor, motion compensation is performed on the positioning coordinates of the epoch. At the same time, lever correction related to the exposure time attitude is applied to obtain the camera projection center coordinates and attitude of each exposure point. S5. Adaptive weighting considering the evolution of lockout: Using the lockout duration calculated from the most recent fixed solution epoch as the independent variable, a variance model is constructed for the observations of each exposure point position. This variance model consists of three superimposed terms: the first term is a differential positioning term that increases with the solution state coefficient and PDOP value; the second term is a random walk term that is directly proportional to the lockout duration and amplified as the visual quality factor decreases; and the third term is a scale drift term that is proportional to the square of the product of the relative mean square error of the visual scale factor, the flight speed, and the lockout duration. The ratio of the baseline variance under the fixed solution state to this variance is used as the weight of the observation of the exposure point position. S6. Image Control-Free Joint Bundle Adjustment: Without using ground image control points, the image point coordinates are the main observations, the camera projection center coordinates of the exposure points obtained in S4 are the absolute observations and participate in the adjustment according to the weights obtained in S5, the relative position increments and relative attitude increments derived from the inter-frame pose increments between adjacent exposure points are the relative observations, and the prior values ​​of the camera self-calibration additional parameters are the virtual observations. A unified weighted least squares objective function is constructed and iteratively solved to obtain the exterior orientation elements of each exposure point and the ground coordinates of the densification points. S7. Accuracy Assessment and Output Production: Statistically measure the internal accuracy of densification points and the mean error of check points, and generate digital surface models, digital elevation models, digital orthophoto maps, and digital line maps based on the internal and external orientation elements obtained from the adjustment.

2. The UAV field mapping method without image control, which integrates real-time differential positioning and visual odometry according to claim 1, is characterized in that, The visual quality factor mentioned in S2 and the visual scale factor mentioned in S3 are calculated using the following formula: In the formula, η k N represents the visual quality factor corresponding to the k-th epoch; k d k r k The values ​​are, in order, the number of tracking feature points at that epoch, the average disparity in pixels, and the average reprojection residual in pixels; N0, d0, and r0 are, in order, the baseline values ​​for the number of tracking feature points, the baseline value for the average disparity, and the baseline value for the average reprojection residual; ŝ is the estimated visual scale factor; and Ω is the scale calibration interval. ω k The weighting coefficient is the ratio of the visual quality factor to the PDOP value in that epoch; Δt k,k+1 b X is a relative translation vector with an undetermined scale between the k-th epoch and the (k+1)-th epoch; G The GNSS positioning coordinates for this epoch.

3. The UAV field mapping method without image control, which integrates real-time differential positioning and visual odometry according to claim 1, is characterized in that... The scale calibration interval described in S3 is taken from a straight flight segment where the change in aircraft attitude angle between adjacent epochs is less than a set attitude stability threshold, so that the difference between the displacement of the GNSS antenna phase center and the displacement of the camera projection center within this interval is negligible; the number of epochs within the scale calibration interval is not less than 20, and the corresponding cumulative track length is not less than 50 meters; after removing outliers from the scale residuals of each epoch within this interval according to the median absolute deviation criterion, the scale is iteratively re-solved until the relative change of the visual scale factor is less than 1×10. −4 The relative mean square error of the visual scale factor is the ratio of the mean square error of the visual scale factor to its estimated value, wherein the mean square error of the visual scale factor is obtained by the unit weight mean square error of the weighted least squares and the weight coefficient corresponding to the visual scale factor according to the error propagation law.

4. The UAV field mapping method without image control, which integrates real-time differential positioning and visual odometry according to claim 1, is characterized in that... The camera projection center coordinates and orientation of the exposure point described in S4 are calculated using the following formula: In the formula, t e For the exposure moment; t p δ is the exposure pulse time; t is the shutter lag correction; k To enclose the differential positioning epoch time of this exposure time; X S X represents the coordinates of the camera projection center in the mapping coordinate system. G R represents the coordinates of the GNSS antenna phase center in the mapping coordinate system. b w The rotation matrix is ​​Δt from the body coordinate system to the mapping coordinate system. k,e b With ΔR k,e b The relative translation vector and relative rotation matrix accumulated within the interval are, in order; b This is the representation of the lever arm vector pointing from the camera projection center to the GNSS antenna phase center in the body coordinate system.

5. The UAV field mapping method without image control, which integrates real-time differential positioning and visual odometry according to claim 4, is characterized in that... In S4, the camera projection center coordinates calculated according to the formula described in claim 4 are subtracted from the camera projection center coordinates obtained by linear interpolation of adjacent epochs. When the difference is greater than the set limit, bidirectional interpolation is performed on both sides of the envelope epoch with the visual quality factor as the weight, and the weighted mean is taken to recalculate. When the difference between the two calculation results is still greater than the set limit, the observation value of the exposure point position is reduced to a weakly constrained observation value. Furthermore, after the adjustment converges for the first time in S6, the exterior azimuth angle elements of each exposure point obtained by adjustment are used to replace the body attitude angles collected in S1. The calculation of the camera projection center coordinates of each exposure point is returned to S4, and the adjustment is repeated in S6. The iteration continues until the difference in the error of the check point obtained by two adjacent adjustments is less than the set convergence limit.

6. The UAV field mapping method without image control, which integrates real-time differential positioning and visual odometry according to claim 1, is characterized in that, The variance model described in S5 is: In the formula, σ k κ represents the standard error of the observation value of the k-th exposure point position; k τ is the solution state coefficient; σ0 is the reference standard error under the fixed solution state; PDOP0 is the reference PDOP value; τ k This is the duration of the lockout calculated from the most recent fixed solution epoch; η k Visual quality factor; ε s denoted as the relative mean square error of the visual scale factor; v is the flight speed; α is the visual random walk coefficient; β is the scale drift coefficient; p k The weight of the observation value at the exposure point location.

7. The UAV field mapping method without image control, which integrates real-time differential positioning and visual odometry according to claim 6, is characterized in that... The solution state coefficients are taken as 1, 3, 10, and 30 respectively when the solution state is identified as a fixed solution, a floating-point solution, a pseudo-range difference decomposition, and a single-point solution; the visual random walk coefficient is taken as 1×10. −4 Up to 1×10 −3 The scale drift coefficient is between 0.5 and 2.0, and the area per square meter per second.

8. The UAV field mapping method without image control, which integrates real-time differential positioning and visual odometry according to claim 1, is characterized in that... S5 further divides the flight strip into strongly constrained segments and weakly constrained segments based on the continuity of the solution state and the visual quality factor. Each end of any weakly constrained segment contains no less than 3 consecutive fixed solution exposure points, so that the strong constrained segments at both ends of the weakly constrained segment provide an absolute reference for the double-end connection. When the weakly constrained segment does not meet this condition, the variance of the observed values ​​of the exposure point positions of the segment is amplified according to the maximum unlock duration of the weakly constrained segment, and a supplementary flight prompt for the segment is output. When the visual quality factor is less than the set visual degradation threshold, the time period is determined to be a visual degradation time period. No relative observation values ​​are formed within the time period, and cross-segment constraints are established with no less than 5 strong constrained exposure points before and after the time period.

9. The UAV field mapping method without image control, which integrates real-time differential positioning and visual odometry according to claim 1, is characterized in that, The unified weighted least squares objective function and the error equation for the relative observations described in S6 are as follows: In the formula, Φ is the objective function; V p V G V R V c The residual vectors are, in order, the image point coordinate observations, the absolute exposure point position observations, the visual relative observations, and the prior virtual observations of the additional parameters; P p P G P R P c The weight matrices are, in order, the four types of observations mentioned above; V R,ij X is the residual vector of the visual relative observations between the i-th exposure point and the j-th exposure point, with the upper half being the relative position residual and the lower half being the relative pose residual; S,i X S,j The coordinates of the camera projection centers for the two exposure points are, respectively; Δt ij b With ΔR ij b These are the cumulative relative translation vector and relative rotation matrix between the two exposure points, respectively; log(·) ∨ This represents the logarithmic mapping operation that maps a rotation matrix to a rotation vector.

10. A UAV field mapping system integrating real-time differential positioning and visual odometry, characterized in that, Includes an onboard data acquisition unit and a ground-based data processing unit for unmanned aerial vehicles (UAVs); The UAV-borne acquisition unit includes: a GNSS differential positioning module for receiving differential correction data and outputting differential positioning epoch data containing positioning coordinates, solution status identifiers, and PDOP values; an attitude measurement module for outputting the aircraft's attitude angle synchronized with the differential positioning epoch; a synchronization time base module for establishing a unified time base using satellite second pulses, triggering camera exposure, and latching the exposure feedback pulse time; a mapping camera with a feedback interface to send exposure pulses back to the synchronization time base module; a visual sensor with a frame rate higher than the output frequency of the differential positioning epoch and rigidly connected to the mapping camera; and an onboard storage module for recording four observation data streams under the same unified time base. The ground data processing unit includes: a visual odometry calculation module, used to calculate the inter-frame relative rotation matrix and the scale-determined relative translation vector, reduce them to the body coordinate system, and calculate the visual quality factor; a scale factor calibration module, used to estimate the visual scale factor based on the GNSS baseline length of a fixed solution epoch interval and simultaneously calculate its relative mean square error; an exposure point calculation module, used to perform visual motion compensation and attitude-related lever correction for exposure time; an adaptive weighting and flight strip segmentation module, used to construct a variance model and weight it by superimposing differential positioning, random walk, and scale drift terms with the unlock duration as the independent variable, and to divide it into strongly constrained and weakly constrained segments; an image control-free joint adjustment module, used to incorporate four types of observations into a unified objective function for iterative calculation; and a results production module, used to generate digital surface models, digital elevation models, digital orthophoto maps, and digital line maps; The system is used to execute the UAV field image-free mapping method that integrates real-time differential positioning and visual odometry as described in any one of claims 1 to 9.

Images