An asymmetric modeling space risk assessment method and system

Abstract

This invention relates to the field of air traffic management and digital airspace safety assessment technology, and provides an asymmetric modeling airspace risk assessment method and system. The method includes: preprocessing raw trajectory data to obtain a trajectory sequence; statistically analyzing trajectory points and extracting the average velocity vector of each voxel, establishing a local trajectory coordinate system, generating an asymmetric three-dimensional Gaussian convolution kernel, and constructing a continuous three-dimensional probability density field; performing a macro-strategic risk quantification assessment on the three-dimensional probability density field to generate the total route risk; based on the relative motion relationship between the aircraft to be assessed and the intruding aircraft, performing a real-time dynamic collision risk quantification assessment on the aircraft to be assessed and the intruding aircraft on a preset route, calculating the peak collision risk index; and outputting a risk report and visualizing it. This invention transforms discrete trajectories into a continuous risk field, realizing a dual-layer risk quantification of macro-strategic and micro-dynamic risks, and improving the continuity and directional sensitivity of airspace risk assessment.

Description

Technical Field

[0001] This invention relates to the field of air traffic management and digital airspace safety assessment technology, and in particular to an asymmetric modeling airspace risk assessment method and system. Background Technology

[0002] With the continuous growth of air traffic, the refined management and utilization of airspace resources has become crucial. Existing airspace traffic modeling and risk assessment methods typically employ two-dimensional raster heatmaps or simple counting statistics of three-dimensional trajectory points. These methods suffer from several drawbacks. First, they struggle to provide consistent spatial resolution and continuous representation in the height dimension, leading to a weakening of traffic differences at different altitudes within the same latitude and longitude location. Second, the discrete statistical results exhibit significant voxel-like discontinuities, making them unsuitable as risk inputs for continuous integration along flight paths. Furthermore, raw flight path data commonly suffers from jumps, missing data, and uneven sampling rates; without systematic preprocessing, these issues will further amplify density estimation biases. Therefore, an asymmetric modeling airspace risk assessment method and system are urgently needed to address the problem of inaccurate airspace traffic risk assessment. Summary of the Invention

[0003] The purpose of this invention is to solve the technical problems of two-dimensional grid heatmaps completely ignoring the height dimension and failing to reflect the differences in safety risks in the vertical direction of the airspace; and three-dimensional trajectory point counting method, although introducing the height dimension, has a discrete statistical model that results in a "blocky" distribution, which violates the continuous motion characteristics of aircraft and is sensitive to data sampling frequency. It is also prone to "zero point" holes when there is noise or uneven sampling, causing abrupt changes in the assessment results and making it difficult to meet the needs of continuous path risk assessment. The invention provides an asymmetric modeling airspace risk assessment method and system.

[0004] To achieve the above-mentioned objectives, the embodiments of the present invention provide the following technical solutions:

[0005] An asymmetric modeling method for airspace risk assessment, the method comprising:

[0006] S1. Collect raw trajectory data, preprocess the raw trajectory data to obtain trajectory sequence;

[0007] S2. Perform spatial statistics on the trajectory sequence using a three-dimensional mesh, and record the average velocity vector of the trajectory points within each voxel; establish a local trajectory coordinate system based on the average velocity vector, generate an asymmetric three-dimensional Gaussian convolution kernel, and construct a three-dimensional probability density field through the asymmetric Gaussian convolution kernel.

[0008] S3. For the preset route input by the user, perform a macro-strategic risk quantification assessment on the three-dimensional probability density field to generate the total route risk;

[0009] S4. Based on the relative motion relationship between the aircraft to be evaluated and the intruding aircraft, perform real-time dynamic collision risk quantification assessment on the aircraft to be evaluated and the intruding aircraft on the preset flight path, and calculate the peak collision risk index.

[0010] S5. Output a risk report based on the total route risk and the peak collision risk index, and display it visually.

[0011] The three-dimensional probability density field constructed in this invention possesses spatial integrity and physical continuity compared to discrete counting models. By transforming discrete points into a continuously distributed "field," it can realistically simulate the spatial probability distribution of aircraft in the airspace, eliminating the "blocky" distortion caused by discrete grids and making airspace safety modeling more consistent with the physical nature of continuous aircraft movement. Based on the characteristics of the continuous density field, this invention transforms risk assessment into a "line integral along the path" problem. This modeling method can provide smooth and quantitative risk indicators for any preset route, not only solving the assessment bias caused by uneven sampling rates but also enabling rigorous horizontal comparability of risks between different routes and time periods.

[0012] Furthermore, the preprocessing in S1 includes:

[0013] Abnormal trajectory points are filtered and removed based on a first preset threshold.

[0014] Grouping by aircraft call sign, segmenting trajectories with time intervals greater than 30 minutes, and performing linear interpolation at 1-second intervals.

[0015] Furthermore, in S2, the three-dimensional probability density field is represented as:

[0016] ;

[0017] in, This represents the probability density of an aircraft appearing per unit volume at spatial location x; k is the voxel index; N k Count the trajectory points within voxel k; x k For the k-th voxel center, R is the local track coordinate transformation matrix. k The transpose of the matrix; is the asymmetric three-dimensional Gaussian convolution kernel corresponding to voxel k.

[0018] Furthermore, in S3, the asymmetric three-dimensional Gaussian convolution kernel is represented as:

[0019] ;

