A multi-dimensional risk map dimension reduction modeling method for complex low-altitude operation

Abstract

The application belongs to the technical field of unmanned aerial vehicle low-altitude airspace operation safety and management, and particularly relates to a multi-dimensional risk map dimension reduction modeling method for complex low-altitude operation; the method is based on constructing an adaptive two-dimensional basic grid and mapping multi-source geographic data; a three-dimensional grid index of dimension reduction is constructed through non-equidistant height layering; in combination with an aerodynamic model and double physical shielding characteristics, dynamic safety risk and dynamic noise risk of nodes are respectively solved; finally, adaptive weighted fusion is carried out based on the ground surface cover semantic to generate a three-dimensional risk map. The application effectively reduces the storage and calculation overhead of spatial data, improves the accuracy of safety and noise evaluation and the scene adaptive ability under complex terrain, innovatively constructs a dimension reduction three-dimensional grid index architecture, and fuses a multi-physical field coupling quantization mechanism, so as to effectively realize high-precision, low-delay risk modeling of complex airspace from the aspects of bottom data structure and quantization mechanism.

Description

Technical Field

[0001] This invention relates to the field of safety and management technology for low-altitude airspace operation of unmanned aerial vehicles (UAVs), specifically to a method for dimensionality reduction modeling of multi-dimensional risk maps for complex low-altitude operations. Background Technology

[0002] In recent years, with the rapid development of the low-altitude economy, the high-density operation of drones in complex urban environments has become increasingly common. To ensure flight safety while also considering the environmental concerns of urban residents, constructing high-fidelity 3D risk maps has become a crucial prerequisite for intelligent drone trajectory planning. Regarding the 3D mapping problem of drones under the dual constraints of safety and noise in complex urban low-altitude environments, existing technologies mainly suffer from the following three limitations: 1. Spatial discretization represents redundancy, resulting in high computational and storage overhead.

[0003] Existing publicly available literature (e.g., Chinese invention patent application CN116126026A, "A Three-Dimensional Path Planning Method for UAVs Based on Risk Maps," and Chinese invention patent application CN121053827A, "Real-Time Planning of Low-Altitude Traffic Flow Airspace") typically uses uniformly spaced three-dimensional grids or voxel structures to divide the task airspace. This approach generates a large number of redundant nodes in the vertical and horizontal directions, leading to a surge in memory consumption of the underlying data and low path retrieval efficiency. Existing three-dimensional voxel modeling methods still fail to address the technical bottleneck of limited airborne computing power in real-time processing of massive grid data in large-scale complex scenes from the underlying data structure perspective.

[0004] 2. Risk quantification models are detached from real physical scenarios, resulting in insufficient assessment accuracy.

[0005] Existing safety risk calculation schemes (e.g., Chinese invention patent application CN120947650A, "An Intelligent Flight Path Planning Method for Unmanned Aerial Vehicles," and Chinese invention patent application CN121167801A, "An Intelligent Site Selection Method and System for UAV Take-off and Landing Points in Urban Building Clusters") are generally based on simplified line-of-sight distance evaluations or rely too heavily on static offline simulations, neglecting the nonlinear aerodynamic drag effect of dynamic fluctuations in air density caused by changes in altitude during actual flight on the aircraft's fall kinetic energy. Furthermore, in terms of noise risk assessment, existing research often ignores the directional shielding effect of the UAV's fuselage structure on sound waves, as well as the physical obstruction and diffraction of sound wave propagation by complex urban building clusters, leading to significant deviations between the risk quantification results and the actual urban physical environment.

[0006] 3. Multi-source risk fusion lacks spatial semantic awareness and has poor scene adaptability.

[0007] When dealing with the cost fusion of multiple objectives such as safety and noise, existing technologies (e.g., Chinese invention patent application CN121187320A, "Adaptive Planning Method for UAV Optimal Path Using Particle Swarm Optimization," and Chinese invention patent application CN120748257A, "Optimization Method and System for Traffic Control in Low-Altitude Economic Zones") largely rely on globally unified fixed weight allocation or purely mathematical multi-objective optimization, neglecting the actual functional attributes of the underlying surface environment. They also fail to realize the differentiated characteristics of different urban functional zones (such as the high noise sensitivity of medical and educational zones and the high collision sensitivity of industrial zones), lacking a dynamic weight adjustment mechanism based on environmental semantics, resulting in generated 3D risk maps that are difficult to conform to the actual norms of refined urban management.

[0008] Therefore, there is an urgent need for a multidimensional risk map construction method that balances computational lightweightness and physical high fidelity, in order to effectively optimize the underlying data index architecture and improve the accuracy of risk assessment and fusion in complex environments. Summary of the Invention

[0009] Faced with the dual constraints of safety and noise in complex low-altitude environments, traditional risk mapping methods struggle to balance efficient airborne computation and high-fidelity assessment. This invention aims to overcome the shortcomings of existing technologies by providing a dimensionality reduction modeling method for multi-dimensional risk maps in complex low-altitude operations. The technical problem addressed by this invention lies in its innovative construction of a dimensionality-reduced three-dimensional grid index architecture, coupled with a multi-physics coupling quantization mechanism. By addressing both the underlying data structure and the quantization mechanism, it effectively achieves high-precision, low-latency risk modeling in complex airspace.

