A non-cooperative UAV flight intention recognition method based on similarity measurement

Abstract

The present invention discloses a non-cooperative UAV flight intention recognition method based on similarity measurement, including: collecting historical trajectory data of multiple types of UAVs to obtain a UAV historical trajectory slice sample set; performing similarity measurement on the trajectory slice sample set; clustering the slice sample set using a clustering algorithm to obtain a clustering result; extracting flight intention features from each cluster of UAV trajectory slice samples; training them separately based on the GRU algorithm to obtain a UAV flight intention recognition model set; and obtaining a UAV flight intention recognition result. The present invention performs similarity measurement and clustering based on the flight behavior characteristics of the UAV, screens out training samples with similar flight characteristics to the non-cooperative UAV trajectory fragments, uses different levels of maneuvering behaviors to characterize flight intentions, improves the recognition accuracy and timeliness of the UAV flight intention recognition model, and enhances its applicability in real-time flight intention recognition tasks.

Description

Technical Field

[0001] The present invention belongs to the technical field of unmanned aerial vehicle traffic management, and in particular relates to a non-cooperative unmanned aerial vehicle flight intention recognition method based on similarity measurement. Background Art

[0002] With the development of new technologies such as 5G communications, unmanned driving, and automatic control, low-altitude drones are developing rapidly, flight demand is increasing day by day, air traffic activities are becoming more frequent, and low-altitude airspace logistics distribution, security patrols, firefighting and disaster relief and other general aviation businesses are becoming increasingly diversified. While low-altitude airspace traffic operations are constantly developing, drone abuse and various unsafe flight behaviors have also emerged. In the face of these dangerous behaviors of non-cooperative drones, if the flight intentions of non-cooperative drones can be identified in real time, it will lay the foundation for drone operation safety control and risk assessment.

[0003] Non-cooperative UAVs have no communication with surrounding aircraft or objects, and cannot share their status information. Their flight intentions are often complex and changeable. When identifying the flight intentions of non-cooperative UAVs, existing AI-based target intention recognition methods often have difficulty obtaining reliable features. The UAV flight intention recognition model based on data mining and machine learning has higher accuracy than the AI-based method. However, due to the small amount of historical sample data of non-cooperative UAVs, it is necessary to solve how to achieve high-precision task recognition in the absence of real historical sample data. Summary of the invention

[0004] In view of the above-mentioned deficiencies in the prior art, the purpose of the present invention is to provide a non-cooperative UAV flight intention recognition method based on similarity measurement, perform similarity measurement and clustering according to the flight behavior characteristics of the UAV, screen out training samples with similar flight characteristics to the non-cooperative UAV trajectory fragments, and use different levels of maneuvering behaviors to characterize flight intentions, thereby improving the recognition accuracy and timeliness of the UAV flight intention recognition model and enhancing its applicability in real-time flight intention recognition tasks.

[0005] To achieve the above object, the technical solution adopted by the present invention is as follows:

[0006] The present invention provides a non-cooperative UAV flight intention recognition method based on similarity measurement, and the steps are as follows:

[0007] 1) Collect historical trajectory data of various types of drones, process the data, and obtain a sample set of drone historical trajectory slices;

[0008] 2) performing similarity measurement on the trajectory slice sample set of the historical trajectory of the UAV obtained in the above step 1);

[0009] 3) Based on the similarity measurement result in step 2), a clustering algorithm is used to cluster the sample set of the drone historical trajectory slices to obtain a clustering result;

[0010] 4) Based on the clustering results in step 3), extract the flight intention features of the UAV trajectory slice samples of each cluster and set the model input and output;

[0011] 5) Based on the flight intention features extracted from the trajectory slice samples of each cluster of drones, they are trained separately based on the GRU algorithm to obtain a set of drone flight intention recognition models;

[0012] 6) Perform intention recognition based on the UAV flight intention recognition model set to obtain the UAV flight intention recognition result.

[0013] Furthermore, in the step 1), historical trajectory data of various types of drones are collected at a sampling frequency of 10 Hz, and the timestamp, identity ID, longitude, latitude, relative flight altitude, speed and attitude angle information of each trajectory point are recorded.

[0014] Furthermore, the data processing in step 1) includes: coordinate conversion, equally spaced sequence processing, trajectory slicing division and screening, and dimensional processing.

[0015] Furthermore, the specific process of step 2) is as follows:

[0016] 21) Select similarity measurement indicators;

[0017] The trajectory slice samples are similarity measured, the dimension is [T in ,f],T in is the time length of the trajectory slice sample, and f represents the feature number of the trajectory slice sample. Based on the motion characteristics of the UAV, it is characterized by the changes in flight speed and turning angle, and similarity measurement indicators are selected, including speed change size, speed change trend, turning angle cosine similarity, and turning angle change amplitude.

[0018] 22) Similarity measurement of speed change;

[0019] The Pearson correlation coefficient is used to represent the similarity of the speed between slices. The Pearson correlation coefficient r between trajectory slice sample i and trajectory slice sample j is ij as follows:

[0020]

[0021] In the formula, s ie 、s je They represent the velocity values ​​corresponding to the trajectory slice samples i and j at the e-th time point, n represents the number of slice samples, l is the slice length of the trajectory slice sample, μi , μ j are the mean speeds of sample i and sample j respectively;

[0022] 23) Similarity measurement of speed change trend;

[0023] Calculate the product of the speed mean ratio and variance ratio of each two trajectory slice samples, and the speed change trend similarity RVC between trajectory slice sample i and trajectory slice sample j ij as follows:

[0024]