[0020] in, This represents a point p in space relative to the voxel center x. k The longitudinal forward offset; This represents a point p in space relative to the voxel center x. k The longitudinal backward offset; u y,k This represents a point p in space relative to the voxel center x. k Lateral offset; u z,k This represents a point p in space relative to the voxel center x. k Vertical offset; This represents a point p in space relative to the voxel center x. k Forward standard deviation; This represents a point p in space relative to the voxel center x. k The backward standard deviation; This represents a point p in space relative to the voxel center x. k Lateral standard deviation; This represents a point p in space relative to the voxel center x. k The vertical standard deviation; Z k This represents the normalization coefficient.

[0021] Furthermore, in S3, the quantitative assessment of macro-strategic risks includes:

[0022] S31. Divide the preset route into multiple preset segments, and extract the risk density value of the starting point of each preset segment in the three-dimensional probability density field by nearest neighbor search.

[0023] S32. Calculate the flight distance for each preset segment;

[0024] S33. Calculate the total route risk based on the risk density value, the flight distance, and the preset risk scaling constant.

[0025] Furthermore, in S4, the real-time dynamic collision risk quantification assessment specifically includes:

[0026] S41. Obtain the relative motion relationship between the aircraft to be evaluated and the intruding aircraft;

[0027] S42. Establish a local coordinate system for relative motion based on the relative motion relationship between the aircraft to be evaluated and the intruding aircraft;

[0028] S43. Combining the three-dimensional collision box model with the Laplace position error distribution, a real-time dynamic collision risk quantification assessment is performed on the aircraft to be evaluated and the intruding aircraft on the preset flight path, and the peak collision risk index is calculated.

[0029] Furthermore, the cumulative distribution function of the Laplace position error distribution is:

[0030]

[0031] in, For a random variable to take values, For position parameters, This is the scale parameter.

[0032] An asymmetric modeling airspace risk assessment system, comprising a preprocessing module, a density field construction module, a collision risk assessment module, a strategic risk assessment module, and a visualization output module:

[0033] The preprocessing module is used to collect raw track data and preprocess the raw track data to obtain a track sequence.

[0034] The density field construction module uses a three-dimensional mesh to perform spatial statistics on the trajectory sequence and records the average velocity vector of trajectory points within each voxel. A local trajectory coordinate system is established based on the average velocity vector, and an asymmetric three-dimensional Gaussian convolution kernel is generated. A three-dimensional probability density field is constructed through the asymmetric Gaussian convolution kernel.

[0035] The strategic risk assessment module is used to perform a macro-strategic risk quantification assessment on the three-dimensional probability density field for a preset route input by the user, and generate the total route risk.

[0036] The collision risk assessment module is used to perform real-time dynamic collision risk quantification assessment of the aircraft to be assessed and the intruding aircraft on a preset flight path based on the relative motion relationship between the aircraft to be assessed and the intruding aircraft, and to calculate the peak collision risk index.

[0037] The visualization output module is used to output a risk report based on the total route risk and the peak collision risk index, and to visualize the report.

[0038] Compared with the prior art, the beneficial effects of the present invention include:

[0039] In terms of airspace modeling: This invention does not merely perform isotropic smoothing on discrete trajectory points, but rather transforms discrete trajectory points into a continuously distributed "field" through a processing chain of "voxel frequency statistics—average velocity extraction—local coordinate system establishment—asymmetric convolution kernel generation," while retaining dominant traffic flow information in the convolution results. This realistically simulates the spatial probability distribution of aircraft in the airspace. It effectively eliminates the "blocky" distortion and "zero-value point" holes caused by discrete grids, making airspace safety modeling more consistent with the physical nature of continuous aircraft movement, and significantly improving the model's noise resistance and accuracy.

[0040] In terms of risk assessment: Based on a continuous density field with directional attributes, this invention can transform macro-strategic risk assessment into a path integral problem along a preset route; simultaneously, at the micro-real-time level, by coupling a three-dimensional collision box model in a local coordinate system of relative motion with a Laplace position error distribution, instantaneous dynamic collision risk indicators can be further obtained. Thus, this invention forms a two-layer quantitative structure of "continuous strategic risk field + micro-dynamic collision risk," enabling smooth, quantitative, and real-time reference-worthy risk assessment results for any preset route. Attached Figure Description

[0041] To more clearly illustrate the technical solutions of the embodiments of the present invention, the accompanying drawings used in the embodiments will be briefly introduced below. It should be understood that the following drawings only show some embodiments of the present invention and should not be regarded as a limitation on the scope. For those skilled in the art, other related drawings can be obtained based on these drawings without creative effort.

[0042] Figure 1 This is a flowchart of the method of the present invention;

[0043] Figure 2(a) is a schematic diagram of the distribution of track points before data preprocessing in this invention;

[0044] Figure 2(b) is a schematic diagram of the distribution of track points after data preprocessing in this invention;

[0045] Figure 3(a) is a schematic diagram of the three-dimensional mesh effect in this invention;

[0046] Figure 3(b) is a schematic diagram of the three-dimensional probability density field in this invention;

[0047] Figure 4(a) is a schematic diagram of the superposition of the three-dimensional probability density field and the preset route in this invention;

[0048] Figure 4(b) is a schematic diagram showing the change of risk accumulation with flight distance in this invention;