[0010] To solve the above-mentioned technical problems, the specific technical solution adopted by the present invention is as follows: A method for dimensionality reduction modeling of multidimensional risk maps for complex low-altitude operations includes the following steps: S1: Construct a semantically adaptive two-dimensional planar basic grid; obtain the preset spatial constraint rules of the target area, adaptively divide the two-dimensional planar basic grid according to the environmental complexity, and map the multi-source geographic data attributes and surface cover semantics to the two-dimensional planar basic grid; S2: Generate a dimensionality-reduced three-dimensional mesh index architecture; perform non-equidistant height layering based on the low-altitude airspace hazard characteristics, and generate a one-dimensional virtual height label array in parallel for the two-dimensional plane base mesh; extract the planar coordinates of the two-dimensional plane base mesh as the spatial index primary key, and associate and map the attribute values ​​of the virtual height label array to construct a dimensionality-reduced three-dimensional mesh index; S3: Solve the nonlinear dynamic safety risk; use the spatial index primary key to retrieve and traverse the virtual altitude labels of each flight-ready node, take the absolute elevation of the ground surface or ground obstacle in the digital elevation model as the physical impact surface, calculate the vertical difference between the virtual altitude label and the physical impact surface as the relative fall height, couple the nonlinear aerodynamic drag model that changes dynamically with altitude with the physical buffering and shielding characteristics of the ground cover, and solve the node dynamic safety risk value. S4: Solve the dynamic noise risk of dual physical occlusion; combine the acoustic spherical geometric divergence model, couple the directional shielding attenuation of the UAV fuselage structure with the ray tracing blocking characteristics of environmental obstacles, and solve the dynamic noise risk value of each node in the reduced-dimensional three-dimensional mesh index; S5: Generate a semantically driven 3D risk map; remove non-flying grids under the preset spatial constraint rules, extract surface cover semantics to construct an adaptive weighting operator, and use the adaptive weighting operator to dynamically weight and fuse dynamic safety risk values ​​and dynamic noise risk values ​​to generate a 3D risk map based on planar coordinates and associated virtual height labels as the underlying addressing basis.

[0011] Furthermore, in step S3, the nonlinear dynamic safety risk is calculated, specifically including: extracting the virtual height label of the current node, using the absolute elevation of the ground surface or ground obstacle in the digital elevation model as the physical impact surface, and calculating the vertical difference between the virtual height label and the physical impact surface as the relative fall height; thereby eliminating the local fall height calculation deviation caused by the unified high-altitude mapping, and accurately reflecting the gravitational potential energy conversion during the UAV fall process; specifically, the calculation is performed through the following mathematical model: Introducing air density that dynamically varies with altitude A nonlinear aerodynamic drag model was constructed to calculate the vertical fall velocity. : ; In the formula, For the quality of drones, It is the acceleration due to gravity. The drag coefficient, For windward area, This refers to the relative fall height. Calculation including horizontal impact velocity The kinetic energy of the terminal impact: ; In the formula, The impact kinetic energy of the terminal falling. The horizontal flight speed of the drone at the time of impact is expressed in m / s. Extract the surface physical buffer shielding coefficient (physical buffer shielding characteristics) mapped to the two-dimensional planar base grid, and calculate the single-person injury probability using the Logistic lethality model. : ; In the formula, This is the injury and death sensitivity adjustment coefficient, with a value ranging from 0.05 to 0.2; This is the surface physical buffering and shielding coefficient, with a value ranging from 0 to 1; The lethal energy threshold; When the population density of the target area is lower than the preset minimum threshold (preset value), the environmental penalty mechanism is triggered to superimpose the terrain collision equipment damage penalty constant, and a dynamic safety risk value is generated through comprehensive calculation.

[0012] Furthermore, in step S4, the dynamic noise risk of dual physical occlusion is calculated, specifically including: Combining the spherical geometric divergence model of acoustic waves, the spatial distance from the UAV node to the ground receiving point is calculated, and the spherical geometric divergence attenuation term of acoustic waves propagating in free space is solved; based on the observation angle between the node center and the ground receiving point of the dimension-reduced three-dimensional grid, the spatial obstruction range of the UAV fuselage structure is matched to establish a directional attenuation function with spatial anisotropy characteristics. By utilizing the ray tracing blocking characteristics of environmental space obstacles, a three-dimensional ray tracing algorithm is used to extract the absolute elevation of surface obstacles on the spatial sound wave propagation path for environmental occlusion determination, calculate the environmental occlusion attenuation, and perform superposition correction to obtain the corrected sound pressure level. The noise sensitivity coefficient obtained based on the semantic mapping of land cover is used, combined with the exponential evaluation function of the environmental background noise tolerance threshold, to convert and solve the modified sound pressure level.