[0025] In the formula, Respectively represent the velocity mean of the i-th and j-th trajectory slice samples, σ i , σ j Respectively represent the velocity variance of the i-th and j-th trajectory slice samples;

[0026] 24) Cosine similarity measurement of turning angle;

[0027] Calculate the cosine similarity of the turning angles of the two trajectory slice samples in the corresponding time. According to the turning angle vector pair between each two time points of the slice sample, calculate the cosine similarity of the turning angles between the corresponding vector pairs of the trajectory. k , and then get the cosine sequence {cos k}, k = 1, 2, ..., l-1, l is the slice length of the trajectory slice sample; the cosine similarity of the turning angle between the trajectory slice sample i and the trajectory slice sample j is cos ij-k as follows:

[0028]

[0029] In the formula, a ik-p 、a jk-p Respectively represent the pth component of the turning angle vector of the i-th and j-th trajectory slice samples in the corresponding time interval k, q represents the number of elements in the turning angle vector coordinates, and the cosine similarity is solved for every two trajectory slice samples to obtain a sequence containing l-1 cosine values. The average value of the sequence, cos ij Perform similarity representation using cosine similarity;

[0030] 25) Similarity measurement of turning angle change amplitude;

[0031] Calculate the position of the maximum difference in turning angle between two trajectory slice samples The ratio of the slice length l of the trajectory slice sample and the initial turning angle difference A of the two trajectory slice samples ij-initial The maximum difference Aij-max The product of the two is used as the similarity R of the turning angle change between the two trajectory slice samples. ij , the similarity R of the turning angle change between trajectory slice sample i and trajectory slice sample j ij as follows:

[0032]

[0033] In the formula, represents the distance from the starting point to the point where the maximum difference in turning angle between trajectory slice sample i and trajectory slice sample j occurs, l is the slice length of the trajectory slice sample, and A ij-initial Represents the initial turning angle difference between two trajectory slice samples, A ij-max is the maximum turning angle difference between two trajectory slice samples;

[0034] 26) Calculate the similarity measure matrix;

[0035] For the drone trajectory slice sample set obtained in step 1), the similarity measurement values ​​between any two trajectory slice samples are calculated according to steps 22) to 25) to obtain the corresponding four similarity measurement matrices.

[0036] Furthermore, the specific process of step 3) is as follows:

[0037] 31) Calculate the comprehensive similarity measure matrix;

[0038] The four similarity measurement matrices obtained in step 26) are weighted averaged according to the average weight η=0.25 to obtain a comprehensive similarity measurement matrix M. The specific calculation formula is as follows:

[0039]

[0040] In the formula, s ij =η×r ij +η×RVC ij +η×cos ij +η×R ij ,s ij represents the similarity value between trajectory slice sample i and trajectory slice sample j in M;

[0041] 32) Calculate the shared similarity samples between any two trajectory slice samples in the comprehensive similarity measurement matrix;

[0042] The similarity measurement threshold is used to screen shared similar samples. The specific calculation formula is as follows:

[0043] SNN(x i ,x j )=τ(xi )∩τ(x j ) (6)

[0044]

[0045] In the formula, τ(x i ) represents the set of similar samples of sample i under the specified threshold; τ(x j ) represents the set of similar samples of sample j under the specified threshold; SNN(x i ,x j ) is given by τ(x i ) and τ(x j ) is obtained by taking the intersection, indicating that samples i and j share a similar sample set; Denotes the set SNN(x i ,x j ) The number of elements in s ip represents the similarity between sample i and sample p in the comprehensive similarity measurement matrix; s jp Represents the similarity between sample j and sample p in the comprehensive similarity measurement matrix; Sim(x i ,x j ) represents the SNN of sample i and sample j (x i ,x j ) Similarity value, which is calculated only when sample i and sample j exist in each other's shared similar samples, otherwise it is zero;

[0046] 33) Calculate the local similarity sample density ρ i , the specific calculation formula is as follows:

[0047]

[0048] In the formula, L(x i ) represents the set of sample points that meet the similarity threshold requirement with sample i under the similarity measurement threshold, Sim(x i ,x j ) indicates that sample i and sample j share similar sample SNN (x i ,x j ) similarity;

[0049] 34) Calculate relative similarity δ i ;

[0050] Calculate the relative similarity δ i , which represents the maximum similarity between sample i and its nearby high-density points. The specific calculation formula is as follows:

[0051]

[0052] 35) Calculate the decision value;

[0053] Use the decision value γ for each sample point xi i To express the possibility of becoming a cluster center, the specific calculation is as follows:

[0054] γ i =ρ i *δ i (10)

[0055] Among them, the decision values ​​γ are sorted in ascending order. i Arrange and get γ i The sequence γ' arranged from small to large;

[0056] 36) Calculate the inflection point position and obtain the number of clusters and cluster centers;

[0057] The cluster center is selected on the sample point with a large decision value. The specific calculation formula is as follows:

[0058]

[0059]

[0060]

[0061] μ i =γ' i+1 -γ' i ,i=n-DN+1,n-DN+2,...,n-1 (14)

[0062] Where DN is the nearest integer value after taking the square root of the number of sample points n; γ' i represents the i-th decision value in the sequence γ', γ' i+1 represents the i+1th decision value in the sequence γ'; μ i Represents the decision value γ' i+1 With γ' i The difference between i+1 Represents the decision value γ' i+2 With γ' i+1 The difference between them; kp is the point with the largest ranking value, which satisfies its corresponding μ i Exceeding a certain threshold, with DN as the starting boundary, the threshold θ is μ in the interval [DN,n] i The average value of

[0063] 37) Determine the clustering results;