[0049] Figure 5 This is a schematic diagram of the system's visual interface in this invention;

[0050] Figure 6 This is a schematic diagram of the system of the present invention. Detailed Implementation

[0051] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. The components of the embodiments of the present invention described and shown in the accompanying drawings can generally be arranged and designed in various different configurations. Therefore, the following detailed description of the embodiments of the present invention provided in the accompanying drawings is not intended to limit the scope of the claimed invention, but merely to illustrate selected embodiments of the invention. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without inventive effort are within the scope of protection of the present invention.

[0052] It should be noted that similar labels and letters in the following figures indicate similar items. Therefore, once an item is defined in one figure, it does not need to be further defined and explained in subsequent figures.

[0053] Example 1:

[0054] This embodiment is achieved through the following technical solution, such as Figure 1 As shown in the figure, this embodiment provides an asymmetric modeling method for airspace risk assessment, the method comprising:

[0055] S1. Collect raw trajectory data, preprocess the raw trajectory data to obtain trajectory sequence.

[0056] The raw flight track data refers to the original record sequence containing aircraft motion information collected from data sources such as ADS-B, including timestamps (MSGTIME), aircraft call signs (CALLSIGN), longitude (LATTD), and flight altitude (FL). This raw flight track data is stored in CSV format and used for subsequent preprocessing, three-dimensional probability density field construction, and risk assessment calculations.

[0057] An aircraft call sign is a code used to uniquely identify an aircraft operating in airspace.

[0058] The preprocessing includes:

[0059] Data points are filtered based on a first preset threshold to remove obviously erroneous data points. In this embodiment, the first preset threshold includes a latitude and longitude jump threshold (e.g., 0.05 degrees, approximately 5.5 kilometers) and an altitude jump threshold (e.g., 100 feet, approximately 30.5 meters). Data points exceeding the first preset threshold are identified as outliers and removed.

[0060] The system groups aircraft callsigns and segments tracks with time intervals greater than 30 minutes, then performs linear interpolation at 1-second intervals. Specifically, the system groups all raw track data with the same aircraft callsign. For each group, it segments the raw track data within that group with time intervals greater than 30 minutes. Further, it uses a resampling method from the Pandas library to perform linear interpolation on each segment of raw track data at 1-second intervals, generating a continuous and smooth track sequence.

[0061] Figures 2(a) and 2(b) provide a schematic diagram comparing the trajectories before and after data preprocessing. It can be seen that the comparison effect of the trajectories before and after preprocessing is obvious: Figure 2(a) is a schematic diagram of the distribution of track points before data preprocessing in this invention. It can be seen that there are obvious abnormal jumps and uneven sampling phenomena in the distribution of original track data points. Some track points are outliers that deviate from the actual track due to data noise; Figure 2(b) is a schematic diagram of the distribution of track points after data preprocessing in this invention. After removing outliers whose latitude, longitude or altitude jumps exceed the threshold, linear interpolation is performed at 1-second intervals to make the track points evenly distributed along the track at equal intervals, thus eliminating the influence of outliers and uneven time intervals.

[0062] S2. Perform spatial statistics on the trajectory sequence using a three-dimensional mesh, and record the average velocity vector of the trajectory points within each voxel; establish a local trajectory coordinate system based on the average velocity vector, generate an asymmetric three-dimensional Gaussian convolution kernel, and construct a three-dimensional probability density field through the asymmetric three-dimensional Gaussian convolution kernel.

[0063] S2 includes:

[0064] Based on the longitude, latitude, and altitude ranges of the target airspace, a three-dimensional mesh is divided using voxels of fixed size, and the average velocity vector within each voxel is recorded simultaneously.

[0065] The target airspace refers to the specific three-dimensional spatial range that requires airspace modeling and risk assessment, and is determined by the longitude interval, latitude interval, and altitude interval.

[0066] The system first sets voxels of fixed dimensions. In this embodiment, the fixed dimensions of the voxels include: a longitude step of 0.00113 degrees (approximately 125 meters), a latitude step of 0.00113 degrees, and the altitude step is dynamically calculated according to the following formula. :

[0067] ;

[0068] Where K is the division index, used to determine the number of layers in the height direction; preferably, K=6; H max H represents the maximum altitude of the target airspace. max =2400 meters; Hmin H represents the minimum altitude of the target airspace. min =200 meters, then the altitude range of the target airspace is [200, 2400] meters. The altitude step size can then be calculated. It is approximately 34.4 meters.

[0069] The Tensor Grid statistical method is used to map the preprocessed trajectory sequence onto a 3D mesh, while simultaneously extracting the instantaneous velocity vector of the trajectory sequence and calculating and recording the voxel center x. k average velocity vector , , , These represent the components of velocity along the x-axis, y-axis, and z-axis, respectively.

[0070] The system establishes a local trajectory coordinate system based on the average velocity vector and generates an asymmetric three-dimensional Gaussian convolution kernel, specifically including:

[0071] With voxel center x k Based on the voxel center x k average velocity vector Transformation matrix R for constructing the local track coordinate system k The transformation matrix R k It is an orthogonal matrix that satisfies T is the matrix transpose. Specifically, the local track coordinate system uses the average velocity vector... The direction is the longitudinal axis, the direction perpendicular to the longitudinal axis and pointing to the right is the lateral axis, and the direction perpendicular to the horizontal plane containing the longitudinal axis and the lateral axis is upward is the vertical axis.