[0013] The calculation involves determining the relative geometric relationship between the node centers of the reduced-dimensional 3D mesh and the ground receiving point to obtain the straight-line distance of noise propagation and the observation angle; and establishing a directional shielding attenuation function by combining the physical shielding characteristics of the UAV's fuselage structure for sound waves. ; By introducing a directional shielding attenuation function Make corrections to ensure the corrected sound pressure level reaches the ground receiving point. The specific calculation model is as follows: ; In the formula, For reference distance The measured sound pressure level of the drone at the location represents the initial noise energy level; This is the spherical geometric divergence attenuation term for sound waves; This refers to the amount of air absorption attenuation. This represents the attenuation caused by environmental obstacles.

[0014] Based on the observed angle Extract the area of ​​sound wave propagation obstruction by the UAV fuselage structure and match the corresponding directional shielding attenuation function. Perform fuselage structure occlusion correction; extract the absolute elevation of surface obstacles along the spatial sound wave propagation path to determine environmental occlusion and perform environmental occlusion overlay correction; obtain the corresponding noise sensitivity coefficient based on the semantic mapping of surface cover in the two-dimensional planar base grid. The noise sensitivity coefficient is used to perform an exponential transformation on the modified sound pressure level to calculate the dynamic noise risk value of the three-dimensional mesh. : ; in, As the amplification factor, This is the environmental background noise tolerance threshold.

[0015] Among them, the amplification factor This is an exponential evaluation function; used to characterize the nonlinear increase in annoyance level among the population as noise energy increases after the sound pressure level exceeds the environmental background noise tolerance threshold; preferably, an amplification factor. The value range is from 0.10 to 0.25, amplification factor. The equivalent sound level penalty is pre-calibrated or dynamically configured based on the acoustic psychology experimental data of the target area.

[0016] Furthermore, in step S4, when calculating the dynamic noise risk of dual physical occlusion, a directional shielding attenuation function with spatial anisotropy characteristics is constructed. The specific logic is as follows: The observation angle between the node center of the reduced 3D mesh and the ground receiving point is calculated. The observation angle is the angle between the UAV propeller disk plane and the line connecting the observation point. A piecewise step function is used to perform directional correction attenuation calculations on the spatial sound field. ; In the formula, when When the region is determined to be an unobstructed radiation zone, the first attenuation constant is matched. ;when The region was determined to be a lateral, mild to moderate structurally obstructed radiation zone, and a second attenuation constant was matched. ;when The region was identified as a severely obstructed radiation zone by an overhead structure, and a third attenuation constant was matched. And satisfy .

[0017] Further, in step S5, a semantically driven 3D risk map is generated, specifically including: using the separating axis theorem to perform 3D rigid body collision detection between UAVs and terrain buildings; performing hard constraint culling on the dimensionality-reduced 3D mesh index using a combined comprehensive constraint potential field; determining a differentiated weighting strategy based on the functional area sensitivity attributes of the land cover semantics, and generating an adaptive weighting operator; using normalization processing to eliminate the dimensional differences between dynamic safety risk values ​​and dynamic noise risk values; and using the adaptive weighting operator to perform weighted summation on the normalized data to generate a comprehensive access cost score.

[0018] Among them, the normalized security risk value is obtained by using an adaptive weighting operator. With noise risk value The data is then integrated to generate a comprehensive passage cost score. : ; in, and The weighting coefficients are adaptively assigned; this is used to construct a three-dimensional risk map for use by the path planning algorithm.

[0019] Compared with the prior art, the present invention has the following beneficial effects: 1. Effectively reduces storage and computing overhead, adapting to limited onboard computing power: The dimension-reduced 3D mesh indexing architecture proposed in this invention improves upon the traditional uniform 3D voxel model, which is prone to generating numerous redundant nodes in large-scale mapping. By employing non-equidistant height layering and eliminating invalid traversals of high-altitude safety zones, the total number of spatially discretized nodes and memory usage are effectively reduced, improving the computational speed of collision detection and spatial retrieval, and meeting the real-time planning requirements of UAVs under limited onboard computing power.

[0020] 2. Reduce errors from idealized assumptions and improve the physical fidelity of risk assessment: To address the common assumptions of free fall and isotropic spherical divergence in traditional risk assessments, this invention introduces a coupling mechanism of nonlinear dynamics and dual acoustic shielding. By calculating the aerodynamic drag caused by changes in air density with altitude, and quantifying the directional shielding of the fuselage structure and the environmental shielding effect of complex urban building complexes, the quantitative results of multidimensional risks more closely reflect the real physical environment of urban operations.

[0021] 3. Improve the singularity of multi-source fusion and achieve refined weighting for scene adaptation: Unlike traditional global fixed weight allocation or single multi-objective optimization, this invention constructs an adaptive weighted operator based on land cover semantics. This operator can dynamically allocate safety and noise fusion weights by combining the differentiated characteristics of different urban functional areas such as medical and educational zones (high noise sensitivity) and industrial zones (high crash sensitivity), making the output three-dimensional risk map more in line with the actual needs of refined urban low-altitude management. Attached Figure Description

[0022] Figure 1 A flowchart of the construction method provided in an embodiment of the present invention; Figure 2 A schematic diagram of traditional true 3D voxel modeling using existing technologies; Figure 3 This is a schematic diagram of a dimension-reduced three-dimensional grid index architecture provided in an embodiment of the present invention; Figure 4 A schematic diagram illustrating the principle of directional physical masking of noise provided in an embodiment of the present invention; Figure 5 A flowchart illustrating the data processing for constructing a multidimensional risk map, as provided in this embodiment of the invention. Detailed Implementation