[0064] According to the cluster centers obtained in step 36), based on the comprehensive similarity measurement matrix, the non-cluster center drone trajectory slice samples are respectively classified into the cluster center with the highest similarity to obtain the clustering result.

[0065] Furthermore, the specific process of step 4) is as follows:

[0066] 41) Extract flight intention features;

[0067] The flight intention features are represented as the UAV's four horizontal flight behaviors, three longitudinal flight behaviors, and speed changes: the horizontal flight behavior includes the heading angle interval between two adjacent moments, which is divided into four intervals at 90° intervals; the longitudinal flight behavior is the pitch angle interval between two adjacent moments, including descent, climb, and level flight; the speed change consists of three situations: acceleration, constant speed, and deceleration between two adjacent moments;

[0068] 42) Calculate flight intention characteristic parameters;

[0069] The first trajectory point of the trajectory slice sample is marked as the origin (0,0,0), then the position coordinates of the i-th trajectory point in the slice sample are expressed as vector P i , P i =(x i ,y i ,z i ), and then the sequence input features are obtained, which are the horizontal heading angle A of the drone between each two adjacent moments in the trajectory slice sample. hi , the vertical pitch angle A vi and acceleration a i , the specific calculation formula is as follows:

[0070]

[0071]

[0072]

[0073] Where Δt is the sampling interval in seconds; |*| is the modulus of the vector; V i+1 and V i are the velocities of the drone at the i+1th and ith trajectory points, respectively, i =(x i ,y i ,z i ) is the position coordinate of the i-th trajectory point in the trajectory slice sample, P i+1 =(x i+1 ,y i+1 ,z i+1 ) is the position coordinate of the i+1th trajectory point in the trajectory slice sample, H i+1 is the vector P i+1 -P i Projection on the horizontal plane, H i+1=(x i+1 -x i ,y i+1 -y i ,0);L i+1 P i+1 -P i The projection in the latitude direction, L i+1 =(0,y i+1 -y i ,0);

[0074] 43) Set input sample;

[0075] The flight intention sequence information is represented as the change process between two adjacent moments, with a time step of T in -1, the flight intention recognition input sample is [b,T inten ], b is the number of batch samples, T inten is the time series intention information of the UAV, and the specific calculation formula is as follows:

[0076]

[0077] 44) Processing output data;

[0078] For the output of flight intention recognition, the behavioral feature values ​​of horizontal flight behavior, longitudinal flight behavior and speed change in step 41) are mapped to a unified interval range, and the output processing of the category label is selected through one-hot encoding. The flight intention of each state quantity at each moment is one-hot encoded as the input of the recognition model, and the category information is retained and converted into a numerical form that can be learned and processed by the algorithm.

[0079] Furthermore, the specific process of step 5) is as follows:

[0080] 51) Determine a training sample set;

[0081] Determine the hybrid trajectory data slice set S, S = {s 1 ,s 2 ,s 3 ,···,s r}, mixed flight intention data slice set D, D = {d 1 ,d 2 ,d 3 ,···,d r}, and a neural network model Q consisting of an input layer, three parallel GRU layers, a parallel fully connected layer, and an output layer, where S and D are divided into training set S in the ratio of 60%, 20%, and 20% respectively. t , D t , validation set S v , Dv , test set S c , D c ,The data inclusion based on the training sample set is as follows;

[0082] Scenario 1: A non-cooperative drone type can be detected, but there is no real historical trajectory of the corresponding type of drone in the training sample set;

[0083] Scenario 2: The type of non-cooperative drone cannot be detected, but the training sample set contains some real historical trajectories of the corresponding type of drone;

[0084] 52) Calculate the slice comprehensive similarity measure matrix;

[0085] Based on the mixed trajectory data slice set S t , according to the similarity measurement method of step 2), the comprehensive similarity measurement matrix M in step 31) is obtained;

[0086] 53) Calculate the number of clusters and cluster centers;

[0087] According to the clustering method in step 3), the number of clusters m and the cluster center c are obtained. i ,i=1,2,...,m;

[0088] 54) Divide the flight intention slice sample set;

[0089] According to the cluster center c i and the similarity between each trajectory sample slice and the cluster center, and the non-cluster center UAV trajectory slice samples are respectively classified into the cluster center with the highest similarity, and then the mixed flight intention data slice set D t Flight intention slice sample sets divided into different categories

[0090] 55) Obtaining a set of UAV flight intention recognition models;

[0091] Slice the mixed flight intention data into a set D t The sample subset for each category in As the input of model Q, the GRU algorithm is used for training to obtain the sample subsets Flight intention recognition model under training i And save it, and then get the UAV flight intention recognition model set O,

[0092] Furthermore, the specific process of step 6) is as follows;

[0093] 61) Determine the corresponding flight intention recognition model of the test sample;

[0094] According to the similarity measurement method in step 2), calculate the test set D c The test samples and the cluster centers c i The similarity value is combined with the flight intention recognition model set O saved in step 55), and the flight intention recognition model O corresponding to the maximum similarity value is selected. i The flight intention recognition model as a test sample;

[0095] 62) Obtaining intent recognition results;

[0096] Input the test sample into the determined flight intention recognition model O i , the flight intention recognition probabilities of the UAV’s horizontal, vertical and speed changes are obtained respectively, and the maneuvering intention features corresponding to the maximum probability are taken as the results of intention recognition.

[0097] Beneficial effects of the present invention:

[0098] The method of the present invention improves the data-driven UAV flight intention recognition task in tasks lacking real training samples, lays a foundation for flight conflict and risk assessment in non-cooperative UAV intrusion scenarios, and provides effective decision-making support for low-altitude airspace UAV operation planning and control.