[0072] For any spatial point p, define spatial point p relative to the voxel center x. k The offset vector is x p -x k Then the spatial point p is relative to the voxel center x k The coordinates u of the offset vector in the local track coordinate system k Represented as:

[0073] ;

[0074] Among them, u x,k u y,k u z,k These represent spatial point p relative to voxel center x. k Longitudinal offset (along the track direction), lateral offset, vertical offset, x p This represents the three-dimensional spatial position of point p.

[0075] Furthermore, to reflect the asymmetry of the aircraft's forward and backward diffusion along the direction of motion, as well as the diffusion characteristics affected by the aircraft's speed and geometry, a spatial point p is defined relative to the voxel center x. k The longitudinal forward offset is Spatial point p relative to voxel center x k The longitudinal backward offset is Simultaneously, the standard deviation parameter for velocity amplitude coupling is set, including:

[0076] ;

[0077] ;

[0078] ;

[0079] ;

[0080] in, This represents a point p in space relative to the voxel center x. k Forward standard deviation; This represents a point p in space relative to the voxel center x. k The backward standard deviation; This represents a point p in space relative to the voxel center x. k Lateral standard deviation; This represents a point p in space relative to the voxel center x. k The vertical standard deviation; The basic diffusion parameter is the forward standard deviation; The basic diffusion parameter is the backward standard deviation; The basic diffusion parameter is the lateral standard deviation; The basic diffusion parameter is the vertical standard deviation; a f The coupling coefficient is the forward standard deviation; a b The coupling coefficient is the backward standard deviation; a y The coupling coefficient is the lateral standard deviation; a z Let a be the coupling coefficient of the vertical standard deviation, satisfying a f >a b ≥0; W represents the equivalent lateral dimension of the aircraft; H represents the vertical dimension of the aircraft; x represents k average velocity vector The magnitude, i.e., the velocity amplitude.

[0081] Based on the spatial point p relative to x k Longitudinal forward offset Longitudinal backward offset Lateral offset u y,k Vertical offset u z,kAnd by setting the standard deviation parameter for velocity amplitude coupling, the voxel center x is obtained. k Corresponding velocity amplitude Coupled asymmetric 3D Gaussian convolution kernel G k for:

[0082] ;

[0083] Among them, Z k This represents the normalization coefficient. Since the asymmetric 3D Gaussian convolution kernel is a piecewise 3D Gaussian convolution kernel with different standard deviations on the front and rear halves, therefore... This can be explicitly guaranteed that the integral of the asymmetric three-dimensional Gaussian convolution kernel function in three-dimensional space is 1.

[0084] Based on the asymmetric 3D Gaussian convolution kernel, a 3D probability density field is constructed, specifically including:

[0085] ;

[0086] Where x is any position in space; N represents the probability density of an aircraft appearing per unit volume at spatial location x; k This is the count of trajectory points within the k-th voxel. The discrete density matrix refers to the three-dimensional discrete array obtained by mapping the preprocessed trajectory sequence to each voxel element after completing the three-dimensional voxel mesh generation, denoted as... .in, This represents the count of trajectory points contained within the voxel cell with index (i,j,m); if no trajectory point falls within this voxel, the corresponding element value is 0. The discrete density matrix... This serves as the input for 3D smoothing operations using an asymmetric 3D Gaussian convolution kernel. The asymmetric 3D Gaussian convolution kernel is then convolved with the discrete density matrix to obtain a continuous 3D probability density field with directional attributes and adaptive velocity diffusion capabilities. The convolution boundary employs nearest-neighbor extrapolation to avoid underestimating the boundary density.

[0087] As shown in Figures 3(a) and 3(b), Figure 3(a) is a schematic diagram of the three-dimensional mesh effect in this invention, used to perform spatial statistics on the trajectory sequence; Figure 3(b) is a schematic diagram of the three-dimensional probability density field in this invention, where the counts of each voxel are diffused and fused along the local trajectory direction by an asymmetric three-dimensional Gaussian convolution kernel to obtain a continuous, smooth, and direction-preserving three-dimensional probability density field.

[0088] It should be noted that in S2, the construction of the three-dimensional probability density field is not a simple combination of average velocity vector statistics and Gaussian convolution. Instead, it first extracts spatial frequency and dominant flow direction simultaneously through tensor grid statistics, and then uses the dominant flow direction to drive the generation of local trajectory coordinate system and asymmetric three-dimensional Gaussian convolution kernel, ultimately forming a three-dimensional probability density field sensitive to local traffic flow direction. Compared to two-dimensional grid heatmaps that completely ignore height dimension information, and traditional three-dimensional trajectory point counting methods that produce "blocky" discontinuous distributions and are highly susceptible to data noise and uneven sampling frequency, the three-dimensional probability density field constructed in this invention combines spatial continuity, orientation preservation capability, and boundary stability.

[0089] It should be noted that the purpose of introducing the asymmetric three-dimensional Gaussian convolution kernel in this embodiment is that when an aircraft appears at a certain coordinate point, its impact on the airspace occupation and potential risk diffusion of the surrounding area is not isotropic, but is significantly constrained by the direction of motion. Therefore, by smoothing the discrete statistical results through the direction-adaptive asymmetric three-dimensional Gaussian convolution kernel, the spatial "holes" between tracks can be filled, and the interference of jumps in the original track data and fluctuations in sampling frequency on the accuracy of the model can be reduced.