[0023] To make the objectives, technical solutions, and advantages of this invention clearer, the execution process of this invention will be described in detail below with reference to specific embodiments.

[0024] This embodiment provides a method for dimensionality reduction modeling of multidimensional risk maps for complex low-altitude operations. Taking a drone logistics delivery task spanning a high-rise commercial center, a medical and health area, and undulating suburban mountains as an example, the application process of this modeling method based on multiphysics coupling and semantic driving is explained in detail. Figure 1 As shown. Specifically, it includes the following steps: S1: Construct a semantically adaptive two-dimensional planar basic grid; obtain the preset spatial constraint rules of the target area, adaptively divide the two-dimensional planar basic grid according to the environmental complexity, and map the multi-source geographic data attributes and surface cover semantics to the two-dimensional planar basic grid; The preset spatial constraint rules specifically include: hard restrictions set for the safety of low-altitude flight of UAVs, policies and regulations or physical environment, including but not limited to at least one of the following: no-fly zone constraints, height-restricted zone constraints, terrain and object rigid body collision constraints, and near-ground restriction layer constraints at a specific altitude.

[0025] Digital elevation model data, building vector surface data, and land use type semantic data from multi-source geographic data attributes are used to extract terrain slope and building density features, which are then normalized and superimposed to generate a planar comprehensive complexity field. .

[0026] Using spatial connectivity algorithms, the parsed population density (regional population distribution density) is obtained. Height of surface obstacles The absolute elevation of the land surface and the physical buffer shielding coefficient extracted by the digital elevation model. The semantics of land cover (such as functional zone attributes like medical zones and industrial zones) are mapped to the corresponding two-dimensional grid attribute fields. Among these, the physical shading coefficient... The data is obtained by looking up the surface cover type in a table. For example, the physical buffering coefficient of the surface cover in dense forest areas is assigned a value of 0.5, and the physical buffering coefficient of the surface cover in hard flat areas is assigned a value of 1.0. This completes the construction of the underlying basic data matrix.

[0027] S2: Generate a dimensionality-reduced 3D mesh index architecture.

[0028] To address the shortcomings of traditional large-scale 3D voxel modeling, which easily generates massive redundant nodes and results in huge data volumes, this step abandons the globally uniform 3D spatial discretization mode. Based on a 2D planar base mesh, and combined with the hazard characteristics of low-altitude airspace, non-equidistant height layering is performed to construct a lightweight dimensionality-reduced 3D mesh index.

[0029] This embodiment abandons the traditional true 3D voxel modeling method, which is prone to generating massive amounts of redundant data, such as... Figure 2 As shown, a dimensionality reduction mapping architecture is adopted, such as... Figure 3 As shown.

[0030] Specifically, this dimensionality reduction mapping architecture uses the coordinates of a two-dimensional planar base grid as the spatial primary key. For example... Figure 3 As shown below, the planar grid G ​​(coordinates) ) and planar grid H (coordinates) This forms the underlying two-dimensional plane index; combined with Figure 3 As shown, the system uses digital elevation model data as the reference surface and sets the study airspace as 0 to 120m above the ground surface. Height parameters are marked along the Z-axis. , , , , This reflects the adaptive hierarchical logic of different regions.

[0031] For different two-dimensional planar base grids, the system generates a one-dimensional virtual height label array in parallel. For example... Figure 3 As shown in the left column, it is a non-equidistant label generated by a two-dimensional planar grid G, where C is the near-ground confinement layer and A is the high-altitude cruise layer; Figure 3 The pillar on the right is a non-equidistant label generated by the two-dimensional planar base mesh H, where D is the near-ground constraint layer, E is the obstacle influence analysis layer, F is the hard constraint obstacle layer, and B is the high-altitude cruise layer.

[0032] Furthermore, the baseline sampling coefficients The value range was determined through offline simulation calibration. A minimum thickness threshold for typical urban suspended obstacles was set to deduce the missed detection step size, thereby determining the baseline sampling coefficient. First, based on the safety upper limit, the minimum step size is calculated by inversely based on the memory allocation threshold of the onboard computing unit (e.g., a maximum allowed number of retrieval nodes of 1 million) and the real-time response time requirement of the algorithm (e.g., <0.5 seconds), and then the baseline sampling coefficient is determined. The lower limit of computing power. In the typical urban high-density building scenario of this embodiment, after multiple simulation optimizations, when the benchmark sampling coefficient... When the system is in the range of 1.5 to 3, it achieves the optimal Pareto front that minimizes the total time spent on global node retrieval, while keeping the obstacle miss rate at 0%.

[0033] Three-dimensional hazard field The value is a normalized superposition of the obstacle collision probability and the ground personnel exposure risk of the corresponding spatial node, with a value range of 0 to 1. In this embodiment, the obstacle collision probability is calculated based on the vertical distance between the node and the obstacle, and the ground personnel exposure risk is obtained by normalizing the population distribution density of the corresponding two-dimensional grid.