[0099] The method of the present invention increases the output of the horizontal, vertical and acceleration maneuvering flight intention probabilities of non-cooperative targets at future prediction moments by dividing different maneuvering flight intention intervals in the data-driven trajectory prediction model, and further expresses the uncertainty of UAV flight. BRIEF DESCRIPTION OF THE DRAWINGS

[0100] Figure 1 Flow chart of the method of the present invention. DETAILED DESCRIPTION

[0101] In order to facilitate the understanding of those skilled in the art, the present invention is further described below in conjunction with embodiments and drawings. The contents mentioned in the implementation modes are not intended to limit the present invention.

[0102] Reference Figure 1 As shown, a non-cooperative UAV flight intention recognition method based on similarity measurement of the present invention comprises the following steps:

[0103] 1) Collect historical trajectory data of various types of drones, process the data, and obtain a sample set of drone historical trajectory slices;

[0104] The historical trajectory data of various types of drones are collected at a sampling frequency of 10HZ, and the timestamp, ID, longitude, latitude, relative flight altitude, speed and attitude angle information of each trajectory point are recorded.

[0105] Data processing includes: coordinate transformation, equally spaced sequence processing, trajectory slicing and screening, and dimensional processing.

[0106] 2) Perform similarity measurement on the trajectory slice sample set of the historical trajectory of the drone obtained in the above step 1); the specific process is as follows:

[0107] 21) Select similarity measurement indicators;

[0108] The trajectory slice samples are similarity measured, the dimension is [T in ,f],T in is the time length of the trajectory slice sample, and f represents the feature number of the trajectory slice sample. Based on the motion characteristics of the UAV, it is characterized by the changes in flight speed and turning angle, and similarity measurement indicators are selected, including speed change size, speed change trend, turning angle cosine similarity, and turning angle change amplitude.

[0109] 22) Similarity measurement of speed change;

[0110] The Pearson correlation coefficient is used to represent the similarity of the speed between slices. The Pearson correlation coefficient r between trajectory slice sample i and trajectory slice sample j is ij as follows:

[0111]

[0112] In the formula, s ie 、s je They represent the velocity values ​​corresponding to the trajectory slice samples i and j at the e-th time point, n represents the number of slice samples, l is the slice length of the trajectory slice sample, μ i , μ j are the mean speeds of sample i and sample j respectively;

[0113] 23) Similarity measurement of speed change trend;

[0114] Calculate the product of the speed mean ratio and variance ratio of each two trajectory slice samples, and the speed change trend similarity RVC between trajectory slice sample i and trajectory slice sample j ij as follows:

[0115]

[0116] In the formula, Respectively represent the velocity mean of the i-th and j-th trajectory slice samples, σ i , σ j Respectively represent the velocity variance of the i-th and j-th trajectory slice samples;

[0117] 24) Cosine similarity measurement of turning angle;

[0118] Calculate the cosine similarity of the turning angles of the two trajectory slice samples in the corresponding time. According to the turning angle vector pair between each two time points of the slice sample, calculate the cosine similarity of the turning angles between the corresponding vector pairs of the trajectory. k , and then get the cosine sequence {cos k}, k = 1, 2, ..., l-1, l is the slice length of the trajectory slice sample; the cosine similarity of the turning angle between the trajectory slice sample i and the trajectory slice sample j is cos ij-k as follows:

[0119]

[0120] In the formula, a ik-p 、a jk-p Respectively represent the pth component of the turning angle vector of the i-th and j-th trajectory slice samples in the corresponding time interval k, q represents the number of elements in the turning angle vector coordinates, and the cosine similarity is solved for every two trajectory slice samples to obtain a sequence containing l-1 cosine values. The average value of the sequence, cos ij Perform similarity representation using cosine similarity;

[0121] 25) Similarity measurement of turning angle change amplitude;

[0122] Calculate the position of the maximum difference in turning angle between two trajectory slice samples The ratio of the slice length l of the trajectory slice sample and the initial turning angle difference A of the two trajectory slice samples ij-initial The maximum difference A ij-max The product of the two is used as the similarity R of the turning angle change between the two trajectory slice samples. ij , the similarity R of the turning angle change between trajectory slice sample i and trajectory slice sample j ij as follows:

[0123]

[0124] In the formula, represents the distance from the starting point to the point where the maximum difference in turning angle between trajectory slice sample i and trajectory slice sample j occurs, l is the slice length of the trajectory slice sample, and A ij-initial Represents the initial turning angle difference between two trajectory slice samples, A ij-max is the maximum turning angle difference between two trajectory slice samples;

[0125] 26) Calculate the similarity measure matrix;

[0126] For the drone trajectory slice sample set obtained in step 1), the similarity measurement values ​​between any two trajectory slice samples are calculated according to steps 22) to 25) to obtain the corresponding four similarity measurement matrices.

[0127] 3) Based on the similarity measurement results in step 2), a clustering algorithm is used to cluster the sample set of drone historical trajectory slices to obtain the clustering results; the specific process is as follows:

[0128] 31) Calculate the comprehensive similarity measure matrix;

[0129] The four similarity measurement matrices obtained in step 26) are weighted averaged according to the average weight η=0.25 to obtain a comprehensive similarity measurement matrix M. The specific calculation formula is as follows:

[0130]

[0131] In the formula, s ij =η×r ij +η×RVC ij +η×cos ij +η×R ij ,s ij represents the similarity value between trajectory slice sample i and trajectory slice sample j in M;

[0132] 32) Calculate the shared similarity samples between any two trajectory slice samples in the comprehensive similarity measurement matrix;

[0133] The similarity measurement threshold is used to screen shared similar samples. The specific calculation formula is as follows:

