Airborne InSAR and LiDAR collaborative planning method and system for complex mountainous terrain mapping

Abstract

This invention relates to the field of collaborative operations in aerial remote sensing mapping and topographic mapping, and specifically discloses an airborne InSAR and LiDAR collaborative planning method and system for complex mountainous terrain mapping. By acquiring digital elevation models and sensor parameters and performing pre-mission observability simulations, this invention enables the pre-determination of overlapping areas, shadow areas, and areas with insufficient coverage, obtaining a quantitative basis for sensor task allocation, thus improving the targeting of mapping operations and the reliability of results. By dividing the survey area into InSAR priority areas, LiDAR priority areas, and transition areas according to observation suitability indicators, and merging them to generate LiDAR flight sub-survey areas, the collaborative relationship between InSAR undertaking the main coverage and LiDAR focusing on key areas is clarified, reducing the cost of using LiDAR for indiscriminate operations across the entire area. By generating InSAR and LiDAR flight paths separately, the inefficient mode of discovering gaps in the results after flight and then conducting supplementary flights is avoided, improving the operational efficiency and flight window utilization rate of aerial mapping in complex mountainous areas.

Description

Technical Field

[0001] This invention relates to the field of collaborative operation technology of aerial remote sensing mapping and topographic mapping, and in particular to an airborne InSAR and LiDAR collaborative planning method and system for complex mountain topographic mapping. Background Technology

[0002] Large-scale topographic mapping in complex mountainous areas typically faces constraints such as significant topographic relief, substantial relative elevation differences, well-developed ridges and valleys, difficult access, complex weather conditions, and short effective flight windows. In such areas, a single sensor solution often cannot simultaneously meet the requirements of rapid coverage, complete data, and detailed representation of local topography.

[0003] Airborne InSAR is an active side-looking radar imaging mapping technology, characterized by its all-weather, all-day coverage capability, high coverage efficiency, and strong mission mobility. Particularly in missions requiring rapid acquisition of large-scale terrain data, airborne InSAR can quickly generate main DEM or DSM results through strip-style side-looking observations. However, due to its distinct side-looking imaging geometry, it is prone to overlay, shadowing, perspective contraction, reduced coherence, and local unwrapping errors in high mountains, deep valleys, steep slopes, and incised ravines, leading to local elevation anomalies or voids in the data.

[0004] Airborne LiDAR acquires 3D point clouds by actively emitting laser pulses, offering advantages in finely depicting complex terrain, extracting ground points in vegetated areas, providing detailed observation of partially obscured areas, and representing high-resolution land surfaces. However, its limitations lie in its reliance on favorable weather conditions and airspace clearance. When undertaking large-scale rapid coverage missions in mountainous areas alone, it may be constrained by weather, flight path density, point cloud density control, and flight costs.

[0005] Therefore, in aerial mapping of complex mountainous areas, airborne InSAR and airborne LiDAR have a clear complementary relationship: InSAR is suitable for large-scale main coverage, while LiDAR is suitable for key data acquisition in areas with overlapping, shadows, complex slopes, and potential voids in the data. The key issue is not simply overlaying the two types of sensors, but rather transforming terrain conditions, observation geometry, and operational constraints into an executable sensor division of labor and flight path planning scheme before the mission is implemented.

[0006] In existing aerial mapping missions in complex mountainous areas, the common practice is to determine InSAR or LiDAR flight paths based on experience, or to arrange local LiDAR operations only after voids or anomalies are discovered during post-flight quality checks. This approach has the following shortcomings:

[0007] (1) Existing solutions lack quantitative identification of overlay, shadows and effective flight band shrinkage caused by airborne InSAR side-looking imaging before the mission, making it difficult to determine which areas are suitable for rapid InSAR coverage and which areas need to be included in the LiDAR operation area for key data acquisition before flight.

[0008] (2) Existing solutions often optimize InSAR or LiDAR routes separately, lacking a method to unify the applicable scope of the two types of sensors into the same task partitioning framework, making it difficult to form clear boundaries between InSAR priority area, transition area and LiDAR priority area.

[0009] (3) Some methods only identify restricted areas by grid or small patches, but fixed-wing aircraft are not suitable for frequent turning operations around scattered patches. If LiDAR routes are directly laid out according to scattered patches, it will lead to broken routes, increased invalid flight distances, complicated airspace approval and reduced operational efficiency.

[0010] (4) Existing mountain route planning usually follows the nominal route width of the plain area or sets a uniform heading and uniform altitude, ignoring the impact of local terrain main direction and elevation differences on InSAR observability, which can easily form continuous overlapping or shadow areas in local mountainous areas.

[0011] (5) For mountainous tasks with strict airspace approval and short weather windows, there is a lack of collaborative mapping task organization methods that take into account observation geometry, weather windows, airspace boundaries, maximum flight path length and flight platform characteristics. Summary of the Invention

[0012] To overcome the aforementioned shortcomings in existing aerial mapping of complex mountainous areas, this invention provides an airborne InSAR and LiDAR collaborative planning method and system for mapping complex mountainous terrain.

[0013] In a first aspect, the present invention provides an airborne InSAR and LiDAR collaborative planning method for mapping complex mountainous terrain, comprising: Acquire the digital elevation model of the area to be measured, as well as the sensor parameters of airborne InSAR and airborne LiDAR; Based on the digital elevation model and the side-looking imaging geometric parameters of the airborne InSAR, a pre-mission observability simulation is performed to determine the observation type of each pixel in the area to be measured. The observation types include normal area, overlapping area, shadow area and insufficient coverage area. Based on the observation type of each pixel, the observation suitability index is statistically analyzed, and the area to be measured is divided into InSAR priority area, LiDAR priority area and transition area according to the observation suitability index. The LiDAR priority area and the transition area are merged into a LiDAR candidate operating range, and the LiDAR candidate operating range is subjected to spatial post-processing to generate at least one LiDAR flight sub-survey area. InSAR flight paths are generated based on the InSAR priority area, and LiDAR flight paths are generated based on the LiDAR flight sub-survey area.