[0034] Based on this, the system performs adaptive height layering for the target airspace: (1) Near-ground constraint layer (0 to 30m): The height range of 0 to 30m above the ground surface is determined as the near-ground constraint layer, and the nodes corresponding to this range are directly marked as non-flying hard constraint attributes; (2) High-risk obstruction zone (30m to 60m): Due to the dense density of building tops in this zone, the risk field height gradient changes drastically, i.e. The sample size is relatively large, and the system automatically derives the dense sampling step size based on the above formula, such as... This generates a dense sequence of virtual height labels containing 35m, 40m, 45m, 50m, 55m, and 60m. (3) High-altitude cruising zone (60m to 120m): This zone has few obstacles and a relatively gentle altitude gradient in the risk field, i.e. When the step size is small, the step size derived from the formula automatically increases, such as... This allows for the sparse generation of virtual height label sequences containing 70m, 80m, and up to 120m.

[0035] The final system uses the underlying two-dimensional planar coordinate key value Using the primary key, the aforementioned non-equidistant virtual height layer sequences are treated as one-dimensional arrays and mapped in parallel to construct a dimensionality-reduced three-dimensional grid index combining two-dimensional key-value pairs and one-dimensional height arrays. This architecture effectively eliminates unnecessary traversal of high-altitude safe zones while preserving the accuracy of three-dimensional topology obstacle avoidance, achieving efficient and lightweight three-dimensional representation of the entire airspace.

[0036] S3: Solving for nonlinear dynamic safety risks. For the flight-safe nodes in the 3D mesh index, this method overcomes the mechanistic distortion of traditional models that treat crashes as idealized free fall, and assesses the physical lethal risk of drone crashes.

[0037] For the current computation node, extract the virtual height label of the current node. Using the absolute elevation of the ground surface or surface obstacle in the digital elevation model as the physical impact surface, calculate the vertical difference between the virtual height label and the physical impact surface as the relative fall height. And based on the standard atmospheric model, the air density at the corresponding altitude was matched. When calculating the collision kinetic energy, calibration is performed based on the physical parameters of a typical logistics drone: the drone's mass is set. It weighs 9.2 kg and has a windward area. for drag coefficient It is 0.8.

[0038] Substituting the above parameters into the nonlinear aerodynamic drag model, the vertical fall velocity is first calculated. and the impact kinetic energy of the terminal falling The calculation formula is as follows: ; ; In the formula, let the horizontal impact speed of the drone during its cruise be denoted as . for .

[0039] Furthermore, the calculated terminal impact kinetic energy Substitute into the logistic fatality rate model to calculate the probability of single-person injury or death. : ; In this embodiment, the lethal energy threshold The value is set at 80 joules, which is based on the international civil aviation authority's general threshold standard for assessing the lethality of ground impact kinetic energy of unmanned aerial vehicles.

[0040] Surface physical shielding coefficient The surface physical shading coefficient is obtained by looking up the semantic type of the land cover mapped in step S1. It is not an arbitrary empirical value, but rather a standardized mapping dictionary based on land cover semantics pre-built by the system. This characterizes the absorption and buffering capacity of different surface physical materials and spatial structures for the impact kinetic energy of a UAV crash, with values ​​ranging from (0,1). For example, on hard, flat ground, since there are no physical obstructions, the impact kinetic energy is not attenuated, so the value is 1.0. In dense forest areas, since the tree canopy can absorb and buffer some of the impact kinetic energy, the value is 0.5.

[0041] Injury and Casualty Sensitivity Adjustment Coefficient The value is 0.12. This coefficient is a conventional empirical fitting parameter for the Logistic fatality rate curve, used to characterize the nonlinear growth rate of the probability of injury or death after the kinetic energy of the fall exceeds the energy threshold.

[0042] Furthermore, the basic failure rate of drone equipment is introduced. (This embodiment takes) (number of flights / flight hours), and population density obtained by combining two-dimensional grid mapping. The dynamic security risk quantification value of the node is obtained through product operation. The complete formula for its calculation is as follows: ; Environmental constraint verification logic is introduced during the security risk quantification process: if the population density of the current node... If the value falls below a preset minimum (e.g., less than 10 people / square kilometer), an environmental penalty mechanism is triggered, forcing [the user to] [remain in a certain position]. The preset terrain and equipment damage penalty constant is superimposed on the basis. To avoid planning dangerous flight paths that are too close to the ground in uninhabited areas and could easily collide with mountains.

[0043] S4: Solving for dynamic noise risk due to dual physical occlusion. For each node in the reduced-dimensional 3D mesh index, the system comprehensively considers geometric divergence attenuation, fuselage structure shielding, and environmental obstacle occlusion to perform noise risk quantification. Specifically, this includes: Combination Figure 4 The diagram shows a schematic of the UAV noise propagation geometric model (spherical acoustic divergence model). The relative geometric relationship between the node centers of the reduced-dimensional 3D mesh and the ground receiving point is calculated. Figure 4 As shown, The center sound power level of the drone's sound source; The sound pressure level at the ground-based sound wave receiver point. This refers to the noise propagation distance. The observation angle (i.e., the spatial angle between the plane of the UAV propeller disk and the line connecting the observation point); This is the attenuation function of the fuselage's physical shielding. The model comprehensively considers the effect of noise propagation with distance. The geometric divergence attenuation and the directional shielding effect brought about by the fuselage structure are used to construct the basic calculation model for low-altitude dynamic noise risk.