[0134] SNN(x i ,x j )=τ(x i )∩τ(x j ) (6)

[0135]

[0136] In the formula, τ(x i ) represents the set of similar samples of sample i under the specified threshold; τ(x j ) represents the set of similar samples of sample j under the specified threshold; SNN(x i ,x j ) is given by τ(x i ) and τ(x j ) is obtained by taking the intersection, indicating that samples i and j share a similar sample set; Denotes the set SNN(x i ,x j ) The number of elements in s iprepresents the similarity between sample i and sample p in the comprehensive similarity measurement matrix; s jp Represents the similarity between sample j and sample p in the comprehensive similarity measurement matrix; Sim(x i ,x j ) represents the SNN of sample i and sample j (x i ,x j ) Similarity value, which is calculated only when sample i and sample j exist in each other's shared similar samples, otherwise it is zero;

[0137] 33) Calculate the local similarity sample density ρ i , the specific calculation formula is as follows:

[0138]

[0139] In the formula, L(x i ) represents the set of sample points that meet the similarity threshold requirement with sample i under the similarity measurement threshold, Sim(x i ,x j ) indicates that sample i and sample j share similar sample SNN (x i ,x j ) similarity;

[0140] 34) Calculate relative similarity δ i ;

[0141] Calculate the relative similarity δ i , which represents the maximum similarity between sample i and its nearby high-density points. The specific calculation formula is as follows:

[0142]

[0143] 35) Calculate the decision value;

[0144] Use the decision value γ for each sample point xi i To express the possibility of becoming a cluster center, the specific calculation is as follows:

[0145] γ i =ρ i *δ i (10)

[0146] Among them, the decision values ​​γ are sorted in ascending order. i Arrange and get γ i The sequence γ' arranged from small to large;

[0147] 36) Calculate the inflection point position and obtain the number of clusters and cluster centers;

[0148] The cluster center is selected on the sample point with a large decision value. The specific calculation formula is as follows:

[0149] kp=max{i||μ i+1 -μ i |≥θ,i=n-DN+1,n-DN+2,...,n-2} (11)

[0150]

[0151]

[0152] μ i =γ' i+1 -γ' i ,i=n-DN+1,n-DN+2,...,n-1 (14)

[0153] Where DN is the nearest integer value after taking the square root of the number of sample points n; γ' i represents the i-th decision value in the sequence γ', γ' i+1 represents the i+1th decision value in the sequence γ'; μ i Represents the decision value γ' i+1 With γ' i The difference between i+1 Represents the decision value γ' i+2 With γ' i+1 The difference between them; kp is the point with the largest ranking value, which satisfies its corresponding μ i Exceeding a certain threshold, with DN as the starting boundary, the threshold θ is μ in the interval [DN,n] i The average value of

[0154] 37) Determine the clustering results;

[0155] According to the cluster centers obtained in step 36), based on the comprehensive similarity measurement matrix, the non-cluster center drone trajectory slice samples are respectively classified into the cluster center with the highest similarity to obtain the clustering result.

[0156] 4) Based on the clustering results in step 3), extract the flight intention features of the UAV trajectory slice samples of each cluster, and set the model input and output; the specific process is as follows:

[0157] 41) Extract flight intention features;

[0158] The flight intention features are represented as the UAV's four horizontal flight behaviors, three longitudinal flight behaviors, and speed changes: the horizontal flight behavior includes the heading angle interval between two adjacent moments, which is divided into four intervals at 90° intervals; the longitudinal flight behavior is the pitch angle interval between two adjacent moments, including descent, climb, and level flight; the speed change consists of three situations: acceleration, constant speed, and deceleration between two adjacent moments;

[0159] 42) Calculate flight intention characteristic parameters;

[0160] The first trajectory point of the trajectory slice sample is marked as the origin (0,0,0), then the position coordinates of the i-th trajectory point in the slice sample are expressed as vector P i , P i =(x i ,y i ,z i ), and then the sequence input features are obtained, which are the horizontal heading angle A of the drone between each two adjacent moments in the trajectory slice sample. hi , the vertical pitch angle A vi and acceleration a i , the specific calculation formula is as follows:

[0161]

[0162]

[0163]

[0164] Where Δt is the sampling interval in seconds; |*| is the modulus of the vector; V i+1 and V i are the velocities of the drone at the i+1th and ith trajectory points, respectively, i =(x i ,y i ,z i ) is the position coordinate of the i-th trajectory point in the trajectory slice sample, P i+1 =(x i+1 ,y i+1 ,z i+1 ) is the position coordinate of the i+1th trajectory point in the trajectory slice sample, H i+1 is the vector P i+1 -P i Projection on the horizontal plane, H i+1 =(x i+1 -x i ,y i+1 -y i ,0);L i+1 P i+1 -Pi The projection in the latitude direction, L i+1 =(0,y i+1 -y i ,0);

[0165] 43) Set input sample;

[0166] The flight intention sequence information is represented as the change process between two adjacent moments, with a time step of T in -1, the flight intention recognition input sample is [b,T inten ], b is the number of batch samples, T inten is the time series intention information of the UAV, and the specific calculation formula is as follows:

[0167]

[0168] 44) Processing output data;

[0169] For the output of flight intention recognition, the behavioral feature values ​​of horizontal flight behavior, longitudinal flight behavior and speed change in step 41) are mapped to a unified interval range, and the output processing of the category label is selected through one-hot encoding. The flight intention of each state quantity at each moment is one-hot encoded as the input of the recognition model, and the category information is retained and converted into a numerical form that can be learned and processed by the algorithm.