[0090] Furthermore, unlike ordinary smoothing of discrete counting results, the above-mentioned spatial gridding method synchronously records spatial frequency and average velocity vector in voxels, so that the direction and diffusion range of subsequent convolution kernels can be jointly determined by the local dominant flow direction.

[0091] S3. For the preset route input by the user, perform a macro-strategic risk quantification assessment in the three-dimensional probability density field to generate the total route risk.

[0092] S3 specifically includes the following steps:

[0093] S31. Divide the preset route into multiple preset segments, and extract the risk density value of the starting point of each preset segment in the three-dimensional probability density field through nearest neighbor search.

[0094] A preset route refers to a flight path input by the user that requires risk assessment. Users input a preset route consisting of multiple waypoints through a web interface. Specifically, the format of a preset route is: [(lon1,lat1,alt1),(lon2,lat2,alt2),...(lon... N ,lat N ,alt N )], where lon i lat i alt i Let i represent the longitude, latitude, and altitude of the i-th waypoint, respectively, and N represent the total number of waypoints, i=1,2,...,N.

[0095] It should be noted that the preset route consists of a series of waypoints (i.e., the starting point of each preset segment), and thus the preset route can be regarded as a continuous spatial curve in three-dimensional space, denoted as C.

[0096] The system divides a preset route into N-1 preset segments based on waypoints. For a preset segment n (n=[1,2,...N-1]), the starting point P of the preset segment n is extracted through nearest neighbor search. n Risk density value in a three-dimensional probability density field :

[0097]

[0098] Among them, M smooth This represents a discrete array storing the smoothed three-dimensional probability density field, where each voxel unit corresponds to a risk density value; [x ndx ,y ndx ,z ndx ] indicates that the starting point P of the preset flight segment n is... n The three-dimensional coordinates are mapped to the voxel index.

[0099] The mapping method is nearest neighbor search, which means finding the starting point P that is n distances from the preset flight segment. n The nearest voxel center, and then the starting point P at a distance n from the preset flight segment. n The index of the most recent voxel center. Specifically, the starting point P of the preset flight segment n. n Match the nearest voxel center to the three-dimensional coordinates and extract the risk density value. .

[0100] S32. Calculate the flight distance for each preset flight segment.

[0101] The system calculates the flight distance d for the preset flight segment n. n The specific formulas include:

[0102] ;

[0103] Where, d horizontal ∆h represents the horizontal distance, calculated using the Haversine formula. n P represents the starting point of the preset flight segment n. n With the endpoint P n+1 The change in vertical height between them.

[0104] S33. Calculate the total route risk based on the risk density value, the flight distance, and the preset risk scaling constant.

[0105] Specifically, the system integrates the three-dimensional probability density field along the preset route C to obtain the total route risk R.total :

[0106] ;

[0107] Where ds is the differential of the three-dimensional arc length of the preset route C; k risk A preset risk scaling constant is used to calculate the original risk accumulation obtained through integration. Mapped to total route risk R at a uniform scale total .

[0108] It should be noted that the total risk R of the above-mentioned routes total The line integral formula is a theoretical expression; therefore, the specific total route risk needs to be obtained through approximate calculation. In a concrete and implementable approach, the system can use the Riemann sum to discretize and approximate the above line integral formula to obtain the total route risk R. total ,include:

[0109] Based on the discretized line integral formula, calculate the risk accumulation R for the preset flight segment n. n : .

[0110] The total route risk R of the preset route is obtained by linearly summing the cumulative risks of all preset segments. total :

[0111] .

[0112] As shown in Figures 4(a) and 4(b), Figure 4(a) is a schematic diagram of the superposition of the three-dimensional probability density field and the preset route in this invention. In the figure, the preset route (orange straight line) is superimposed on the three-dimensional probability density field. Figure 4(b) is a schematic diagram of the change of risk accumulation with flight distance in this invention. The total route risk of the entire preset route is obtained by summing the risk accumulation of all preset segments.

[0113] S4. Based on the relative motion relationship between the aircraft to be evaluated and the intruding aircraft, perform real-time dynamic collision risk quantification assessment on the aircraft to be evaluated and the intruding aircraft on the preset flight path, and calculate the peak collision risk index.

[0114] It should be noted that S4 is used for real-time microscopic dynamic collision risk assessment and does not rely on the three-dimensional probability density field constructed in S2.

[0115] Before conducting real-time dynamic collision risk quantification assessment, the system needs to perform vertical screening of aircraft to identify intruding aircraft and synchronize data simultaneously, including:

[0116] The extrapolated position of the aircraft to be evaluated and the surrounding traffic situation are obtained based on time series data.

[0117] The extrapolated position refers to the estimated position obtained by the system based on the latest trajectory data (including position, speed, acceleration, etc.) of the aircraft to be evaluated during the real-time dynamic collision risk assessment process. This is done in order to achieve synchronous analysis of the aircraft to be evaluated and surrounding aircraft under the same time reference. The system uses linear extrapolation or uniform acceleration motion model to extrapolate the position of the aircraft to be evaluated to the current unified analysis time. This is done because there are deviations or delays in the reporting time of ADS-B data of different aircraft.

[0118] Surrounding traffic situation refers to the real-time or near-real-time movement and spatial distribution of other aircraft surrounding the aircraft being evaluated within the airspace where the aircraft being evaluated is located.