[0044] In order to overcome the physical distortion of the traditional spherical isotropic divergent model when calculating dynamic noise risk, the system introduces a directional shielding attenuation function with spatial anisotropy characteristics. .

[0045] Considering the actual physical structure of a typical multi-rotor logistics drone, sound wave propagation experiences significant directional attenuation due to the unobstructed area below the rotor, while being blocked laterally and upward by the battery compartment and the top structure of the fuselage, respectively. In this embodiment, the directional shielding attenuation function... It is constructed using a piecewise step function, the specific expression of which is as follows: ; In the formula, when At that time, the observation point was located directly below or slightly below the drone's propeller disk, and the acoustic radiation was minimally blocked by the fuselage. The baseline attenuation was then used. ;when At that time, the observation point was located to the side, and the propagation of sound waves was physically blocked by the middle of the fuselage and the battery compartment (such as a typical battery compartment size of 30cm×20cm×15cm), resulting in additional structural attenuation. ;when At that time, the observation point was located above the drone, and the sound waves were severely blocked by the main body on the top of the drone. The extreme value of attenuation was taken. .

[0046] By introducing this The function is modified to adjust the sound pressure level reaching the ground receiving point. The specific calculation model is as follows: ; In the formula, For reference distance The measured sound pressure level of the drone at the location represents the initial noise energy level, such as when the DJI Matrice 350 RTK is hovering. , ; This is the spherical geometric divergence attenuation term for sound waves; This represents the air absorption attenuation; it characterizes the energy consumed by sound waves in overcoming air viscosity. This embodiment adopts the ISO 9613-1 outdoor sound propagation attenuation standard published by the International Organization for Standardization. The measurement was performed at 20 degrees Celsius and relative humidity... Under typical urban climate conditions, the air absorption attenuation coefficient in the mid-frequency band is set to Then the system dynamically solves .

[0047] In this embodiment, is the directional shielding attenuation function of the fuselage, used to characterize the spatial shielding noise reduction effect of the fuselage structure; This represents the attenuation caused by environmental obstacles. The system employs a 3D ray tracing algorithm, emitting virtual sound waves from the UAV node to the ground receiving point. If the ray intersects with a building polygon in the 3D geographic information model (i.e., line-of-sight occlusion exists), the Fresnel half-wave zone diffraction model is further invoked to dynamically calculate the attenuation based on the path difference of the sound wave diffracting around the building edge. This attenuation value typically fluctuates dynamically between 5 dB(A) and 20 dB(A). If the rays do not intersect, then .

[0048] Call Mapping the dynamic noise risk value of the 3D mesh .

[0049] Noise sensitivity coefficient To achieve core adjustment coefficients for underlying spatial semantic perception, the system pre-constructs a mapping dictionary between surface cover semantics and sensitivity coefficients. In this embodiment, for highly sensitive semantic grids such as hospitals and schools, values ​​are assigned... For regular residential areas, assign values For low-sensitivity semantic meshes such as factories and main roads, assign values... This coefficient directly determines the severity of risk amplification once the tolerance threshold is exceeded.

[0050] In this actual calculation process, the amplification factor in the exponential function This plays a crucial role in determining the nonlinear growth of noise risk penalty weights with exceeding sound pressure levels. To ensure the generated cost map accurately adapts to the refined management needs of different time periods, this embodiment... Context-based dynamic value assignment: When drones perform routine urban logistics delivery tasks during the daytime, the system has a relatively high tolerance for new noise due to the high background noise level in the daytime environment. A value of 0.15 is taken as an empirical reference. If the drone is performing a mission at night, the background noise level is extremely low, and people are extremely sensitive to even the slightest noise. In this case, the system will... It was dynamically adjusted to 0.23.

[0051] Environmental background noise tolerance threshold The red line for triggering exponential noise risk penalties is dynamically obtained through table lookup based on relevant standards such as the "Environmental Quality Standard for Noise" (GB3096-2008) of the People's Republic of China. For example, when a drone flies over a Class 1 noise environment functional zone (such as residential areas or medical and health areas), the daytime... The daytime calibrator is set at 55 dB(A), and the nighttime calibrator is set at 45 dB(A); when flying over Class 3 acoustic environment functional zones (such as industrial production areas), the daytime calibrator is set at 45 dB(A). It is calibrated to 65 dB(A).

[0052] This dynamic parameter adjustment allows for the correction of the sound pressure level reaching the ground. Once the tolerance threshold is exceeded The calculated dynamic noise risk value This will increase exponentially, forcing subsequent trajectory planning algorithms to actively avoid noise-sensitive areas such as medical facilities and residential areas with stricter hard constraints.

[0053] S5: Generate a semantically driven 3D risk map. In the cost fusion stage, unlike traditional global fixed weight allocation or single mathematical multi-objective optimization, this method introduces an adaptive weighting operator based on land cover semantics.