[0170] 5) Based on the flight intention features extracted from the trajectory slice samples of each cluster of drones, they are trained separately based on the GRU algorithm to obtain a set of drone flight intention recognition models; the specific process is as follows:

[0171] 51) Determine a training sample set;

[0172] Determine the hybrid trajectory data slice set S, S = {s 1 ,s 2 ,s 3 ,···,s r}, mixed flight intention data slice set D, D = {d 1 ,d 2 ,d 3 ,···,d r}, and a neural network model Q consisting of an input layer, three parallel GRU layers, a parallel fully connected layer, and an output layer, where S and D are divided into training set S in the ratio of 60%, 20%, and 20% respectively. t , D t , validation set S v , D v , test set S c , D c,The data inclusion based on the training sample set is as follows;

[0173] Scenario 1: A non-cooperative drone type can be detected, but there is no real historical trajectory of the corresponding type of drone in the training sample set;

[0174] Scenario 2: The type of non-cooperative drone cannot be detected, but the training sample set contains some real historical trajectories of the corresponding type of drone;

[0175] 52) Calculate the slice comprehensive similarity measure matrix;

[0176] Based on the mixed trajectory data slice set S t , according to the similarity measurement method of step 2), the comprehensive similarity measurement matrix M in step 31) is obtained;

[0177] 53) Calculate the number of clusters and cluster centers;

[0178] According to the clustering method in step 3), the number of clusters m and the cluster center c are obtained. i ,i=1,2,...,m;

[0179] 54) Divide the flight intention slice sample set;

[0180] According to the cluster center c i and the similarity between each trajectory sample slice and the cluster center, and the non-cluster center UAV trajectory slice samples are respectively classified into the cluster center with the highest similarity, and then the mixed flight intention data slice set D t Flight intention slice sample sets divided into different categories

[0181] 55) Obtaining a set of UAV flight intention recognition models;

[0182] Slice the mixed flight intention data into a set D t The sample subset for each category in As the input of model Q, the GRU algorithm is used for training to obtain the sample subsets Flight intention recognition model under training i And save it, and then get the UAV flight intention recognition model set O,

[0183] 6) Perform intention recognition based on the UAV flight intention recognition model set to obtain the UAV flight intention recognition result; the specific process is as follows:

[0184] 61) Determine the corresponding flight intention recognition model of the test sample;

[0185] According to the similarity measurement method in step 2), calculate the test set D cThe test samples and the cluster centers c i The similarity value is combined with the flight intention recognition model set O saved in step 55), and the flight intention recognition model O corresponding to the maximum similarity value is selected. i The flight intention recognition model as a test sample;

[0186] 62) Obtaining intent recognition results;

[0187] Input the test sample into the determined flight intention recognition model O i , the flight intention recognition probabilities of the UAV’s horizontal, vertical and speed changes are obtained respectively, and the maneuvering intention features corresponding to the maximum probability are taken as the results of intention recognition.

[0188] Taking Mini 2 and Mavic 2 as non-cooperative drones as examples, the flight intention recognition accuracy Train with the mixed trajectory slice dataset as the training sample is calculated based on the intention recognition results of the test samples in the prediction duration of 5s, 10s, 15s, and 20s in scenario 1 and scenario 2. a1 With Train a2 , the flight intention recognition accuracy Train using the trajectory slice dataset after similarity measurement and clustering screening as the training sample c1 With Train c2 In general, the flight intention recognition accuracy result using the corresponding type of drone trajectory slice dataset as the training sample is the best. s As a comparison item, the UAV flight intention recognition results are evaluated. As shown in Table 1, the flight intention recognition results based on the GRU model at different time lengths, the method proposed in the present invention has excellent performance in the flight intention recognition method for UAVs that lack sufficient real historical trajectories compared to the flight intention recognition model without sample screening; Table 1 is as follows:

[0189] Table 1

[0190]

[0191]

[0192] The present invention has many specific application paths. The above is only a preferred embodiment of the present invention. It should be pointed out that for ordinary technicians in this technical field, several improvements can be made without departing from the principle of the present invention. These improvements should also be regarded as the protection scope of the present invention.

Claims