[0014] According to one specific implementation, the above planning method includes performing pre-task observability simulation, specifically comprising: Establish the aircraft's flight reference point based on the candidate heading and flight altitude; The pixel center of the digital elevation model is taken as the ground point to be determined; Based on the spatial relationship between the aircraft flight reference point and the ground point to be determined, calculate the lateral distance, elevation difference, and nadir angle of the ground point relative to the aircraft. If the ground point to be determined is located on the radar observation side, and the carrier aircraft is higher than the point, and the nadir angle is within the preset observation angle range, then the point is determined to have entered the effective observation domain. Within the effective observation domain, azimuth profiles are organized according to the distance along the flight path. Within each profile, the profiles are sorted by lateral distance, and overlapping pixels and shadow pixels are identified by using the monotonicity of slant distance and the line-of-sight occlusion relationship.

[0015] According to a specific implementation method, in the above planning method, the observation suitability index includes the percentage of overlapping area, the percentage of shaded area, the percentage of insufficiently covered area, the percentage of observation-restricted connected patches, and a comprehensive suitability index. The comprehensive suitability index Q is calculated using the following formula: Q = 0.35p L +0.45p S +0.20p E p L p represents the percentage of overlapping area. S p represents the proportion of the shaded area. E The percentage of areas not covered; Based on the observation suitability index, the area to be measured is divided into an InSAR priority area, a LiDAR priority area, and a transition area, specifically including: The area to be tested is divided into several grid units according to a preset spatial scale; For each grid cell, the distribution of observation types of all pixels within that grid cell is statistically analyzed, and a comprehensive suitability index is calculated. When p S <0.012, p L When Q < 0.035 and Q < 0.030, the grid cell is designated as an InSAR priority area; When p S≥0.040, or p L +p S ≥0.080, or Q≥0.050, or p C When p ≥ 0.100, the grid cell is designated as a LiDAR priority region, where p C To observe the proportion of restricted connected patch area; The remaining grid cells are designated as transition zones.

[0016] According to a specific implementation method, the above planning method involves spatial post-processing of the LiDAR candidate operation area to generate at least one LiDAR flight sub-survey area, specifically including: Merge LiDAR candidate operation patches that are spatially adjacent or whose distance is less than a preset threshold. Hole filling is performed on the InSAR priority area surrounded by the LiDAR candidate operation range; Remove isolated task fragments with an area smaller than a preset threshold; Based on the major axis direction of the smallest circumscribed rotating rectangle of the candidate operational area after merging, filling, and eliminating, the terrain elevation zone, or the airspace boundary, the continuous candidate operational areas are divided into one or more LiDAR flight sub-areas.

[0017] According to one specific implementation, the above planning method, which generates InSAR flight routes based on the InSAR priority area, further includes: The nominal ground cover width W is calculated based on the average elevation, maximum observed nadir angle, and minimum observed nadir angle of the InSAR sub-region. I The calculation formula is: W I =(H f - (tanθ) max -tanθ min ), where H f For the first designed flight altitude, θ represents the average elevation of the sub-survey area. max and θ min These are the maximum and minimum observed nadir angles, respectively. According to the overlapping area ratio p in the sub-test area L and the proportion of shadow area p S The nominal ground coverage width is corrected to obtain the effective flight strip width W. I,eff The calculation formula is: W I,eff =W I (1-p L -p S ); Based on the effective flight strip width and safety factor β I Determine the InSAR line spacing D IThe calculation formula is: D I =β I W I,eff The safety factor β I The value ranges from 0.85 to 0.95.

[0018] According to a specific implementation, the above planning method, which generates LiDAR flight paths based on the LiDAR flight sub-survey areas, further includes: For the k-th LiDAR flight sub-region, based on the highest point z in the digital elevation model within that sub-region... max,k and minimum safe relative flight altitude h min Determine the second design flight altitude H L,k The calculation formula is: H L,k =z max,k +h min ; Based on the LiDAR's field of view φ and lateral overlap rate o L Determine the conservative flight band width W L and the distance between the routes D L The calculation formula is: W L =2h AGL tan(φ / 2), D L =W L (1-o L ); The direction of the long side of the smallest circumscribed rotating rectangle of the LiDAR flight sub-region is taken as the direction of the LiDAR flight path, and the path is extended outward by a preset buffer distance after intersecting with the actual boundary.

[0019] According to one specific implementation, in the above planning method, after generating InSAR and LiDAR flight paths, the method further includes: Hole checks and accuracy checks are performed on the InSAR results obtained using the InSAR route and the LiDAR results obtained using the LiDAR route. If a void or elevation anomaly appears in the transition zone or the InSAR priority zone, the area will be added to the LiDAR operation range. Based on the post-flight inspection results, the threshold values ​​for the observation suitability index are corrected for subsequent mission parameter optimization.

[0020] According to one specific implementation, in the above planning method, after obtaining the digital elevation model and sensor parameters of the area to be measured, the method further includes: The digital elevation model is subjected to projection unification, outlier removal, and boundary clipping. The area to be measured is divided into multiple InSAR sub-regions according to topographic relief, airspace boundaries, or rectangular grids; For each InSAR sub-area, the main elevation increase direction is determined based on the gradient statistics of the digital elevation model within that sub-area, and the direction perpendicular to this main elevation increase direction is taken as the recommended heading.