[0119] Furthermore, the system sets a vertical height filtering threshold to eliminate irrelevant aircraft whose absolute height difference is greater than the vertical height filtering threshold, and regards the remaining surrounding aircraft as intruders.

[0120] In one embodiment, the vertical altitude filtering threshold can be set to 600 meters to quickly exclude aircraft with excessively large altitude differences from the aircraft to be evaluated that do not pose a spatial conflict possibility in a short period of time; in other embodiments, the vertical altitude filtering threshold can also be adjusted according to the type of operating airspace, aircraft flight performance, and risk sensitivity requirements.

[0121] The specific steps for the system to calculate the peak collision risk index include:

[0122] S41. Obtain the relative motion relationship between the aircraft to be evaluated and the intruding aircraft.

[0123] The aircraft to be evaluated refers to the aircraft that requires real-time dynamic collision risk quantification assessment.

[0124] Relative motion relationships include relative position vectors With relative velocity vector .

[0125] Define the relative position vectors of the aircraft to be evaluated and the intruding aircraft i at time t. for:

[0126] ;

[0127] Among them, P own (t) represents the three-dimensional spatial position of the aircraft to be evaluated at time t; p i (t) represents the three-dimensional spatial position of the intruder i at time t.

[0128] The system calculates the velocity vectors of the aircraft to be evaluated and the intruding aircraft i at time t, including:

[0129] ;

[0130] ;

[0131] The relative velocity vectors of the aircraft to be evaluated and the intruder i at time t are then evaluated. for:

[0132] ;

[0133] Among them, v own (t) represents the velocity vector of the aircraft to be evaluated at time t; Indicates the aircraft to be evaluated at time The three-dimensional spatial position; Indicates the aircraft to be evaluated at time The three-dimensional spatial position; v i (t) represents the velocity vector of the intruding machine i at time t; Indicates that the intrusion machine i is at time i The three-dimensional spatial position; Indicates that the intrusion machine i is at time i The three-dimensional spatial position; Indicates the time step.

[0134] S42. Establish a local coordinate system for relative motion based on the relative motion relationship between the aircraft to be evaluated and the intruding aircraft.

[0135] It should be noted that the local trajectory coordinate system in S2 is constructed based on the average velocity vector of the historical trajectory points within the voxel and fixed at the voxel center, and is used for the directional adaptive diffusion of the asymmetric three-dimensional Gaussian convolution kernel; while the relative motion local coordinate system in S42 is constructed based on the real-time relative motion dynamics between the aircraft to be evaluated and the intruding aircraft and is updated with the aircraft's motion, and is used for real-time dynamic collision risk assessment coupled with the three-dimensional collision box model and the Laplacian error.

[0136] The system establishes a local coordinate system based on the relative motion relationship and utilizes the coordinate transformation matrix R. i (t) The relative position vector Projecting these coordinates onto the local coordinate system of relative motion, we obtain the longitudinal, lateral, and vertical distances between the aircraft to be evaluated and the intruding aircraft i:

[0137] ;

[0138] in, This indicates the longitudinal distance between the aircraft to be evaluated and the intruding aircraft i; This indicates the lateral distance between the aircraft to be evaluated and the intruding aircraft i; This indicates the vertical distance between the aircraft to be evaluated and the intruding aircraft i.

[0139] S43. Combining the three-dimensional collision box model with the Laplace position error distribution, a real-time dynamic collision risk quantification assessment is performed on the aircraft to be evaluated and the intruding aircraft on the preset flight path, and the peak collision risk index is calculated.

[0140] It should be noted that the 3D collision box model can treat the aircraft to be evaluated and the intruding aircraft as cubes with certain geometric dimensions. Therefore, the collision conflict domain of the aircraft to be evaluated and the intruding aircraft can be equivalent to a joint 3D rectangular box, which can be used to quantify the critical conditions for the spatial overlap of the aircraft to be evaluated and the intruding aircraft.

[0141] The system combines a 3D collision box model with the Laplace position error distribution to perform real-time dynamic collision risk quantification assessment of the aircraft to be evaluated and the intruding aircraft on a preset flight path, and calculates the peak collision risk index, specifically including:

[0142] Define the equivalent length, width, and height of the aircraft to be evaluated as L. own W own H own The equivalent length, width, and height of the intrusion machine i are L, ... i W i H i The joint conflict domain of the collision between the aircraft to be evaluated and the intruding aircraft i can be equivalent to a three-dimensional collision box model, wherein half-length b of the three-dimensional collision box model is... x,i Half-width b y,i With half-height b z,i They are respectively:

[0143] ;

[0144] The Laplace position error distribution is a heavy-tailed probability distribution used to describe the positioning error of an aircraft.

[0145] Suppose that the longitudinal, lateral, and vertical positions each follow a Laplace distribution, then the cumulative distribution function of the Laplace position error distribution is... for:

[0146]

[0147] in, For a random variable to take values, For position parameters, Let be the scaling parameter, and exp be an exponential function with base e. In one embodiment, assuming the navigation error is symmetrically distributed around zero, the position parameter... Scale parameters The degree of dispersion of position error can be statistically estimated from historical navigation and positioning error samples, or set according to the accuracy of the navigation system used.