[0054] Hard constraint culling is performed: the Separated Axis Theorem (SAT) is used to conduct 3D rigid body collision detection between the UAV's 3D bounding box and the building boundary. For meshes located within the near-ground constraint layer at an altitude of 0 to 30m above the ground surface or meshes experiencing rigid body interference, the system utilizes a comprehensive constraint potential field. Either eliminate it directly or set the cost to infinity.

[0055] Perform constraint filtering judgment: Obtain the absolute elevation of the grid and compare it with the sum of the height of surface obstacles, the relative ground elevation, and the preset safety margin. If a collision condition is triggered, directly set the passage cost of the node to the preset infinite penalty value.

[0056] Adaptive weighting is performed: First, Min-Max normalization is applied to the dynamically calculated security risk values. With dynamic noise risk value Mapping to a unified dimension interval [0,1] yields the normalized safety risk value. With noise risk value .

[0057] Based on the semantic labels of land cover obtained from the underlying two-dimensional planar grid mapping, a corresponding weighting strategy is matched. For example, when the semantic label is identified as medical and health land, a noise suppression priority allocation strategy is implemented, and noise weights are configured. The safety weight is 0.7. The value is 0.3; when the semantic label is identified as industrial land, an anti-collision priority allocation strategy is implemented, and noise weights are configured. The safety weight is 0.2. It is 0.8.

[0058] Finally, the adaptive weighting operator is used to perform a weighted fusion of the normalized safety risk value and the noise risk value to generate a comprehensive passage cost score for the 3D mesh node. The calculation formula is as follows:

[0059] like Figure 5 The data processing flow for constructing a multi-dimensional risk map is shown. The system constructs a three-dimensional matrix from the cost scores of all nodes in the region, thus completing the modeling of the three-dimensional risk map. The underlying data matrix, which includes multi-dimensional penalty weights, can be directly used by path planning algorithms such as multi-objective particle swarm optimization to plan the optimal low-altitude flight path that balances physical safety and noise pollution.

[0060] To quantify and verify the beneficial technical effects of this embodiment in reducing computational and storage overhead, real urban geographic environment data was extracted for software simulation testing. The selected target area covered an area of ​​2km x 2km, with terrain features including dense commercial areas and undulating terrain, and an airspace height limit of 120m. The test used the controlled variable method to compare the traditional uniform three-dimensional voxel modeling method (using a global 5m x 5m x 5m equidistant resolution) with the adaptive dimensionality reduction three-dimensional mesh indexing method.

[0061] The hardware environment used for testing was: Intel Core i7 processor, 16GB RAM, and NVIDIA RTX 3050 high-performance graphics processor; the software environment was: Python 3.11. After data processing and algorithm execution, the comparative test results are shown in Table 1. Table 1: Comparison of mapping and computational costs for 3D risk maps; Total number of spatial nodes generated Approximately 3,840,000 Approximately 910,000 The number of nodes decreased by approximately 76.3%. Underlying data memory usage Approximately 14.65MB Approximately 3.47MB Memory overhead was reduced by approximately 76.3%. Time taken for 1000 ray-tracing collision detections Approximately 24.00 seconds Approximately 8.00 seconds The processing speed has increased by approximately 200%. As shown in Table 1, the test data demonstrates that this method achieves significant lightweighting while maintaining near-ground obstacle avoidance accuracy. It effectively eliminates redundant high-altitude nodes, reducing the total number of generated nodes and underlying memory usage by approximately 76.3%, directly overcoming the airborne storage bottleneck in complex, large-scale mapping. By limiting spatial retrieval to the effective influence layer, it eliminates invalid traversal of high-altitude areas, accelerating dynamic calculations by approximately 200%, and significantly reducing the time required for collision detection such as 3D ray tracing, thus meeting the low-latency requirements of UAV low-altitude real-time trajectory planning.

[0062] The above description is merely a preferred embodiment of the present invention and is not intended to limit the scope of protection of the present invention. Those skilled in the art can make various modifications and variations to the present invention without departing from its principles. Any modifications, equivalent substitutions, and improvements made within the principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A method for dimensionality reduction modeling of multidimensional risk maps for complex low-altitude operations, characterized in that, Includes the following steps: S1: Construct a semantically adaptive two-dimensional planar basic mesh; Obtain the preset spatial constraint rules of the target area, adaptively divide the two-dimensional plane basic grid according to the environmental complexity, and map the multi-source geographic data attributes and surface cover semantics to the two-dimensional plane basic grid; S2: Generate a dimensionality-reduced three-dimensional mesh index architecture; perform non-equidistant height layering based on the low-altitude airspace hazard characteristics, and generate a one-dimensional virtual height label array in parallel for the two-dimensional plane base mesh; extract the planar coordinates of the two-dimensional plane base mesh as the spatial index primary key, and associate and map the attribute values ​​of the virtual height label array to construct a dimensionality-reduced three-dimensional mesh index; S3: Solve the nonlinear dynamic safety risk; use the spatial index primary key to retrieve and traverse the virtual altitude labels of each flight-ready node, take the absolute elevation of the ground surface or ground obstacle in the digital elevation model as the physical impact surface, calculate the vertical difference between the virtual altitude label and the physical impact surface as the relative fall height, couple the nonlinear aerodynamic drag model that changes dynamically with altitude with the physical buffering and shielding characteristics of the ground cover, and solve the node dynamic safety risk value. S4: Solve the dynamic noise risk of dual physical occlusion; combine the acoustic spherical geometric divergence model, couple the directional shielding attenuation of the UAV fuselage structure with the ray tracing blocking characteristics of environmental obstacles, and solve the dynamic noise risk value of each node in the reduced-dimensional three-dimensional mesh index; S5: Generate a semantically driven 3D risk map; remove non-flying grids under the preset spatial constraint rules, extract surface cover semantics to construct an adaptive weighting operator, and use the adaptive weighting operator to dynamically weight and fuse dynamic safety risk values ​​and dynamic noise risk values ​​to generate a 3D risk map based on planar coordinates and associated virtual height labels as the underlying addressing basis.