[0021] According to a specific implementation, in the above planning method, the sensor parameters include: the observation side-looking direction, candidate heading angle, minimum nadir angle, maximum nadir angle, design relative flight altitude, and strip overlap requirements for airborne InSAR; the field of view, lateral overlap rate, minimum relative flight altitude, point cloud density requirements, and maximum allowable single route length for airborne LiDAR; and the minimum safe flight altitude, turning buffer distance, available airspace boundaries, and operating windows of the flight platform.

[0022] Secondly, this invention provides an airborne InSAR and LiDAR collaborative planning system for mapping complex mountainous terrain, comprising: The data acquisition module is used to acquire the digital elevation model of the area to be measured, as well as the sensor parameters of the airborne InSAR and airborne LiDAR. The observability analysis module is used to perform pre-mission observability simulation based on the digital elevation model and the side-looking imaging geometric parameters of the airborne InSAR, in order to determine the observation type of each pixel in the area to be measured. The observation type includes normal area, overlapping area, shadow area and insufficient coverage area. The task partitioning module is used to calculate the observation suitability index according to the observation type of each pixel, and divide the area to be measured into an InSAR priority area, a LiDAR priority area, and a transition area according to the observation suitability index. The job area generation module is used to merge the LiDAR priority area and the transition area into a LiDAR candidate job area, and to perform spatial post-processing on the LiDAR candidate job area to generate at least one LiDAR flight sub-survey area. The flight path generation module is used to generate InSAR flight paths based on the InSAR priority area and to generate LiDAR flight paths based on the LiDAR flight sub-survey area.

[0023] Compared with the prior art, the beneficial effects of the present invention are as follows: This invention, by acquiring digital elevation models and sensor parameters and performing pre-mission observability simulations, enables the pre-determination of overlapping areas, shadow areas, and areas with insufficient coverage, providing a quantitative basis for sensor task allocation and improving the targeting of surveying operations and the reliability of results quality. By dividing the survey area into InSAR priority areas, LiDAR priority areas, and transition areas based on observation suitability indicators, and merging them to generate LiDAR flight sub-survey areas, the collaborative relationship between InSAR undertaking the main coverage and LiDAR focusing on key areas is clarified, reducing the cost of using LiDAR for indiscriminate operations across the entire area. By generating InSAR and LiDAR flight paths separately, the inefficient mode of discovering gaps in the results after flight and then conducting supplementary flights is avoided, improving the operational efficiency of aerial surveying in complex mountainous areas and the utilization rate of flight windows. Attached Figure Description

[0024] Figure 1 This is a flowchart illustrating an airborne InSAR and LiDAR collaborative planning method for mapping complex mountainous terrain, provided as an embodiment of the present invention. Detailed Implementation

[0025] The present invention will now be described in further detail with reference to specific embodiments. However, this should not be construed as limiting the scope of the present invention to the following embodiments; all technologies implemented based on the content of the present invention fall within the scope of the present invention.

[0026] Unless otherwise specified, in the description of specific embodiments of the present invention, "several", "more than", or "a number of" represent at least two. The number can be any number, such as two, three, four, five, six, seven, eight, or nine, and can even exceed nine.

[0027] This invention is primarily applied in the field of aerial remote sensing and mapping, particularly for mapping complex mountainous terrain characterized by dramatic topographic relief, deep valleys, and rapid slope changes. In the description of this invention's embodiments, airborne InSAR refers to an interferometric synthetic aperture radar mounted on a flight platform, which transmits and receives radar waves via a side-looking approach, generating a digital elevation model (DEM) or digital surface model through interferometric processing of complex images acquired by two or more antennas. Airborne LiDAR refers to a lidar system mounted on a flight platform, which acquires high-density three-dimensional point cloud data by actively emitting laser pulses and receiving echoes. Overlay refers to the phenomenon where the slant range order of radar beams reaching ground targets is reversed compared to the actual spatial order of the ground targets, typically occurring on the side of a hillside facing the radar. Shadowing refers to terrain areas behind the radar beams, typically occurring on steep slopes or the back of mountains facing away from the radar. A digital elevation model (DEM) is ground elevation data stored in a regular grid format. The technical solution of this invention allows for observability analysis based on prior DEMs and sensor parameters before the mission, thereby optimizing sensor allocation and flight path planning.

[0028] The technical solution of this invention will be described in detail below with reference to a specific complex mountain terrain mapping task. In this embodiment, the survey area is a high mountain canyon region, covering an area of ​​approximately 500 square kilometers, with a relative elevation difference exceeding 1500 meters. It is planned to use a fixed-wing UAV equipped with InSAR and LiDAR systems respectively to carry out a large-scale terrain mapping task.

[0029] Please refer to Figure 1 This illustration shows a flowchart of an airborne InSAR and LiDAR collaborative planning method for mapping complex mountainous terrain, provided by an embodiment of the present invention. The method includes: Step S1: Obtain the digital elevation model of the area to be measured, as well as the sensor parameters of the airborne InSAR and airborne LiDAR.