[0148] Furthermore, the longitudinal probability P of the Laplace position error distribution of the intruding machine i falling into the interval of the three-dimensional collision box model at time t is... x,i (t), lateral probability P y,i (t) and vertical probability P z,i (t) are respectively:

[0149] ;

[0150] ;

[0151] ;

[0152] Among them, F x F y F z These are the cumulative distribution functions of the Laplace distributions followed by the longitudinal, lateral, and vertical position errors, respectively.

[0153] Therefore, the instantaneous collision probability P at time t col,1 (t) is:

[0154] .

[0155] Furthermore, to avoid calculation errors caused by relying solely on the extreme values ​​at a single sampling moment, this embodiment introduces a preset short time window. Define the cumulative collision probability P of the intruding machine i within a preset short time window starting at time t. win,t (t) is:

[0156] ;

[0157] in, Indicates the index of the discrete time step, with values ​​ranging from time t to... Each sampling time between; Indicates a preset short time window; Indicates that the intrusion machine i is at time i The probability that no collision occurs; This represents the approximate probability that no collision occurs at any of the sampling times within the preset short time window, under the condition of a preset short time window discrete approximation, when the absence of collision at each time moment is considered a conditionally independent event.

[0158] Ultimately, the system outputs the peak collision risk index P. peak for:

[0159] ;

[0160] in, The peak collision risk index represents the cumulative peak collision probability of all intruding machines at each time t within the assessment period. This means taking the maximum value for the time dimension and iterating through all time indices t; This means taking the maximum value for the intrusion machine dimension and iterating through all monitored intrusion machine indices i.

[0161] It should be noted that S3 is used for microscopic real-time dynamic collision risk assessment. It does not simply take the maximum value of the collision probability at a single moment. Instead, under the constraint of relative motion relationship, it couples the Laplace position error distribution with the collision domain of the three-dimensional collision box model, and further constructs the peak collision risk index through a preset short time window accumulation mechanism.

[0162] S5. Output a risk report based on the total route risk and the peak collision risk index, and display it visually.

[0163] In one embodiment, the system uses a 3D voxel rendering method to display the 3D probability density field, specifically including: projecting the smoothed voxel density matrix into a set of voxel surfaces according to a second preset threshold or continuous color mapping, and providing multi-view rendering such as top view, front view, back view, left view, and right view; simultaneously overlaying and displaying a preset flight path and waypoint markers in the 3D view. The second preset threshold is obtained based on prior experience.

[0164] The risk report outputs the start and end points, segment distances, and the starting point P of the preset route n in tabular form for each preset segment. n Risk density value in a three-dimensional probability density field The cumulative risk R of the preset flight segment n n The system then summarizes and outputs the total risk of the flight route. It can further map the total risk to risk levels based on preset classification thresholds and display these levels on the interface using color or icons. These preset classification thresholds can be obtained based on prior experience.

[0165] In another embodiment, the visualization front-end can use WebGL / Three.js to implement volume drawing or point cloud rendering, and the back-end provides density field and risk data through a REST API.

[0166] like Figure 5 As shown in the diagram, the system visualization interface illustrates the multi-view rendering of the three-dimensional probability density field and the comprehensive display of risk assessment results. Figure 5The central area uses a 3D voxel rendering method to present a smoothed 3D probability density field. Different voxels are colored according to their density values ​​using continuous color mapping (such as red-yellow-green gradient) or a second preset threshold, which intuitively reflects the airspace risk distribution. A multi-view switching button is provided in the upper left corner of the interface, which makes it easy for users to observe the spatial structure of the density field from different angles. The view also displays the preset route polyline trajectory and its waypoint markers entered by the user. The risk report is output in tabular form on the right side of the interface, including the start and end points of each segment, segment distance, density sampling value, segment risk value, and a summary of the total route risk. The total risk is mapped to a risk level (such as low / medium / high) according to a preset classification threshold and indicated by color icons.

[0167] Example 2:

[0168] This embodiment is a further optimization based on Embodiment 1 described above, such as... Figure 6 As shown, an asymmetric modeling airspace risk assessment system is provided. The system includes: a preprocessing module, a density field construction module, a strategic risk assessment module, a collision risk assessment module, and a visualization output module.

[0169] The preprocessing module is used to collect raw track data and preprocess the raw track data to obtain a track sequence.

[0170] The density field construction module uses a three-dimensional mesh to perform spatial statistics on the trajectory sequence and records the average velocity vector of trajectory points within each voxel. A local trajectory coordinate system is established based on the average velocity vector, and an asymmetric three-dimensional Gaussian convolution kernel is generated. A three-dimensional probability density field is constructed through the asymmetric three-dimensional Gaussian convolution kernel.

[0171] The strategic risk assessment module is used to perform a macro-level strategic risk quantification assessment on the three-dimensional probability density field for the preset route input by the user, and generate the total route risk.

[0172] The collision risk assessment module is used to perform real-time dynamic collision risk quantification assessment of the aircraft to be assessed and the intruding aircraft on a preset flight path based on the relative motion relationship between the aircraft to be assessed and the intruding aircraft, and to calculate the peak collision risk index.

[0173] The visualization output module is used to output a risk report based on the total route risk and the peak collision risk index, and to visualize the report.

[0174] The other parts of this embodiment are the same as those in the above embodiments, so they will not be described again.

[0175] The above description is merely a specific embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the technical scope disclosed in the present invention should be included within the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be determined by the scope of the claims.