2. The method for dimensionality reduction modeling of multidimensional risk maps for complex low-altitude operations according to claim 1, characterized in that, Step S3 involves calculating the nonlinear dynamic safety risk, specifically including: extracting the virtual height label of the current node; using the absolute elevation of the ground surface or surface obstacle in the digital elevation model as the physical impact surface; and calculating the vertical difference between the virtual height label and the physical impact surface as the relative fall height. Specifically, the solution is obtained through the following mathematical model: Introducing air density that dynamically varies with altitude A nonlinear aerodynamic drag model was constructed to calculate the vertical fall velocity. : ； In the formula, For the quality of drones, It is the acceleration due to gravity. The drag coefficient, For windward area, This refers to the relative fall height. Calculation including horizontal impact velocity The kinetic energy of the terminal impact: ； In the formula, The impact kinetic energy of the terminal falling. The horizontal flight speed of the drone at the time of impact is expressed in m / s. Extract the surface physical buffer shielding coefficient mapped to the two-dimensional planar base grid, and calculate the single-person injury probability using the logistic fatality rate model. : ； In the formula, This is the injury and death sensitivity adjustment coefficient, with a value ranging from 0.05 to 0.2; This is the surface physical buffering and shielding coefficient, with a value ranging from 0 to 1; The lethal energy threshold; When the population density of the target area is lower than the preset minimum threshold, an environmental penalty mechanism is triggered to superimpose a terrain collision equipment damage penalty constant, and a dynamic safety risk value is generated through comprehensive calculation.

3. The method for dimensionality reduction modeling of multidimensional risk maps for complex low-altitude operations according to claim 1, characterized in that, Step S4 calculates the dynamic noise risk of dual physical occlusion, specifically including: By combining the spherical geometric divergence model of sound waves, the spatial distance from the UAV node to the ground receiving point is calculated, and the spherical geometric divergence attenuation term of sound waves propagating in free space is solved. Based on the observation angle between the node center and the ground receiving point of the dimensionality-reduced 3D grid, the spatial occlusion range of the UAV fuselage structure is matched to establish a directional attenuation function with spatial anisotropy characteristics. By utilizing the ray tracing blocking characteristics of environmental space obstacles, a three-dimensional ray tracing algorithm is used to extract the absolute elevation of surface obstacles on the spatial sound wave propagation path for environmental occlusion determination, calculate the environmental occlusion attenuation, and perform superposition correction to obtain the corrected sound pressure level. The noise sensitivity coefficient obtained based on the semantic mapping of land cover is used, combined with the exponential evaluation function of the environmental background noise tolerance threshold, to convert and solve the modified sound pressure level.

4. The method for dimensionality reduction modeling of multidimensional risk maps for complex low-altitude operations according to claim 3, characterized in that, In step S4, when calculating the dynamic noise risk of dual physical blockage, a directional shielding attenuation function with spatial anisotropy characteristics is constructed. The specific logic is as follows: Calculate the observation angle between the node center of the reduced-dimensional 3D mesh and the ground receiving point. Observe the included angle The angle between the UAV propeller disk plane and the line connecting the observation point; a piecewise step function is used to calculate the directional correction attenuation of the spatial sound field: ； In the formula, when When the region is determined to be an unobstructed radiation zone, the first attenuation constant is matched. ;when The region was determined to be a lateral, mild to moderate structurally obstructed radiation zone, and a second attenuation constant was matched. ;when The region was identified as a severely obstructed radiation zone by an overhead structure, and a third attenuation constant was matched. And satisfy .

5. The method for dimensionality reduction modeling of multidimensional risk maps for complex low-altitude operations according to claim 1, characterized in that, Step S5 generates a semantically driven 3D risk map, specifically including: The separation axis theorem is used to perform 3D rigid body collision detection between UAVs and terrain structures, and a combined integrated constraint potential field is used to perform hard constraint culling on the dimensionality-reduced 3D mesh index. Based on the functional zone sensitivity attributes of land cover semantics, a differentiated weighting strategy is determined, and an adaptive weighting operator is generated; Normalization is used to eliminate the dimensional differences between dynamic safety risk values ​​and dynamic noise risk values. An adaptive weighting operator is used to perform weighted summation on the normalized data to generate a comprehensive passage cost score.

Images