[0030] Specifically, before performing the surveying task, the boundaries, prior DEM, and existing terrain or airspace constraint data of the area to be surveyed are acquired. For example, vector boundaries and high-resolution DEM data of the high mountain and canyon area are extracted from a Geographic Information System (GIS) database. The grid spacing of the DEM data is 5 meters. Simultaneously, airborne InSAR parameters are acquired, including: observation side-looking direction (e.g., set to right-side look), candidate heading angles (e.g., one candidate direction every 15° from 0° to 360°), minimum nadir angle (e.g., 30°), maximum nadir angle (e.g., 60°), design relative flight altitude (e.g., 3000 meters), and strip overlap requirements (e.g., lateral overlap of no less than 20%). Airborne LiDAR parameters are also acquired, including: field of view (e.g., 70°), lateral overlap rate (e.g., 30%), minimum relative flight altitude (e.g., 800 meters), point cloud density requirements (e.g., no less than 4 points per square meter), and maximum permissible single flight path length (e.g., limited by UAV communication range, set to 50 kilometers). Simultaneously, flight platform parameters are acquired, including: minimum safe flight altitude (e.g., no less than 500 meters above the highest point on the ground), turning buffer distance (e.g., 500 meters), applicable airspace boundaries (e.g., temporary flight airspace areas approved by the civil aviation administration), and operating windows (e.g., daily from 8:00 AM to 4:00 PM). Obtaining these parameters provides a complete data foundation for subsequent observability analysis and route planning, thereby improving the accuracy and feasibility of the planning results.

[0031] Optionally, after acquiring the above data, a preprocessing step for the survey area is also included. Specifically, the DEM undergoes projection unification, outlier removal, and boundary clipping. For example, the DEM's coordinate system is converted to the WGS84 projected coordinate system, and outlier grid points with elevation values ​​significantly exceeding the normal range (e.g., points with negative elevations or exceeding the actual highest peak by more than 200 meters) are removed. Then, the DEM is clipped according to the survey area boundaries. The area to be surveyed is divided into multiple InSAR sub-survey areas based on terrain undulation, airspace boundaries, or rectangular grids. For example, based on the slope changes calculated from the DEM, continuous areas with gentle slope changes and elevation differences less than 300 meters are divided into a sub-survey area. For each InSAR sub-survey area, the main direction of elevation increase is determined based on the gradient statistics of the DEM within that sub-survey area, and the direction perpendicular to this main direction of elevation increase is used as the recommended heading. Specifically, the slope direction is calculated for each grid point within the sub-survey area, and the direction with the highest frequency is statistically analyzed as the main direction of elevation increase. If the main direction of elevation increase is north-south (i.e., elevation increases from south to north), then an east-west heading is recommended. This ensures that the InSAR's side-view direction is roughly perpendicular to the terrain contour lines, thereby reducing the probability of overlay and shadowing. This terrain-feature-based heading recommendation method can adaptively adjust to the local terrain of different sub-areas, avoiding the localized geometric degradation caused by a uniform heading across the entire area.

[0032] Step S2: Based on the digital elevation model and the side-looking imaging geometric parameters of the airborne InSAR, perform pre-mission observability simulation to determine the observation type of each pixel in the area to be measured. The observation type includes normal area, overlapping area, shadow area and insufficient coverage area.

[0033] For example, for each candidate heading and altitude, a flight reference point for the aircraft is established along the flight path. Specifically, a set of flight paths is set within the sub-survey area, and an aircraft position point P is taken at regular intervals (e.g., 100 meters) along each flight path. f P f Planar coordinates (x) including the carrier aircraft ff ,y ff The center of each pixel in the DEM within the sub-area is taken as the ground point P to be determined. P contains the plane coordinates (x, y) and elevation value z of the pixel.

[0034] According to P f Given the relative spatial position of point P, calculate the planar offset of the ground point relative to the aircraft. Decompose this offset into a distance A along the heading (i.e., P). f The projected length of P in the flight direction and the lateral distance C perpendicular to the flight direction (i.e., P) f (The distance to P perpendicular to the flight path). Simultaneously calculate the elevation difference Δh = Hz between the aircraft and the ground point, the slant range R = sqrt(C² + Δh²), and the nadir angle θ = arctan(C / Δh).

[0035] A preliminary assessment of pixel observability is performed: if a pixel is located on the radar observation side (e.g., when looking from the right, C>0), and the aircraft is higher than the pixel (Δh>0), and the nadir angle θ is within a preset observation angle range (e.g., 30°≤θ≤60°), then the pixel is considered to have entered the effective observation domain. Otherwise, it is considered not to participate in the effective observation of the current flight path. Through the above processing, the effective observation range under different headings and flight altitudes can be obtained before the mission is carried out. The advantage of this simulation method is that it can predict the radar beam illumination of the ground without actual flight, providing input for the subsequent accurate identification of overlay and shadows.

[0036] Within the effective observation domain, effective pixels are organized into azimuth profiles based on their distance A along the flight path. Within each profile, pixels are sorted by lateral distance C from near to far (i.e., from directly below the flight path to furthest from the flight path). The monotonicity of slant range is used for discrimination: within the same azimuth profile, if two ground points P1 and P2 have the same slant range R, but P2's lateral distance C2 is greater than C1, it indicates that the far and near points are compressed into the same slant range cell, and P2 is determined to be an overlay pixel. Line-of-sight occlusion is used for discrimination: if P2's slant range R2 is greater than P1's slant range R1, but P2's nadir angle is less than P1's nadir angle, it indicates that P1 occludes P2, and P2 is determined to be a shadow pixel. If both overlay and occlusion characteristics exist simultaneously, it is determined to be a coexisting region. Finally, a unified coded mask is output for the normal region, overlay region, shadow region, and coexisting region. The coded value of each pixel in the mask represents the observation type of that pixel. Areas at the edges of the flight strip, partially uncovered areas, or gaps caused by flight strip contraction are further marked as insufficiently covered areas. This fine-grained, pixel-level discrimination allows for the acquisition of observation type labels for each pixel within the survey area, providing a data foundation for subsequent suitability statistics based on spatial units.

[0037] Step S3: Based on the observation type of each pixel, statistical observation suitability indicators are calculated, and the area to be measured is divided into InSAR priority area, LiDAR priority area and transition area according to the observation suitability indicators.