Claims

1. An asymmetric modeling spatial risk assessment method characterized by, The method comprises: S1, collecting original flight path data, preprocessing the original flight path data to obtain a trajectory sequence; S2, using a three-dimensional grid to statistically analyze the trajectory sequence, and recording the average velocity vector of the trajectory points in each voxel; establishing a local flight path coordinate system according to the average velocity vector, generating an asymmetric three-dimensional Gaussian convolution kernel, and constructing a three-dimensional probability density field through the asymmetric three-dimensional Gaussian convolution kernel; S3, for a preset flight path input by a user, quantitatively evaluating the macroscopic strategic risk of the three-dimensional probability density field, and generating a total flight path risk; S4, based on the relative motion relationship between the aircraft to be evaluated and the intruder, quantitatively evaluating the real-time dynamic collision risk of the aircraft to be evaluated and the intruder on the preset flight path, and calculating a peak collision risk index; S5, outputting a risk report according to the total flight path risk and the peak collision risk index, and visually displaying the risk report.

2. The asymmetric modeling spatial risk assessment method of claim 1, wherein, The preprocessing in S1 comprises: screening and removing abnormal trajectory points based on a first preset threshold value; the first preset threshold value comprises a latitude jump threshold value and a height jump threshold value; grouping according to the aircraft call sign, segmenting the trajectory with a time interval greater than 30 minutes, and performing linear interpolation at an interval of 1 second.

3. The asymmetric modeling spatial risk assessment method of claim 1, wherein, In S2, the three-dimensional probability density field is represented as: ； wherein, represents the probability density of the aircraft appearance in a unit volume at spatial position x; k is the voxel index; N k is the trajectory point count within voxel k; x k is the center of the kth voxel, is the transpose matrix of the local trajectory coordinate transformation matrix R k ; is the asymmetric three-dimensional Gaussian convolution kernel corresponding to voxel k.

4. The asymmetric modeling spatial risk assessment method of claim 1 or 3, wherein, The asymmetric three-dimensional Gaussian convolution kernel is represented as: ； wherein, denotes the longitudinal forward offset of the spatial point p relative to the voxel center x k ; denotes the longitudinal backward offset of the spatial point p relative to the voxel center x k ; y,k denotes the lateral offset of the spatial point p relative to the voxel center x k ; z,k denotes the vertical offset of the spatial point p relative to the voxel center x k ; denotes the forward standard deviation of the spatial point p relative to the voxel center x k ; denotes the backward standard deviation of the spatial point p relative to the voxel center x k ; denotes the lateral standard deviation of the spatial point p relative to the voxel center x k ; denotes the vertical standard deviation of the spatial point p relative to the voxel center x k ; k denotes a normalization coefficient.

5. The asymmetric modeling spatial risk assessment method of claim 1, wherein, In S3, the macroscopic strategic risk quantitative evaluation comprises: S31, dividing the preset flight path into multiple preset flight segments, and extracting the risk density value of each preset flight segment in the three-dimensional probability density field through nearest neighbor search; S32, calculating the flight distance of each preset flight segment; S33, calculating the total flight path risk based on the risk density value, the flight distance and a preset risk scaling constant.

6. The asymmetric modeling spatial risk assessment method of claim 1, wherein, In S4, the real-time dynamic collision risk quantitative evaluation specifically comprises: S41, obtaining the relative motion relationship between the aircraft to be evaluated and the intruder; S42, establishing a relative motion local coordinate system based on the relative motion relationship between the aircraft to be evaluated and the intruder; S43, combining a three-dimensional collision box model and a Laplace position error distribution to quantitatively evaluate the real-time dynamic collision risk of the aircraft to be evaluated and the intruder on the preset path, and calculate the peak collision risk index.

7. The asymmetric modeling spatial risk assessment method of claim 1, wherein, The cumulative distribution function of the Laplace position error distribution is: wherein is a random variable, is a location parameter, is a scale parameter.

8. An asymmetrically modelled airspace risk assessment system for implementing the asymmetrically modelled airspace risk assessment method of any one of claims 1 to 7, characterised in that, The system comprises a preprocessing module, a density field construction module, a collision risk evaluation module, a strategic risk evaluation module and a visualization output module: The preprocessing module is used to collect original flight path data, preprocess the original flight path data to obtain a trajectory sequence; The density field construction module uses a three-dimensional grid to statistically analyze the trajectory sequence, and records the average velocity vector of the trajectory points in each voxel; a local flight path coordinate system is established according to the average velocity vector, an asymmetric three-dimensional Gaussian convolution kernel is generated, and a three-dimensional probability density field is constructed through the asymmetric three-dimensional Gaussian convolution kernel; The strategic risk evaluation module is used to quantitatively evaluate the macroscopic strategic risk of the three-dimensional probability density field for a preset flight path input by a user, and generate a total flight path risk; The collision risk assessment module is configured to establish a relative motion local coordinate system based on the relative motion relationship between the aircraft to be evaluated and the intruder aircraft, and to combine a three-dimensional collision box model and a Laplace position error distribution to perform real-time dynamic collision risk quantitative evaluation on the aircraft to be evaluated and the intruder aircraft on the preset flight path, and to calculate a peak collision risk index. The visualization output module is configured to output a risk report according to the total risk of the flight path and the peak collision risk index, and to perform visual display.

Images