1. A non-cooperative UAV flight intention recognition method based on similarity measurement, characterized in that: Here are the steps: 1) Collect historical trajectory data of various types of drones, process the data, and obtain a sample set of drone historical trajectory slices; 2) performing similarity measurement on the trajectory slice sample set of the historical trajectory of the UAV obtained in the above step 1); 3) Based on the similarity measurement result in step 2), a clustering algorithm is used to cluster the sample set of the drone historical trajectory slices to obtain a clustering result; 4) Based on the clustering results in step 3), extract the flight intention features of the UAV trajectory slice samples of each cluster and set the model input and output; 5) Based on the flight intention features extracted from the trajectory slice samples of each cluster of drones, they are trained separately based on the GRU algorithm to obtain a set of drone flight intention recognition models; 6) Performing intention recognition based on the UAV flight intention recognition model set to obtain the UAV flight intention recognition result; The specific process of step 2) is as follows: 21) Select similarity measurement indicators; The trajectory slice samples are similarity measured, the dimension is [T in ,f],T in is the time length of the trajectory slice sample, and f represents the feature number of the trajectory slice sample. Based on the motion characteristics of the UAV, it is characterized by the changes in flight speed and turning angle, and similarity measurement indicators are selected, including speed change size, speed change trend, turning angle cosine similarity, and turning angle change amplitude. 22) Similarity measurement of speed change; The Pearson correlation coefficient is used to represent the similarity of the speed between slices. The Pearson correlation coefficient r between trajectory slice sample i and trajectory slice sample j is ij as follows: In the formula, s ie 、s je They represent the velocity values ​​corresponding to the trajectory slice samples i and j at the e-th time point, n represents the number of slice samples, l is the slice length of the trajectory slice sample, μ i , μ j are the mean speeds of sample i and sample j respectively; 23) Speed ​​change trend similarity measurement; Calculate the product of the speed mean ratio and variance ratio of each two trajectory slice samples, and the speed change trend similarity RVC between trajectory slice sample i and trajectory slice sample j ij as follows: In the formula, Respectively represent the velocity mean of the i-th and j-th trajectory slice samples, σ i , σ j Respectively represent the velocity variance of the i-th and j-th trajectory slice samples; 24) Cosine similarity measurement of turning angle; Calculate the cosine similarity of the turning angles of the two trajectory slice samples in the corresponding time. According to the turning angle vector pair between each two time points of the slice sample, calculate the cosine similarity of the turning angles between the corresponding vector pairs of the trajectory. k , and then get the cosine sequence {cos k }, k = 1, 2, ..., l–1; l is the slice length of the trajectory slice sample; the cosine similarity of the turning angle between trajectory slice sample i and trajectory slice sample j is cos ij-k as follows: In the formula, a ik-p 、a jk-p Respectively represent the pth component of the turning angle vector of the i-th and j-th trajectory slice samples in the corresponding time interval k, q represents the number of elements in the turning angle vector coordinates, and the cosine similarity is solved for every two trajectory slice samples to obtain a sequence containing l-1 cosine values. The average value of the sequence, cos ij Perform similarity representation using cosine similarity; 25) Similarity measurement of turning angle change amplitude; Calculate the position of the maximum difference in turning angle between two trajectory slice samples The ratio of the slice length l of the trajectory slice sample and the initial turning angle difference A of the two trajectory slice samples ij-initial The maximum difference A ij-max The product of the two is used as the similarity R of the turning angle change between the two trajectory slice samples. ij , the similarity R of the turning angle change between trajectory slice sample i and trajectory slice sample j ij as follows: In the formula, represents the distance from the starting point to the point where the maximum difference in turning angle between trajectory slice sample i and trajectory slice sample j occurs, l is the slice length of the trajectory slice sample, and A ij-initial Represents the initial turning angle difference between two trajectory slice samples, A ij-max is the maximum turning angle difference between two trajectory slice samples; 26) Calculate the similarity measure matrix; For the drone trajectory slice sample set obtained in step 1), the similarity measurement values ​​between any two trajectory slice samples are calculated according to steps 22) to 25) to obtain the corresponding four similarity measurement matrices; The specific process of step 3) is as follows: 31) Calculate the comprehensive similarity measure matrix; The four similarity measurement matrices obtained in step 26) are weighted averaged according to the average weight η=0.25 to obtain a comprehensive similarity measurement matrix M. The specific calculation formula is as follows: In the formula, s ij =η×r ij +η×RVC ij +η×cos ij +η×R ij ,s ij represents the similarity value between trajectory slice sample i and trajectory slice sample j in M; 32) Calculate the shared similarity samples between any two trajectory slice samples in the comprehensive similarity measurement matrix; The similarity measurement threshold is used to screen shared similar samples. The specific calculation formula is as follows: SNN(x i ,x j )=τ(x i )∩τ(x j ) (6) In the formula, τ(x i ) represents the set of similar samples of sample i under the specified threshold; τ(x j ) represents the set of similar samples of sample j under the specified threshold; SNN(x i ,x j ) is given by τ(x i ) and τ(x j ) is obtained by taking the intersection, indicating that samples i and j share a similar sample set; Denotes the set SNN(x i ,x j ) The number of elements in s ip represents the similarity between sample i and sample p in the comprehensive similarity measurement matrix; s jp Represents the similarity between sample j and sample p in the comprehensive similarity measurement matrix; Sim(x i ,x j ) represents the SNN of sample i and sample j (x i ,x j ) Similarity value, which is calculated only when sample i and sample j exist in each other's shared similar samples, otherwise it is zero; 33) Calculate the local similarity sample density ρ i , the specific calculation formula is as follows: In the formula, L(x i ) represents the set of sample points that meet the similarity threshold requirement with sample i under the similarity measurement threshold, Sim(x i ,x j ) indicates that sample i and sample j share similar sample SNN (x i ,x j ) similarity; 34) Calculate relative similarity δ i ; Calculate the relative similarity δ i , which represents the maximum similarity between sample i and its nearby high-density points. The specific calculation formula is as follows: 35) Calculate the decision value; Use each sample point x i The decision value of γ i To express the possibility of becoming a cluster center, the specific calculation is as follows: c i =ρ i *d i (10) Among them, the decision values ​​γ are sorted in ascending order. i Arrange and get γ i The sequence γ' arranged from small to large; 36) Calculate the inflection point position and obtain the number of clusters and cluster centers; The cluster center is selected on the sample point with a large decision value. The specific calculation formula is as follows: kp=max{iμ i+1 -m i ≥θ,i=n-DN+1,n-DN+2,...,n-2} (11) m i =c' i+1 -c' i ,i=n-DN+1,n-DN+2,...,n-1 (14) Where DN is the nearest integer value after taking the square root of the number of sample points n; γ' i represents the i-th decision value in the sequence γ', γ' i+1 represents the i+1th decision value in the sequence γ'; μ i Represents the decision value γ' i+1 With γ' i The difference between i+1 Represents the decision value γ' i+2 With γ' i+1 The difference between them; kp is the point with the largest ranking value, which satisfies its corresponding μ i Exceeding a certain threshold, with DN as the starting boundary, the threshold θ is μ in the interval [DN,n] i The average value of 37) Determine the clustering results; According to the cluster centers obtained in step 36), the non-cluster center drone trajectory slice samples are respectively classified into the cluster center with the highest similarity according to the comprehensive similarity measurement matrix to obtain the clustering result; The specific process of step 4) is as follows: 41) Extract flight intention features; The flight intention features are represented as the UAV's four horizontal flight behaviors, three longitudinal flight behaviors, and speed changes: the horizontal flight behavior includes the heading angle interval between two adjacent moments, which is divided into four intervals at 90° intervals; the longitudinal flight behavior is the pitch angle interval between two adjacent moments, including descent, climb, and level flight; the speed change consists of three situations: acceleration, constant speed, and deceleration between two adjacent moments; 42) Calculate flight intention characteristic parameters; The first trajectory point of the trajectory slice sample is marked as the origin (0,0,0), then the position coordinates of the i-th trajectory point in the slice sample are expressed as vector P i , P i =(x i ,y i ,z i ), and then the sequence input features are obtained, which are the horizontal heading angle A of the drone between each two adjacent moments in the trajectory slice sample. hi , the vertical pitch angle A vi and acceleration a i , the specific calculation formula is as follows: Where Δt is the sampling interval in seconds; |*| is the modulus of the vector; V i+1 and V i are the velocities of the drone at the i+1th and ith trajectory points, respectively, i =(x i ,y i ,z i ) is the position coordinate of the i-th trajectory point in the trajectory slice sample, P i+1 =(x i+1 ,y i+1 ,z i+1 ) is the position coordinate of the i+1th trajectory point in the trajectory slice sample, H i+1 is the vector P i+1 -P i Projection on the horizontal plane, H i+1 =(x i+1 -x i ,y i+1 -y i ,0);L i+1 P i+1 -P i The projection in the latitude direction, L i+1 =(0,y i+1 -y i ,0); 43) Set input sample; The flight intention sequence information is represented as the change process between two adjacent moments, with a time step of T in -1, the input sample of flight intention recognition is [b,T inten ], b is the number of batch samples, T inten is the time series intention information of the UAV, and the specific calculation formula is as follows: 44) Processing output data; For the output of flight intention recognition, the behavior characteristic values ​​of horizontal flight behavior, longitudinal flight behavior and speed change in step 41) are mapped to a unified interval range, and the output processing of the category label is selected through one-hot encoding. The flight intention of each state quantity at each moment is respectively one-hot encoded as the input of the recognition model, and the category information is retained and converted into a numerical form that can be learned and processed by the algorithm; The specific process of step 5) is as follows: 51) Determine a training sample set; Determine the hybrid trajectory data slice set S, S = {s1, s2, s3, ···, s r }, mixed flight intention data slice set D, D = {d1, d2, d3, ···, d r }, and a neural network model Q consisting of an input layer, three parallel GRU layers, a parallel fully connected layer, and an output layer, where S and D are divided into training set S in the ratio of 60%, 20%, and 20% respectively. t , D t , validation set S v , D v , test set S c , D c ,The data inclusion based on the training sample set is as follows; Scenario 1: A non-cooperative drone type can be detected, but there is no real historical trajectory of the corresponding type of drone in the training sample set; Scenario 2: The type of non-cooperative drone cannot be detected, but the training sample set contains some real historical trajectories of the corresponding type of drone; 52) Calculate the slice comprehensive similarity measure matrix; Based on the mixed trajectory data slice set S t , according to the similarity measurement method of step 2), the comprehensive similarity measurement matrix M in step 31) is obtained; 53) Calculate the number of clusters and cluster centers; According to the clustering method in step 3), the number of clusters m and the cluster center c are obtained. i ,i=1,2,...,m; 54) Divide the flight intention slice sample set; According to the cluster center c i and the similarity between each trajectory sample slice and the cluster center, and the non-cluster center UAV trajectory slice samples are respectively classified into the cluster center with the highest similarity, and then the mixed flight intention data slice set D t Flight intention slice sample sets divided into different categories 55) Obtaining a set of UAV flight intention recognition models; Slice the mixed flight intention data into a set D t The sample subset for each category in As the input of model Q, the GRU algorithm is used for training to obtain the sample subsets Flight intention recognition model under training i And save it, and then get the UAV flight intention recognition model set O, 2. The non-cooperative UAV flight intention recognition method based on similarity measurement according to claim 1 is characterized in that: In the step 1), the historical trajectory data of various types of drones are collected at a sampling frequency of 10 Hz, and the timestamp, identity ID, longitude, latitude, relative flight altitude, speed and attitude angle information of each trajectory point are recorded.

3. The non-cooperative UAV flight intention recognition method based on similarity measurement according to claim 1 is characterized in that: The data processing in step 1) includes: coordinate conversion, equal-interval sequence processing, trajectory slicing and screening, and dimensional processing.

4. The non-cooperative UAV flight intention recognition method based on similarity measurement according to claim 1 is characterized in that: The specific process of step 6) is as follows; 61) Determine the corresponding flight intention recognition model of the test sample; According to the similarity measurement method in step 2), calculate the test set D c The test samples and the cluster centers c i The similarity value is combined with the flight intention recognition model set O saved in step 55), and the flight intention recognition model O corresponding to the maximum similarity value is selected. i The flight intention recognition model as a test sample; 62) Obtaining intent recognition results; Input the test sample into the determined flight intention recognition model O i , the flight intention recognition probabilities of the UAV’s horizontal, vertical and speed changes are obtained respectively, and the maneuvering intention features corresponding to the maximum probability are taken as the results of intention recognition.

Images