[0038] Specifically, after the processing in step S2, each pixel in the survey area is assigned an independent observation type label. Based on this, the area to be surveyed is divided into several grid units according to a preset spatial scale, and the pixel-level discrimination results are aggregated into unit-level suitability indicators. For example, the area to be surveyed is divided into several grid units using a 1 km × 1 km regular grid. For each grid unit, the distribution of observation types of all pixels within that grid unit is statistically analyzed, and the following calculations are performed using a uniform grid unit (e.g., a 1 km × 1 km grid) as the statistical object, calculating the overlap area ratio p. L (Ratio of overlapping and coexisting areas to the total area), percentage of shaded area p S (Percentage of shaded and coexisting areas), percentage of under-covered areas p E The proportion of connected patch area p under observation limitations C (That is, the proportion of the area of ​​connected patches in overlapping, shaded, and undercovered areas to the total area of ​​the survey area). The overall suitability index is defined as Q = 0.35p. L +0.45p S +0.20p E The weighting for shadowing is higher than that for overlaying because shadowed areas typically cannot yield any usable data using InSAR technology, while overlaying areas can sometimes still provide some elevation information. The weighting for insufficient coverage is used to reflect the risk of continuity issues at strip edges.

[0039] In this embodiment, the three-level partitioning rules are set as follows: when p S <0.012, p L When p < 0.035 and Q < 0.030, the corresponding grid cell is designated as an InSAR priority area. S ≥0.040, or p L +p S ≥0.080, or Q≥0.050, or p C When the threshold is ≥0.100, the corresponding grid cell is designated as a LiDAR priority zone. The remaining grid cells are designated as transition zones. These thresholds were adjusted based on quality feedback from historical mapping projects in the mountainous area and the sensor parameters for this mission. For example, in previous InSAR data for this region, when the shadowed area ratio exceeded 4%, the hole rate in the results increased significantly; therefore, the shadow threshold was set to 0.04. Through a three-level partitioning system, the survey area is divided into three types: InSAR priority zones are suitable for rapid and continuous coverage using airborne InSAR; LiDAR priority zones have a higher risk of observational limitations and require priority for LiDAR operations; and transition zones are included in post-flight quality checks. This partitioning method quantifies the division of sensor tasks, avoiding the uncertainty caused by experience-based decision-making.

[0040] Furthermore, for InSAR priority areas, airborne InSAR is prioritized for rapid and continuous coverage to acquire the main DEM or DSM results. For LiDAR priority areas, airborne LiDAR is prioritized for acquiring key data. For transitional areas, they are included in the post-flight quality inspection scope, and whether to add LiDAR operations is determined based on post-flight results porosity, elevation anomalies, phase unwrapping quality, or point cloud coverage. This collaborative approach clarifies the task boundaries of InSAR undertaking large-scale main coverage and LiDAR focusing on complex obstructed areas, thereby reducing the indiscriminate operation cost of using LiDAR across the entire region.

[0041] Step S4: Merge the LiDAR priority area and the transition area into a LiDAR candidate operating range, and perform spatial post-processing on the LiDAR candidate operating range to generate at least one LiDAR flight sub-survey area.

[0042] Specifically, the LiDAR priority area, transition area, and manually designated key data acquisition area are merged into a LiDAR candidate operation range. Since these areas may be small, scattered patches in space, and fixed-wing UAVs are not suitable for frequent turns, spatial post-processing is required.

[0043] Specifically, LiDAR candidate operational patches that are spatially adjacent or within a preset threshold (e.g., 200 meters) are merged. For example, if there is a narrow gap between two patches, but the gap width is less than 200 meters, the gap is filled and the two patches are merged into a single continuous area. Small InSAR priority areas surrounded by candidate operational areas are filled with holes, as the aircraft cannot skip these small areas during continuous flight. Isolated fragments with an area smaller than a preset threshold (e.g., 0.1 square kilometers) are removed, as these fragments are extremely inefficient to manage individually and can be merged into adjacent InSAR operational areas or processed through post-flight verification. Then, based on the major axis direction of the smallest circumscribed rotation rectangle of the candidate operational area, the terrain elevation zone, and the airspace boundary, the continuous candidate operational area is divided into several LiDAR flight sub-areas. For example, if the major axis of the smallest bounding rectangle of a continuous operating area is 30 degrees east of north, and the terrain elevation within this area is divided into a low-elevation zone (800-1200 meters) and a high-elevation zone (1200-1600 meters), then it can be divided into two sub-areas along the boundary with the gentler elevation change. The flight altitude within each sub-area can be optimized separately. Through this post-processing, the gridded restricted area is transformed into an executable continuous operating area, reducing unnecessary turns and detours by the flight platform and improving operational efficiency.

[0044] Step S5: Generate InSAR flight paths based on the InSAR priority area and generate LiDAR flight paths based on the LiDAR flight sub-survey area.

[0045] Specifically, for an InSAR sub-survey area, the nominal ground coverage width is first calculated. Let the first design flight altitude be H. f The average elevation of the sub-survey area is Maximum observed nadir angle θ max (e.g., 60°), minimum observed nadir angle θ min (e.g., 30°), then the nominal ground cover width W I =(H f - (tanθ) max -tanθ min For example, H f =3000 meters, =1000 meters, then H f - =2000 meters, tan60°≈1.732, tan30°≈0.577, difference 1.155, W I ≈2310 meters.

[0046] Considering the coverage shrinkage caused by overlay and shading, adjust the effective flight strip width W. I,eff =W I (1-pL -p S If the average p in this sub-region L =0.02, p S =0.01, then 1-p L -p S =0.97, W I,eff ≈2240 meters. Lane spacing: D I =β I W I,eff , where β I For continuous coverage safety factor, the value ranges from 0.85 to 0.95. For example, if we take 0.9, then D... I ≈2016 meters. The InSAR flight path generated in this way can ensure effective coverage even with a certain proportion of observational limitations, avoiding gaps in coverage between flight strips.

[0047] For a LiDAR flight sub-region, let z be the highest point of the DEM in the k-th sub-region. max,k (For example, 1600 meters), the minimum safe relative altitude is h. min (For example, 800 meters), then the second design flight altitude H L,k =z max,k +h min =2400 meters. If the LiDAR field of view φ=70°, the lateral overlap rate o L =30%, then the conservative flight band width W L =2h AGL tan(φ / 2)=2 800 0.7 = 1120 meters, D L =W L (1-o L ) = 784 meters. h AGL h is taken as the relative flight altitude used for strip width estimation. AGL =h min Conservative estimations are used to avoid overestimating the flight path width due to sudden changes in local elevation. The LiDAR flight path direction is preferentially determined by the longer side of the smallest circumscribed rotating rectangle corresponding to the sub-survey area. If this longer side is 30 degrees east of north, the flight path is set parallel to this direction. After the flight path intersects the actual boundary, a preset buffer distance (e.g., 500 meters) is added at both ends to allow the aircraft to enter a stable flight state earlier. If the length of a single flight path exceeds the maximum allowable length (e.g., 50 kilometers), it is segmented according to airspace boundaries or length thresholds. For example, if a flight path is 70 kilometers long, a segmentation point is set at 50 kilometers, dividing it into two independent flight paths. Flight paths generated using the above method comply with the operational constraints of fixed-wing platforms and can be directly used for flight mission planning.

[0048] Furthermore, after data acquisition, this embodiment of the invention also includes performing hole checks, elevation anomaly checks, DEM / DSM consistency checks, and field elevation checkpoint accuracy statistics on the InSAR and LiDAR results. For example, the DEM generated by InSAR is compared with the elevations of pre-deployed ground control points. If the elevation error in a certain area exceeds 5 meters, it is marked as an elevation anomaly area. If localized holes, abnormal elevations, or quality indicators that do not meet requirements appear in the transition zone or InSAR priority area, the area is added to the LiDAR operation range, and a supplementary flight is arranged. Simultaneously, the post-flight inspection results are fed back to p. L p S p E p C The threshold correction process for Q. For example, in this task, it was found that at p S The region with a threshold of 0.035 (less than the original threshold of 0.04) also showed many data holes, indicating that the original shadow threshold was too high. Therefore, in subsequent similar mountainous tasks, this threshold was lowered to 0.03. Through this feedback mechanism, this method can adaptively optimize according to different terrain types and sensor parameters, and has good engineering adaptability.

[0049] Based on the above technical solutions, this invention achieves the pre-determination of overlapping areas, shadow areas, and areas with insufficient coverage by acquiring digital elevation models and sensor parameters and performing pre-task observability simulations. This provides a quantitative basis for the division of sensor tasks, improving the targeting of surveying operations and the reliability of results quality. By dividing the survey area into InSAR priority areas, LiDAR priority areas, and transition areas according to observation suitability indicators, and merging them to generate LiDAR flight sub-survey areas, the collaborative relationship between InSAR undertaking the main coverage and LiDAR focusing on key areas is clarified, reducing the cost of using LiDAR for indiscriminate operations across the entire area. By generating InSAR and LiDAR flight paths separately, the inefficient mode of discovering gaps in the results after flight and then conducting supplementary flights is avoided, improving the operational efficiency of aerial surveying in complex mountainous areas and the utilization rate of flight windows.

[0050] The beneficial effects of this invention are as follows: (1) It can identify potential overlapping, shadow and insufficient coverage areas in complex mountainous areas before flight, and transform the observation geometry problem that is difficult to judge directly into spatial partitioning results.

[0051] (2) It can clearly define the mission boundaries of InSAR and LiDAR, enabling InSAR to undertake large-scale rapid main coverage, and enabling LiDAR to focus on complex obstruction and result void risk areas, thereby reducing the cost of LiDAR operation across the entire domain.

[0052] (3) By connecting patches, filling holes and dividing the fixed wing area, the gridded restricted area can be transformed into an executable LiDAR operation area, reducing frequent turns, detours and invalid flight paths.

[0053] (4) By determining the recommended InSAR heading and local flight altitude by sub-survey area, the local observation geometric deviation caused by uniform heading or uniform flight altitude throughout the area can be reduced.

[0054] (5) The threshold is updated by post-flight quality feedback, which enables the method to be corrected according to different terrain types, sensor parameters and scale requirements, and has good engineering adaptability.

[0055] On the other hand, this invention also provides an airborne InSAR and LiDAR collaborative planning system for mapping complex mountainous terrain. This system includes the following functional modules: a data acquisition module, an observability analysis module, a task partitioning module, a work area generation module, and a flight path generation module.

[0056] The data acquisition module is used to acquire the digital elevation model of the area to be measured, as well as the sensor parameters of airborne InSAR and airborne LiDAR. For example, this module can be a data interface that reads the DEM file and sensor configuration file from an external storage device or a cloud database, as described in the previous embodiments.

[0057] The observability analysis module is used to perform pre-mission observability simulation based on the digital elevation model and the side-looking imaging geometric parameters of the airborne InSAR, in order to determine the observation type of each pixel in the area to be measured. The observation type includes normal area, overlapping area, shadow area and insufficient coverage area.

[0058] The task partitioning module is used to calculate the observation suitability index according to the observation type of each pixel, and divide the area to be measured into an InSAR priority area, a LiDAR priority area, and a transition area according to the observation suitability index.

[0059] The job area generation module is used to merge the LiDAR priority area and the transition area into a LiDAR candidate job area, and to perform spatial post-processing on the LiDAR candidate job area to generate at least one LiDAR flight sub-survey area.

[0060] The flight path generation module is used to generate InSAR flight paths based on the InSAR priority area and to generate LiDAR flight paths based on the LiDAR flight sub-survey area.

[0061] The specific working principles and implementation details of each of the above modules have been described in detail in the aforementioned method embodiments, and will not be repeated here.

[0062] The modules described above can be deployed in a ground planning workstation. This workstation includes a processor, memory, input / output interfaces, and a display. The processor can be a central processing unit (CPU), graphics processing unit (GPU), or application-specific integrated circuit (ASIC). The memory stores computer program instructions, which, when executed by the processor, implement the various steps of the aforementioned method. The input / output interfaces are used to connect a keyboard, mouse, and external data sources. The display is used to show the zoning results and flight path planning maps. Optionally, the system can also be implemented using a cloud platform, where planners submit survey area data and parameters via a web interface, and the cloud server completes the calculations and distributes the flight path files to the flight mission management terminal.

[0063] Furthermore, the present invention also provides an electronic device, including a memory, a processor, and a computer program stored in the memory and executable on the processor. When the processor executes the program, it implements the technical solution described in any of the methods above. Exemplarily, this electronic device may be a portable ruggedized laptop computer for use by field surveying engineers.

[0064] The present invention also provides a computer-readable storage medium having a computer program stored thereon. When the program is executed by a processor, it implements the technical solution described in any of the above methods. Exemplarily, the storage medium may be a USB flash drive, a portable hard drive, a read-only memory (ROM), a random access memory (RAM), a magnetic disk, or an optical disk, etc. Engineers can load the program from the storage medium into a planning computer to execute the collaborative planning method of the present invention.

[0065] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A collaborative planning method for airborne InSAR and LiDAR in complex mountainous terrain mapping, characterized in that, include: Acquire the digital elevation model of the area to be measured, as well as the sensor parameters of airborne InSAR and airborne LiDAR; Based on the digital elevation model and the side-looking imaging geometric parameters of the airborne InSAR, a pre-mission observability simulation is performed to determine the observation type of each pixel in the area to be measured. The observation types include normal area, overlapping area, shadow area and insufficient coverage area. Based on the observation type of each pixel, the observation suitability index is statistically analyzed, and the area to be measured is divided into InSAR priority area, LiDAR priority area and transition area according to the observation suitability index. The LiDAR priority area and the transition area are merged into a LiDAR candidate operating range, and the LiDAR candidate operating range is subjected to spatial post-processing to generate at least one LiDAR flight sub-survey area. InSAR flight paths are generated based on the InSAR priority area, and LiDAR flight paths are generated based on the LiDAR flight sub-survey area. The observation suitability indicators include the percentage of overlapping area, the percentage of shaded area, the percentage of insufficiently covered area, the percentage of connected patches with limited observation area, and a comprehensive suitability indicator. The comprehensive suitability index Q is calculated using the formula: Q = 0.35p L +0.45p S +0.20p E p L p represents the percentage of overlapping area. S p represents the proportion of the shaded area. E The percentage of areas not covered; Based on the observation suitability index, the area to be measured is divided into an InSAR priority area, a LiDAR priority area, and a transition area, specifically including: The area to be tested is divided into several grid units according to a preset spatial scale; For each grid cell, the distribution of observation types of all pixels within that grid cell is statistically analyzed, and a comprehensive suitability index is calculated. When p S <0.012, p L When Q < 0.035 and Q < 0.030, the grid cell is designated as an InSAR priority area; When p S ≥0.040, or p L +p S ≥0.080, or Q≥0.050, or p C When p ≥ 0.100, the grid cell is designated as a LiDAR priority region, where p C To observe the proportion of restricted connected patch area; The remaining grid cells are designated as transition zones.

2. The airborne InSAR and LiDAR collaborative planning method for complex mountain terrain mapping according to claim 1, characterized in that, Pre-mission observability simulations include: Establish the aircraft's flight reference point based on the candidate heading and flight altitude; The pixel center of the digital elevation model is taken as the ground point to be determined; Based on the spatial relationship between the aircraft flight reference point and the ground point to be determined, calculate the lateral distance, elevation difference, and nadir angle of the ground point relative to the aircraft. If the ground point to be determined is located on the radar observation side, and the carrier aircraft is higher than the point, and the nadir angle is within the preset observation angle range, then the point is determined to have entered the effective observation domain. Within the effective observation domain, azimuth profiles are organized according to the distance along the flight path. Within each profile, the profiles are sorted by lateral distance, and overlapping pixels and shadow pixels are identified by using the monotonicity of slant distance and the line-of-sight occlusion relationship.

3. The airborne InSAR and LiDAR collaborative planning method for complex mountain terrain mapping according to claim 1, characterized in that, Spatial post-processing is performed on the LiDAR candidate operational area to generate at least one LiDAR flight sub-area, specifically including: Merge LiDAR candidate operation patches that are spatially adjacent or whose distance is less than a preset threshold. Hole filling is performed on the InSAR priority area surrounded by the LiDAR candidate operation range; Remove isolated task fragments with an area smaller than a preset threshold; Based on the major axis direction of the smallest circumscribed rotating rectangle of the candidate operational area after merging, filling, and eliminating, the terrain elevation zone, or the airspace boundary, the continuous candidate operational areas are divided into one or more LiDAR flight sub-areas.

4. The airborne InSAR and LiDAR collaborative planning method for complex mountain terrain mapping according to claim 3, characterized in that, Generating InSAR flight paths based on the InSAR priority area further includes: The nominal ground cover width W is calculated based on the average elevation, maximum observed nadir angle, and minimum observed nadir angle of the InSAR sub-region. I The calculation formula is: W I =(H f - (tanθ) max -tanθ min ), where H f For the first designed flight altitude, θ represents the average elevation of the sub-survey area. max and θ min These are the maximum and minimum observed nadir angles, respectively. According to the overlapping area ratio p in the sub-test area L and the proportion of shadow area p S The nominal ground coverage width is corrected to obtain the effective flight strip width W. I,eff The calculation formula is: W I,eff =W I (1-p L -p S ); Based on the effective flight strip width and safety factor β I Determine the InSAR flight path spacing D I The calculation formula is: D I =β I W I,eff The safety factor β I The value ranges from 0.85 to 0.

95.

5. The airborne InSAR and LiDAR collaborative planning method for complex mountain terrain mapping according to claim 4, characterized in that, Generating LiDAR flight paths based on the aforementioned LiDAR flight sub-areas further includes: For the k-th LiDAR flight sub-region, based on the highest point z in the digital elevation model within that sub-region... max,k and minimum safe relative flight altitude h min Determine the second design flight altitude H L,k The calculation formula is: H L,k =z max,k +h min ; Based on the LiDAR's field of view φ and lateral overlap rate o L Determine the conservative flight band width W L and the distance between the routes D L The calculation formula is: W L =2h AGL tan(φ / 2), D L =W L (1-o L ); The direction of the long side of the smallest circumscribed rotating rectangle of the LiDAR flight sub-region is taken as the direction of the LiDAR flight path, and the path is extended outward by a preset buffer distance after intersecting with the actual boundary.

6. The airborne InSAR and LiDAR collaborative planning method for complex mountain terrain mapping according to claim 1, characterized in that, After generating InSAR and LiDAR flight paths, the method further includes: Hole checks and accuracy checks are performed on the InSAR results obtained using the InSAR route and the LiDAR results obtained using the LiDAR route. If a void or elevation anomaly appears in the transition zone or the InSAR priority zone, the area will be added to the LiDAR operation range. Based on the post-flight inspection results, the threshold values ​​for the observation suitability index are corrected for subsequent mission parameter optimization.

7. The airborne InSAR and LiDAR collaborative planning method for complex mountain terrain mapping according to claim 1, characterized in that, After acquiring the digital elevation model and sensor parameters of the area to be measured, the method further includes: The digital elevation model is subjected to projection unification, outlier removal, and boundary clipping. The area to be measured is divided into multiple InSAR sub-regions according to topographic relief, airspace boundaries, or rectangular grids; For each InSAR sub-area, the main elevation increase direction is determined based on the gradient statistics of the digital elevation model within that sub-area, and the direction perpendicular to this main elevation increase direction is taken as the recommended heading.

8. The airborne InSAR and LiDAR collaborative planning method for complex mountain terrain mapping according to claim 1, characterized in that, The sensor parameters include: the observation side-looking direction, candidate heading angle, minimum nadir angle, maximum nadir angle, design relative flight altitude, and strip overlap requirements for airborne InSAR; the field of view, lateral overlap rate, minimum relative flight altitude, point cloud density requirements, and maximum allowable single flight path length for airborne LiDAR; and the minimum safe flight altitude, turning buffer distance, available airspace boundaries, and operating windows for the flight platform.

9. An airborne InSAR and LiDAR collaborative planning system for mapping complex mountainous terrain, characterized in that, include: The data acquisition module is used to acquire the digital elevation model of the area to be measured, as well as the sensor parameters of the airborne InSAR and airborne LiDAR. The observability analysis module is used to perform pre-mission observability simulation based on the digital elevation model and the side-looking imaging geometric parameters of the airborne InSAR, in order to determine the observation type of each pixel in the area to be measured. The observation type includes normal area, overlapping area, shadow area and insufficient coverage area. The task partitioning module is used to calculate the observation suitability index according to the observation type of each pixel, and divide the area to be measured into an InSAR priority area, a LiDAR priority area, and a transition area according to the observation suitability index. The job area generation module is used to merge the LiDAR priority area and the transition area into a LiDAR candidate job area, and to perform spatial post-processing on the LiDAR candidate job area to generate at least one LiDAR flight sub-survey area. The flight path generation module is used to generate InSAR flight paths based on the InSAR priority area and to generate LiDAR flight paths based on the LiDAR flight sub-survey area. The observation suitability indicators include the percentage of overlapping area, the percentage of shaded area, the percentage of insufficiently covered area, the percentage of connected patches with limited observation area, and a comprehensive suitability indicator. The comprehensive suitability index Q is calculated using the following formula: Q = 0.35p L +0.45p S +0.20p E p L p represents the percentage of overlapping area. S p represents the proportion of the shaded area. E The percentage of areas not covered; Based on the observation suitability index, the area to be measured is divided into an InSAR priority area, a LiDAR priority area, and a transition area, specifically including: The area to be tested is divided into several grid units according to a preset spatial scale; For each grid cell, the distribution of observation types of all pixels within that grid cell is statistically analyzed, and a comprehensive suitability index is calculated. When p S <0.012, p L When Q < 0.035 and Q < 0.030, the grid cell is designated as an InSAR priority area; When p S ≥0.040, or p L +p S ≥0.080, or Q≥0.050, or p C When p ≥ 0.100, the grid cell is designated as a LiDAR priority region, where p C To observe the proportion of restricted connected patch area; The remaining grid cells are designated as transition zones.

Images