A low-altitude unmanned aerial vehicle signal analysis method and system based on multi-source features

Abstract

The present application relates to the technical field of signal analysis, and particularly relates to a low-altitude unmanned aerial vehicle signal analysis method and system based on multi-source characteristics. The method comprises the following steps: controlling a multi-source sensor array to perform full polarization scanning on a low-altitude target airspace, collecting an original radio frequency signal stream, and generating a multi-channel original signal tensor; performing nonlinear decoupling and channel estimation on the multi-channel original signal tensor and constructing a space-time-frequency three-dimensional feature map to generate a space-time-frequency feature cube; inputting the space-time-frequency feature cube into an intention perception engine to generate a behavior intention candidate set; activating a corresponding interference strategy module according to the behavior intention candidate set and generating a directional countermeasure signal waveform, while performing pulse shaping and time slot alignment to generate a high-hit-rate interference pulse sequence; and sending the high-hit-rate interference pulse sequence into a radio frequency power amplifier, radiating through a directional antenna array, and recording a countermeasure operation log. The present application can greatly improve the discrimination efficiency of unmanned aerial vehicle flight intentions.

Description

Technical Field

[0001] This invention relates to the field of signal analysis technology, specifically to a method and system for low-altitude unmanned aerial vehicle (UAV) signal analysis based on multi-source features. Background Technology

[0002] With the widespread adoption of consumer and industrial-grade small drones, incidents of unauthorized drones entering airports, important industrial parks, and key protected areas are frequent, leading to a continuous increase in demand for low-altitude security and low-altitude control in key areas. The low-altitude open airspace environment is complex, containing a large number of electromagnetic signals from broadcast communications, mobile communications, and various industrial wireless devices, creating a dense and complex electromagnetic clutter environment that poses a significant challenge to drone target detection and signal identification. Traditional low-altitude drone control systems mostly use single-type sensors for signal acquisition, relying solely on single-dimensional radio frequency characteristics for drone identification. This limited detection method has weak anti-interference capabilities and is highly susceptible to interference from various clutter signals in the airspace, resulting in problems such as missed target detection and high false alarm rates.

[0003] Current mainstream methods for detecting low-altitude drones include electronic reconnaissance, millimeter-wave radar, lidar, optoelectronic monitoring, and acoustic signature detection. Each type of sensor has its own advantages and disadvantages: electronic reconnaissance excels at capturing radio frequency signals, radar has ranging and velocity measurement capabilities, and optoelectronic devices can achieve visual target confirmation. Most existing prevention and control solutions deploy and analyze data from various sensors independently, lacking a unified mechanism for collaborative acquisition of multi-source signals and joint feature extraction, thus failing to achieve multi-sensor signal fusion analysis. Single-sensor identification modes are insufficient for extracting weak drone radio frequency signals in complex electromagnetic environments, cannot distinguish between civilian wireless clutter and drone remote control and image transmission signals, and struggle to accurately determine the drone's target identity and flight intent, resulting in an inability to conduct targeted countermeasures. Summary of the Invention

[0004] Therefore, it is necessary for the present invention to provide a low-altitude UAV signal analysis method and system based on multi-source features to solve at least one of the above-mentioned technical problems.

[0005] To achieve the above objectives, a low-altitude UAV signal analysis method based on multi-source features includes the following steps: Step S1: Control the multi-source sensor array to perform full polarization scanning of the low-altitude target airspace, acquire the original radio frequency signal stream containing amplitude, phase and time-frequency characteristics, and generate multi-channel original signal tensors; Step S2: Perform nonlinear decoupling and channel estimation on the multi-channel original signal tensor, extract the channel state information vector corresponding to each channel, and construct a space-time-frequency three-dimensional feature map based on the channel state information vector to generate a space-time-frequency feature cube; Step S3: Input the spatiotemporal feature cube into the intent perception engine, call the preset behavior prototype in the UAV behavior pattern library, extract the UAV identity fingerprint tag and behavior intent vector corresponding to the behavior prototype, and generate a behavior intent candidate set; Step S4: Activate the corresponding interference strategy module according to the candidate set of behavioral intentions, dynamically calculate the optimal transmission power and radiation pattern of the interference signal in combination with the background noise spectrum of the current electromagnetic environment, generate the directional countermeasure signal waveform, and perform pulse shaping and time slot alignment on the directional countermeasure signal waveform to generate a high hit rate interference pulse sequence. Step S5: The high-hit-rate jamming pulse sequence is sent to the radio frequency power amplifier and radiated to the target airspace through the directional antenna array. At the same time, the transmission time, power, azimuth angle and frequency band information are recorded to generate a countermeasure operation log.

[0006] Furthermore, the present invention also provides a low-altitude unmanned aerial vehicle (UAV) signal analysis system based on multi-source features, including a processor, a memory, and a computer program stored in the memory and executable on the processor, for performing the low-altitude UAV signal analysis method based on multi-source features as described above.

[0007] The beneficial effects of this invention are: The low-altitude UAV signal analysis method and system based on multi-source features proposed in this invention have the following advantages over existing technologies: By controlling a multi-source sensor array to complete a full polarization scan of the low-altitude target airspace, the original radio frequency signal stream containing amplitude, phase, and time-frequency characteristics is collected, and a multi-channel original signal tensor is constructed. This step relies on the array to carry out a full polarization collaborative scan, simultaneously collecting multi-dimensional radio frequency basic information. The original features of various electromagnetic signals in the airspace are completely preserved through the multi-channel tensor form, overcoming the information limitations of single-channel acquisition, distinguishing the inherent feature differences between UAV target signals and environmental clutter, and providing complete and homogeneous original signal data support for subsequent channel decoupling and multi-dimensional feature construction. By performing nonlinear decoupling and channel estimation on the multi-channel raw signal tensors, channel state information vectors for each channel are extracted. Based on these vectors, a three-dimensional feature map of space-time-frequency is constructed, generating a space-time-frequency feature cube. This step decouples the multi-channel signals, eliminating feature distortion caused by mutual channel interference. It mines the inherent attributes of the signal from the spatial, temporal, and frequency domains, forming a highly discriminative three-dimensional feature structure. This suppresses feature interference from environmental electromagnetic noise, improves the recognition of weak UAV signals, and provides highly reliable feature input for UAV identity and behavior identification. The space-time-frequency feature cube is then fed into the intent perception engine, which calls preset behavior prototypes from the UAV behavior pattern library to extract corresponding identity fingerprint tags and behavior intent vectors, outputting a candidate set of behavior intents. This step, based on the three-dimensional features, matches the target identity with the behavior pattern, establishes the correlation between signal features and flight behavior, and outputs multiple sets of candidate intent results. This upgrades from simple target detection to behavior cognition, providing a decision-making basis for subsequent dynamic matching of targeted interference strategies. Based on the candidate set of behavioral intentions, the corresponding interference strategy module is activated. The transmission power and radiation pattern of the interference signal are dynamically calculated using the background noise spectrum of the electromagnetic environment, generating a directional countermeasure waveform. This waveform is then processed through pulse shaping and time slot alignment to obtain an interference pulse sequence. This step dynamically adjusts the interference parameters in conjunction with the real-time electromagnetic environment of the airspace, using the directional waveform to constrain the radiation range. Pulse timing optimization further enhances the signal superposition effect, strengthening the suppression of the target UAV's communication link and reducing negative impacts on surrounding legitimate wireless devices. The interference pulse sequence is then fed into an RF power amplifier and radiated into the target airspace via a directional antenna array. Simultaneously, key transmission parameters are recorded to generate a countermeasure operation log. This step completes the closed-loop execution of the directional radiation of the interference signal, synchronously and completely archiving all countermeasure operation parameters to form a long-term storable operation record. This improves the traceability capability of low-altitude security response processes and supports the review, operation and maintenance management, and compliance control of airspace events in key protection areas. Attached Figure Description

[0008] Other features, objects, and advantages of the invention will become more apparent from the following detailed description of non-limiting embodiments with reference to the accompanying drawings: Figure 1 This is a schematic diagram of the steps of the low-altitude UAV signal analysis method based on multi-source features of the present invention; Figure 2 for Figure 1 A detailed flowchart of step S2. Detailed Implementation

[0009] The technical method of the present invention will now be clearly and completely described with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of the present 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.

[0010] To achieve the above objectives, please refer to Figures 1 to 2 This invention provides a method for analyzing low-altitude unmanned aerial vehicle (UAV) signals based on multi-source features. In the embodiments of this invention, please refer to... Figure 1 The diagram shown is a block diagram of the low-altitude UAV signal analysis system based on multi-source features according to the present invention. In this example, the low-altitude UAV signal analysis method based on multi-source features includes the following steps: Step S1: Control the multi-source sensor array to perform full polarization scanning of the low-altitude target airspace, acquire the original radio frequency signal stream containing amplitude, phase and time-frequency characteristics, and generate multi-channel original signal tensors; In this embodiment, when controlling the multi-source sensor array to perform a full-polarization scan of the low-altitude target airspace, a uniform circular array consisting of eight dual-polarized antenna elements is first deployed. The array radius is 15 cm, and the operating frequency covers two ISM bands: 2.4 GHz to 2.5 GHz and 5.7 GHz to 5.8 GHz. Each antenna element simultaneously receives radio frequency signals from both horizontal and vertical polarization modes. The scan control module issues a trigger command, driving the array to perform a mechanical rotation scan at a rate of 10 times per second, with each rotation covering a 360-degree azimuth angle and a 90-degree elevation angle range. At each scan position, the receiving link digitally acquires the radio frequency signal at a sampling rate of 100 MHz, with each acquisition lasting 10 milliseconds, obtaining 2048 sampling points. The acquired data includes in-phase components. and orthogonal components The two signals are merged into a single complex sequence. For 8 antenna elements, 2 polarizations, and 2048 sampling points, each scan generates data with a dimension of 8 × 2 × 2048. After 100 consecutive scans, all data are stacked along the time dimension to generate a fifth-order tensor. Its dimensions are 8×2×2048×100. This tensor contains the original information of amplitude, phase, and time-frequency characteristics. The amplitude is calculated... The square root is obtained, and the phase is calculated. The time-frequency characteristics are obtained by performing a short-time Fourier transform on the complex sequence of each channel. The final generated multi-channel original signal tensor... This serves as the input data for subsequent nonlinear decoupling and channel estimation.

[0011] Step S2: Perform nonlinear decoupling and channel estimation on the multi-channel original signal tensor, extract the channel state information vector corresponding to each channel, and construct a space-time-frequency three-dimensional feature map based on the channel state information vector to generate a space-time-frequency feature cube; In this embodiment, the multi-channel original signal tensor When performing nonlinear decoupling and channel estimation, the fifth-order tensor is first... Remodeled into a fourth-order tensor By merging the polarization and antenna dimensions, a third-order tensor with dimensions of 16 × 2048 × 100 is obtained. A Hilbert transform is performed on the time-domain signal of each channel to construct an analytic signal, separating the real and imaginary components. The instantaneous amplitude of each channel is then calculated based on the analytic signal. and instantaneous phase An amplitude-phase joint histogram is constructed, with the phase axis divided into 360 grids and the amplitude axis into 256 grids. The number of points in each grid is counted to generate an amplitude-phase feature map. A short-time Fourier transform is performed on the amplitude-phase feature map of each frame using a Hamming window with a window length of 256 points and an overlap rate of 50%, yielding the time-frequency energy spectrum. Channel estimation is performed based on the time-frequency energy spectrum. The energy values ​​of each subcarrier are extracted in the frequency domain, and the channel response is calculated using the least squares criterion. ,in The received pilot subcarrier energy, The pilot symbol is transmitted. The channel responses of all subcarriers are combined into a Channel State Information (CSI) vector of length K×2, where K is the total number of subcarriers set to 64, including fading factors and phase rotation factors. The CSI is normalized, the spatial correlation matrix is ​​calculated, and principal component analysis is performed to extract the first three principal component directions. The principal component vector is then multiplied by the time and frequency dimensions to generate a high-dimensional feature tensor. The high-dimensional feature tensor is compressed to 128 dimensions using an autoencoder, and then expanded to N×T×K dimensions through a fully connected layer, where N=16, T=100, and K=64. Finally, it is reshaped into a three-dimensional tensor G with dimensions of 16×100×64, which is the space-time-frequency feature cube.

[0012] Step S3: Input the spatiotemporal feature cube into the intent perception engine, call the preset behavior prototype in the UAV behavior pattern library, extract the UAV identity fingerprint tag and behavior intent vector corresponding to the behavior prototype, and generate a behavior intent candidate set; In this embodiment, when the spatiotemporal feature cube G is input to the intent perception engine, a pre-trained UAV behavior pattern library is first loaded. This library contains three-dimensional template tensors for five behavior prototypes, each template having the same dimensions as G, which is 16×100×64. The intent perception engine performs a three-dimensional cross-correlation operation on G and each behavior prototype template, and the cross-correlation response map is generated. The calculation formula is: Where G is the spacetime-frequency feature cube, The behavior prototype template tensor is set with a stride of 1 and zero-padded boundaries. A 16×100×64 response map is obtained for each behavior prototype, which is then stacked to form a 5×16×100×64 five-channel output tensor. Phase jitter patterns are extracted from each behavior prototype template. Spectrum broadening trend and frame interval stability parameters Construct a behavioral prototype feature dictionary The spatiotemporal feature cube is divided into sliding windows of 16 frames in length and 4 frames in step along the time axis. Local feature vectors are extracted from each window, and cosine similarity is calculated between these windows and feature vectors in the dictionary to obtain the local matching degree score matrix. The dimension is 22×5. A temporal behavior chain is constructed based on the score matrix. The behavior type sequence is obtained by taking the maximum value index of each row, and the state transition probability matrix is ​​calculated. Based on path continuity and stability indices, the current action category is output. and the probability distribution of the next action prediction Query from the preset database Corresponding model Firmware version and encryption protocol types Combined into an identity fingerprint tag .Will and Package and generate a set of behavioral intent candidates. .

[0013] Step S4: Activate the corresponding interference strategy module according to the candidate set of behavioral intentions, dynamically calculate the optimal transmission power and radiation pattern of the interference signal in combination with the background noise spectrum of the current electromagnetic environment, generate the directional countermeasure signal waveform, and perform pulse shaping and time slot alignment on the directional countermeasure signal waveform to generate a high hit rate interference pulse sequence. In this embodiment, based on the behavioral intent candidate set When activating the corresponding interference strategy module, first from Extract The system retrieves the corresponding interference template from the interference template library. The template contains frequency band avoidance strategies. and Power ramp-up curve Modulation mode switching rules. According to... Adjust the interference parameters during hovering behavior. Reduce to 30% and extend Up to twice the speed, the return voyage will be Increase to 150% and shorten Up to 0.5 times. The current electromagnetic environment background noise spectrum is collected, and a frequency sweep measurement is performed on the 20MHz to 6GHz band using a spectrum monitoring module, with a frequency step of 100kHz and a dwell time of 10ms. The signal-to-noise ratio and interference tolerance for each frequency band are calculated. Combined with the target UAV's receiver sensitivity curve, the minimum effective power required for interference in each frequency band is calculated. ,in Set to 6dB. Generate a dynamic power adjustment curve based on the minimum effective power and the preset power ramp-up rule. By combining the target UAV's historical trajectory, the spatial position at the next moment is predicted using a Kalman filter, the azimuth and pitch angles are calculated, and an optimized spatial pointing vector is generated. .Will Slice into 100-millisecond time windows. Discretize the problem at a resolution of 1 degree and solve for the interference efficiency for each time window. A maximized power-pointing pairing scheme is employed. A continuously varying time series of interference parameters is generated through cubic spline interpolation, which is then converted into an executable command sequence for the RF front-end. This command sequence is input into a waveform synthesizer to generate a frequency-domain spectral template. Power weighting and a steering vector are applied, and an inverse Fourier transform is performed to obtain a time-domain baseband waveform. This waveform is then up-converted to the target carrier frequency to generate a directional counter-interference signal waveform. This waveform is then subjected to raised cosine pulse shaping with a roll-off factor of 0.25. The time offset is calculated using preamble correlation detection, and time slot alignment is performed. The waveform is then segmented into a pulse sequence, with each pulse lasting 10% of the receiving window period and spaced 90% apart, generating a high-hit-rate interference pulse sequence. .

[0014] Step S5: The high-hit-rate jamming pulse sequence is sent to the radio frequency power amplifier and radiated to the target airspace through the directional antenna array. At the same time, the transmission time, power, azimuth angle and frequency band information are recorded to generate a countermeasure operation log.

[0015] In this embodiment, a high-hit-rate interference pulse sequence is used. When fed into the RF power amplifier, the pulse waveform of the digital baseband is first converted into an analog signal by a digital-to-analog converter (DAC), with a sampling rate of 1 GHz and a quantization bit depth of 14 bits. The analog signal then passes through an up-converter, converting it from the baseband frequency to the RF carrier frequency. The upconversion uses quadrature modulation to convert the baseband... The two signals are multiplied by the cosine and sine components of the local oscillator respectively, and then added together to obtain the radio frequency signal. The radio frequency signal enters the power amplifier, and the amplifier's gain... Set the gain value according to the optimal power calculated in step S4. ,in This refers to the input power of the power amplifier. To achieve optimal transmission power, the amplified signal is transmitted to the directional antenna array via a feeder. The beam pointing of the antenna array is controlled by the spatial pointing optimization vector optimized in step S4. The phase of each array element is adjusted by a phase shifter to achieve beam pointing towards the predicted target position. During the radiation process, the monitoring module records key parameters for each transmission in real time: transmission timestamp. Acquired via system clock, with an accuracy of 1 microsecond; transmit power The actual output power value is obtained by sampling and detection through a coupler; direction angle Information obtained from encoder readings of the antenna array, with an accuracy of 0.1 degrees; frequency band information. Record the start and end frequencies of the currently used interference band. Write all entries to the countermeasures operation log in chronological order. Each log entry is formatted as follows: The log file is stored in binary format, with each record having a fixed length of 64 bytes, facilitating high-speed writing and subsequent retrieval. This log provides complete and traceable data for evaluating the effectiveness of countermeasures and auditing operations.

[0016] Furthermore, step S2 includes the following steps: Step S21: Perform complex domain decomposition on the multi-channel original signal tensor, separate the real and imaginary components, and generate a complex signal matrix; In this embodiment, when performing complex domain decomposition on the multi-channel original signal tensor, the time-domain discrete signal sequences of the acquired N channels are first... These are combined into a third-order tensor S with dimensions N×T×1, where T represents the number of sampling points in a single channel. An Euler transform is performed on each element of this tensor, mapping it from the real domain to the complex domain, forming a complex signal matrix C. The matrix elements... The core of the decomposition operation lies in constructing an analytic signal using the Hilbert transform, for each channel's real signal. The spectrum of the signal is calculated using a Fast Fourier Transform (FFT), the negative frequency components are set to zero, and then an Inverse Fourier Transform is performed to obtain the complex analytic signal. In this process, the real part 'a' is the original signal. The imaginary part b is The Hilbert transform result. The final generated complex signal matrix C has a dimension of N×T, where each row represents a complex analytic signal sequence of one channel, and the matrix elements... The magnitude reflects the instantaneous envelope of the signal, while the amplitude contains the phase information of the signal, providing a complete complex domain foundation for subsequent amplitude and phase analysis.

[0017] Step S22: Calculate the instantaneous amplitude and instantaneous phase of each channel based on the complex signal matrix, construct the amplitude-phase joint distribution map, and generate the amplitude-phase feature map; In this embodiment, when calculating the instantaneous amplitude and instantaneous phase of each channel based on the complex signal matrix C, for each element in the matrix... Directly calculate the modulus As the instantaneous amplitude, its argument is also calculated. As the instantaneous phase. Where, Im(c n (t) represents the element The imaginary part, Re(c n (t) represents the element The real part of the amplitude-phase joint distribution map is as follows: The instantaneous amplitudes of all N channels at the same time t are... and instantaneous phase Consider it as a point on a two-dimensional plane, with the x-coordinate of the point representing the phase value. The vertical axis represents the amplitude value. Traversing all T sampling times, all N×T points are plotted in a single two-dimensional scatter plot. To generate an amplitude-phase feature map, this scatter plot is further divided into M×P grid cells, where M is the number of grid cells in the phase direction, set to 360, with each cell corresponding to 1 degree; P is the number of grid cells in the amplitude direction, set to 256 based on the signal's dynamic range. The number of points falling into each grid cell is counted, and the point density is represented by color intensity, thus generating a grayscale image. The pixel values ​​of this image represent the count results within that grid cell. This amplitude-phase feature map visually presents the statistical distribution characteristics of the signal under different amplitude and phase combinations.

[0018] Step S23: Perform time-frequency transformation on the amplitude-phase feature map, extract the time-varying spectral energy distribution, and generate a time-frequency energy spectrum; In this embodiment, when performing time-frequency transformation on the amplitude-phase feature map, the grayscale image is first considered as a two-dimensional matrix F with size M×P. Since the image itself does not contain time-dimensional information, a time variable needs to be reintroduced. In practice, steps S21 and S22 are repeated for each signal segment within a time window to generate a series of amplitude-phase feature map sequences that vary over time. ,in This represents the index of the time window. Then, for each frame feature map in the sequence... A short-time Fourier transform is applied, with a Hamming window as the kernel function. The window length is set to 256 sampling points, and the overlap rate is 50%. For each frame of feature map, a one-dimensional Fourier transform is performed along its row direction (i.e., the phase axis) to obtain the spectral coefficient matrix. These spectral coefficients are stacked in time order to form a three-dimensional tensor E, whose three dimensions are time, time, and phase. Frequency f and phase shift k. Elements in tensor E Indicates time window Within the range, the energy density at a frequency component f and a phase shift k is calculated. Finally, the time-frequency energy spectrum is obtained by summing over all phase shifts k. The spectrum clearly shows the distribution of signal energy in both time and frequency dimensions.

[0019] Step S24: Perform channel estimation based on the time-frequency energy spectrum, calculate the channel response vector of each subcarrier, and generate the channel state information vector; In this embodiment, based on time-frequency energy spectrum When performing channel estimation, it is assumed that the system uses orthogonal frequency division multiplexing (OFDM) technology and the total number of subcarriers is K. For each time window... K energy values ​​corresponding to K subcarrier frequency positions are extracted from the time-frequency energy spectrum to form a vector. The core of channel estimation is to solve for the channel response using known pilot symbols. At the transmitting end, pilot symbols with known amplitude and phase are pre-inserted on specific subcarriers. The receiving end from The energy values ​​on these pilot subcarriers are extracted, and combined with the known values ​​of the pilot symbols, the channel response of each subcarrier is calculated using the least squares criterion. The calculation formula is: ,in It is the received pilot subcarrier energy. This refers to the transmitted pilot symbols. For non-pilot subcarriers, their channel response is obtained through interpolation. The interpolation method used is linear interpolation, which estimates the response using the response values ​​of adjacent pilot subcarriers. The resulting channel response vector is... This is the channel state information vector, which describes the current time window. The channel attenuation and phase rotation experienced by each subcarrier are shown below.

[0020] Step S25: Map the channel state information vector to a space-time-frequency three-dimensional grid to construct a space-time-frequency feature cube, where the spatial dimension corresponds to the antenna array index, the time dimension corresponds to the sampling time, and the frequency dimension corresponds to the subcarrier index.

[0021] In this embodiment, when mapping the channel state information vector to a space-time-frequency three-dimensional grid, the three axes of the three-dimensional coordinate system are first defined: the spatial axis corresponds to the antenna array index, ranging from 1 to N; the time axis corresponds to the sampling time index, ranging from 1 to T; and the frequency axis corresponds to the subcarrier index, ranging from 1 to K. For each time window... A channel state information vector obtained from step S24 Only corresponds to a single moment To fill the entire 3D mesh, steps S21 to S24 need to be repeated for all T time windows to obtain a set of channel state information vectors {H1, H2, ..., H...}. T}. Then, a three-dimensional empty array G of size N×T×K is constructed. For each antenna channel n, the corresponding time series is... Assign to grid ,in Indicates the first The channel response vector of the nth antenna within a time window. After this operation, any point in grid G... The value represents the channel response of the nth antenna on the kth subcarrier at time t. This space-time-frequency feature cube G completely characterizes the propagation characteristics of low-altitude UAV signals in a multi-antenna system as time and frequency change, providing structured multi-dimensional feature inputs for subsequent target detection and recognition.

[0022] Furthermore, step S24 includes the following steps: Step S241: Extract the energy peak value corresponding to each subcarrier and its time window position from the time-frequency energy spectrum to generate a set of energy peak position pairs; In this embodiment, from the time-frequency energy spectrum When extracting the energy peak corresponding to each subcarrier, first traverse all K subcarrier indices k. For each fixed subcarrier k, along the time axis... Directional scan energy value sequence Define a sliding window of length L=5, with its center point at... Calculate the mean of all energy values ​​within the window. and standard deviation When the center point Energy value at the location Greater than When the value is zero, mark that point as a candidate peak. Perform non-maximum suppression on all candidate peaks, retaining local maxima points, i.e., if... If the energy value at a given point is greater than the energy values ​​at the two adjacent time points before and after it, then that point is confirmed as an energy peak. Record the energy value of this peak. and their corresponding time window positions , where m represents the m-th peak detected on the subcarrier. For each subcarrier k, the positions of all detected peaks are... Form a set Finally, the peak position pairs of all subcarriers are merged into a global set. Each element in this set uniquely corresponds to a subcarrier index and a specific time window position, and completely records the spatiotemporal distribution information of energy bursts on each subcarrier.

[0023] Step S242: Construct an energy propagation path graph between subcarriers based on the set of energy peak positions. The nodes in the path graph are subcarrier indices, and the weight of the edges is the time difference of the energy peak of adjacent subcarriers, generating an energy propagation topology graph. In this embodiment, when constructing the inter-carrier energy propagation path graph based on the set Q of energy peak locations, the K subcarrier indices are first defined as nodes of the graph. For any two adjacent subcarriers and Extract the peak position sequences of each element in the time dimension from set Q. and The Dynamic Time Warping (VTW) algorithm is used to calculate the optimal alignment path between the two sequences. The core of this algorithm is to recursively find the alignment with the minimum cumulative distance, defined as follows: On the alignment path, calculate the average time difference of all matching point pairs. The average time difference is the connection node. and weight of the edge If there are no matching peak pairs between a pair of adjacent subcarriers, the edge weight is set to infinity, indicating that there is no direct energy propagation relationship. This constructs a weighted directed graph. Its adjacency matrix elements All other elements are zero. This energy propagation topology diagram visually depicts the time delay characteristics of energy transfer between different subcarriers.

[0024] Step S243: Perform graph neural network inference on the energy propagation topology graph, extract the neighbor aggregation features of each node, and generate the topology embedding vector; In this embodiment, the energy propagation topology is described. When performing graph neural network inference, a two-layer graph convolutional network is first constructed. The input layer is the node feature matrix X, with dimensions K×D, where D is set to 64, and each node... The initial characteristics are determined by the normalized frequency index of its subcarrier. Encoded into f kFor the k-th subcarrier frequency, The sampling rate is used. The definition of graph convolution is... ,in Let I be the normalized adjacency matrix, and I be the identity matrix. for The degree matrix, , and The weight matrices are trainable, with dimensions of 64×128 and 128×32, respectively. The first layer of graph convolution outputs a feature matrix of dimension K×128, which, after passing through the ReLU activation function, is input into the second layer of graph convolution, outputting a feature matrix of dimension K×32. .right Global average pooling is performed along the node dimension to obtain a one-dimensional vector of length 32. This vector is then mapped to an output space of length 16 through a fully connected layer. This output vector is the topology embedding vector Z, with a dimension of 1×16. This embedding vector aggregates the neighbor information of all subcarrier nodes in the graph structure, including the energy propagation timing relationships and topological connectivity between nodes, compressing the high-dimensional graph structure information into a compact low-dimensional feature representation.

[0025] Step S244: Channel splicing of the topology embedding vector with the original time-frequency energy spectrum to generate an enhanced energy spectrum that fuses topology information; In this embodiment, the topological embedding vector Z is compared with the original time-frequency energy spectrum. When performing channel stitching, the dimensions of both need to be unified first. The original time-frequency energy spectrum P has a size of T×K, where T is the total number of time windows and K is the total number of subcarriers. The topology embedding vector Z has a length of 16 and is a one-dimensional vector. To perform channel stitching, Z is copied and expanded into a tensor with the same spatial dimensions as P. Specifically, Z is repeated T times along the time axis and then K times along the frequency axis, generating an embedding feature map M with dimensions T×K×16. At the same time, the original time-frequency energy spectrum P is expanded from dimension T×K to a single-channel feature map P' with dimensions T×K×1. Then, P' is stitched with M in the channel dimension to obtain a multi-channel feature map with dimensions T×K×17. This stitching operation integrates graph structure information into the original time-frequency energy spectrum as an additional channel, allowing subsequent processing to utilize both numerical energy distribution information and topological correlation information between subcarriers. The resulting enhanced energy spectrum... Each spatial location Each has 17 feature channels. The first channel retains the original energy value, while the remaining 16 channels carry global topological information from the graph neural network.

[0026] Step S245: Perform complex exponential fitting on the enhanced energy spectrum for each subcarrier, extract the channel fading factor and phase rotation factor corresponding to each subcarrier, and generate the channel state information vector.

[0027] In this embodiment, the enhanced energy spectrum is... When performing a subcarrier-by-subcarrier complex exponential fitting, for each subcarrier k, from Extract the 17-dimensional feature vector sequence of the subcarrier across all time windows. ,in The goal of fitting to T is to find a complex exponential function. To make it consistent with the observation sequence The square root of the first channel (i.e., the original energy value) is closest in the least squares sense. In practice, the complex exponential function to be fitted is rewritten in linear form. After taking the natural logarithm of the observed sequence, a hypothesis about the independent variable is constructed. Solve the system of linear equations using the least squares method to find the parameters. and The real and imaginary parts of the integer. The solution obtained... This is the channel fading factor, which reflects the degree of amplitude attenuation of the signal on that subcarrier; This is the frequency offset used to calculate the phase rotation factor. ,in The sampling interval is given. Finally, for each subcarrier k, the channel state information vector... By decay factor and phase rotation factor Composed of two scalars. The sum of all K subcarriers. When concatenated, they form a one-dimensional vector of length 2K, which is the final Channel State Information (CSI) vector. This vector fully characterizes the propagation characteristics of the signal in a frequency-selective fading channel.

[0028] Furthermore, step S243 includes the following steps: Each node in the energy propagation topology is initialized with the average energy value of that subcarrier, and an initial feature vector of the node is generated. In this embodiment, when initializing each node in the energy propagation topology to the average energy value of that subcarrier, the first step is to start from the time-frequency energy spectrum. Extract the energy time series corresponding to each subcarrier k For each subcarrier k, calculate its arithmetic mean energy over all T time windows. The average value This reflects the average energy level of the subcarrier throughout the entire observation period. This scalar value... As a node The initial features are then constructed. To build an initial feature vector with uniform dimensions, the scalar features of each node are expanded into a fixed-length vector. The dimension D of the feature vector is set to 64, specifically by creating a 64-dimensional vector of all zeros. Assign its first element a value The remaining 63 elements are kept to zero. This sparse coding method preserves the average energy information of each subcarrier, while reserving sufficient dimensional space for subsequent message passing and feature fusion. Therefore, each node... A 64-dimensional initial feature vector was obtained. The initial feature vectors of all nodes constitute a node initial feature matrix of dimension K×64. , where K is the total number of subcarriers, and this matrix serves as the starting point for the inference process of the graph neural network.

[0029] Perform message passing operations on each node, aggregate the feature vectors of its first-order neighbor nodes and sum them with weights to generate neighbor aggregate features; In this embodiment, when performing message passing operations on each node, each node is first determined according to the adjacency matrix A of the energy propagation topology graph. First-order neighbor node set For nodes Its first-order neighbors are defined as all neighbors that satisfy... or nodes The message passing process consists of three steps. The first step is to calculate the value of each edge. message vectors on ,in It is a learnable weight matrix with dimensions 64×64. It is a node The feature vector from the previous iteration. The second step is to process the nodes... The received messages from neighboring nodes are summed using a weighted average, with the weights determined by the edge weights. Decide, That is, the reciprocal of the calculated peak energy time difference between adjacent subcarriers, i.e. ,in To prevent division by zero. The aggregation formula is: The third step is to aggregate the results. The final neighbor aggregation feature vector is obtained through a linear transformation layer. ,in Its dimensions are 64×64. This is a 64-dimensional bias vector. The neighbor aggregation feature... Integrating nodes Information on all first-order neighbors in the surrounding topology, with the inverse of the energy propagation time difference as the importance weight.

[0030] The initial feature vector of a node is fused with the aggregated features of its neighbors using a gating system. The ratio between the two is dynamically adjusted through a learnable gating unit to generate a gating update feature vector. In this embodiment, the initial feature vector of the node is... Features of aggregation with neighbors When performing gated fusion, a learnable gate unit is first designed. The input of this gate unit is the initial features of the nodes. and neighbor aggregation features spliced ​​vector Its dimension is 128. The gating unit consists of a single-layer fully connected network and a sigmoid activation function, outputting a 64-dimensional gating vector. ,in Its dimensions are 64×128. A 64-dimensional bias vector. Gated vector. Each element takes a value between 0 and 1, controlling the mixing ratio of the initial feature and neighboring features during fusion. The formula for gating fusion is as follows: ,in This indicates element-wise multiplication. The formula means that when the gate value is close to 1, the fusion result tends to retain the node's own initial features; when the gate value is close to 0, the fusion result adopts neighbor aggregation features more. The resulting gated update feature vector... It dynamically balances local node information with global topology information, enabling the model to adaptively adjust the importance of information sources according to the characteristics of different subcarriers.

[0031] A nonlinear activation function is applied to the gated update feature vector, and residual connections are introduced to preserve the original information, thereby generating a deep graph feature vector; In this embodiment, the gating update feature vector is performed. When applying a nonlinear activation function and introducing residual connections, firstly... The layer is passed through a fully connected layer with the LeakyReLU activation function, where the negative slope parameter α of LeakyReLU is set to 0.01. The weight matrix of this fully connected layer is shown below. The dimension is 64×64, and the bias vector is... With a dimension of 64, intermediate features were calculated. Subsequently, residual connections are introduced to gate the updated feature vectors. With intermediate features Adding them together yields the feature vector of the deep image. The introduction of residual connections ensures that the original gated fusion information can be directly passed to subsequent layers, avoiding the gradient vanishing problem. To further enhance the expressive power of features, [further details are needed]. Apply another batch normalization layer; the batch normalization operation applies to... Calculate the mean within a mini-batch for each feature dimension. and variance Then, normalization is performed: ,in Finally, multiply by the learnable scaling factor γ and add the offset term β, where both γ and β have dimensions of 64. The result after batch normalization is... This is the final deep graph feature vector, which contains both high-order semantic information after nonlinear transformation and retains the gated fusion features of the underlying layer through residual connections.

[0032] Stack the deep graph feature vectors of all nodes in subcarrier index order to generate a topology embedding vector.

[0033] In this embodiment, when stacking the deep graph feature vectors of all nodes in subcarrier index order, the processed node is first obtained. Deep graph feature vectors The vector has a dimension of 64. K subcarriers correspond to K nodes, and each node corresponds to a subcarrier index k, ranging from 1 to K. The deep graph feature vector of each node is extracted sequentially according to the ascending order of subcarrier index k. Then, these vectors are concatenated in the row direction. Specifically, an empty two-dimensional array Z is created, with K rows and 64 columns. For k ranging from 1 to K, the vectors are concatenated... Assign the value to the k-th row of Z, i.e. After assigning values ​​to all K rows, a two-dimensional matrix Z with dimensions K×64 is obtained. Each row of this matrix corresponds to a deep map feature of a subcarrier, and each column corresponds to a feature dimension. To generate the final topology embedding vector, matrix Z is flattened into a one-dimensional vector. The flattening operation is performed in row-major order, that is, the first 64 elements are taken first, followed by the second 64 elements, and so on, up to the Kth row 64 elements. The final one-dimensional vector has a length L = K×64, and this vector is the topology embedding vector. This embedding vector fully encodes the deep feature information of all subcarrier nodes after inference through a multi-layer graph neural network, compressing the structural characteristics of the entire energy propagation topology into a compact one-dimensional representation.

[0034] Furthermore, step S25 includes the following steps: Step S251: Normalize the channel state information vector so that its mean is zero and its standard deviation is one, and generate a standardized channel vector; In this embodiment, when normalizing the Channel State Information (CSI) vector, the final CSI vector is first obtained. This vector has a length of 2K, where the first K elements are the fading factors of each subcarrier. The last K elements are the phase rotation factors of each subcarrier. Calculate the arithmetic mean of the entire CSI vector. At the same time, calculate its standard deviation. The normalization operation performs a linear transformation on each element of a vector: This transformation makes the mean of the output vector CSI' zero and the standard deviation one. In actual calculations, to avoid division by zero due to a zero standard deviation, a very small constant is added to the denominator. That is, the actual operation is Normalization eliminates the inconsistency in numerical scales caused by power differences or different channel conditions among different subcarriers, enabling subsequent spatial correlation analysis and principal component extraction to be performed within a uniform numerical range. The standardized channel vector CSI' retains the relative change pattern of the original channel state information but removes the overall offset and scale effects, providing a numerically stable input for subsequent feature extraction.

[0035] Step S252: Calculate the spatial correlation matrix based on the standardized channel vector, extract the principal component directions in the spatial dimension, and generate the spatial principal component vector; In this embodiment, when calculating the spatial correlation matrix based on the normalized channel vector CSI', CSI' first needs to be reconstructed into a two-dimensional matrix. Since CSI' contains the fading factors and phase rotation factors of K subcarriers, it is split into two subvectors of length K: A' = [A'1, A'2, ..., A'2]. K ] and θ'=[θ'‌1,θ'‌2,...,,θ'‌ K Using these two sub-vectors as two rows, construct a matrix R of dimension 2×K, where the first row is A' and the second row is... The spatial correlation matrix is ​​calculated as the transpose of matrix R. The covariance matrix, i.e. ,in The dimension of R is K×2, and the dimension of R is 2×K. Their product yields a covariance matrix C of dimension K×K. The elements of matrix C... Let C represent the joint covariance of the i-th subcarrier and the j-th subcarrier in terms of fading factor and phase rotation factor. Perform eigenvalue decomposition on the covariance matrix C and solve the characteristic equation. We obtain K eigenvalues ​​λ1≥λ2≥...≥λ k and its corresponding eigenvectors Select the eigenvectors corresponding to the P largest eigenvalues, and set P to 3, to form the projection matrix. Its dimension is K×3. The spatial principal component vector is calculated as follows: ,in This is the transpose of the normalized channel vector, with a dimension of 2K×1. Because... The dimension is K×3, so CSI' needs to be reshaped into a K×2 matrix before it can be combined with... Matrix multiplication yields a one-dimensional vector of length 3×2=6, which is the spatial principal component vector SP. This vector captures the main variation direction of channel state information in the spatial dimension.

[0036] Step S253: Perform tensor product operation between the spatial principal component vector and the time and frequency dimensions to generate a high-dimensional feature tensor; In this embodiment, when performing a tensor product operation on the spatial principal component vector SP with the time and frequency dimensions, the scale of the three dimensions is first defined. The spatial principal component vector SP has a length of 6, and is considered as a first-order tensor S with a dimension of 6. The time dimension corresponds to the total number of sampling times T, and the frequency dimension corresponds to the total number of subcarriers K. A time basis vector t=[1,2,...,T] of length T and a frequency basis vector f=[1,2,...,K] of length K are defined. The core of the tensor product operation is to perform an outer product expansion on these three first-order tensors to generate a third-order tensor H. The specific operation is as follows: for any spatial index... Time Index and frequency index The value of tensor H at this position is determined by the product of the three corresponding elements, i.e. This operation is equivalent to combining each element of the spatial principal component vector SP with the product of the time index and the frequency index, forming a third-order tensor H of dimension 6×T×K. For example, when hour, This tensor product operation integrates the originally independent spatial, temporal, and frequency information into a unified high-dimensional space, allowing each spatial principal component to expand along the time and frequency axes, thereby generating a high-dimensional feature tensor rich in structural information.

[0037] Step S254: Dimensionality reduction and compression are performed on the high-dimensional feature tensor, and low-dimensional hidden layer representation is extracted using an autoencoder to generate compact spatiotemporal frequency feature representation; In this embodiment, a four-layer autoencoder network is used to perform dimensionality reduction and compression on the high-dimensional feature tensor H. The encoder part of the autoencoder consists of two fully connected layers. The first layer maps the flattened dimension of the input tensor (6×T×K) to an intermediate dimension (512), and the second layer maps from 512 to the target low-dimensional hidden layer dimension (128). The decoder part is symmetrically configured, with the first layer mapping from 128 back to 512, and the second layer mapping from 512 back to the original dimension (6×T×K). The training objective of the autoencoder is to make the reconstructed output as close as possible to the original input, and the loss function is defined as the mean squared error. ,in To reconstruct the tensor, N is the total number of tensor elements. In the actual inference phase, only the encoder part is used. First, the high-dimensional feature tensor H is flattened into a one-dimensional vector. Its length is 6×T×K. Then... The input encoder undergoes the first linear transformation. ,in The dimension is 512×(6TK), It is a 512-dimensional bias vector. Then it undergoes a second linear transformation. ,in Its dimensions are 128×512. Given a 128-dimensional bias vector, a low-dimensional hidden layer representation z with a length of 128 is obtained. This compact hidden layer representation z significantly reduces the data dimensionality, eliminates redundancy and noise, and provides an efficient and expressive feature vector for subsequent generation of the spatiotemporal feature cube.

[0038] Step S255: Reshape the compact spatiotemporal feature representation into a spatiotemporal three-dimensional mesh structure to generate a spatiotemporal feature cube.

[0039] In this embodiment, when reshaping the compact spatiotemporal feature representation z into a spatiotemporal three-dimensional mesh structure, the dimensions of the target three-dimensional mesh are first defined. The spatial dimension of the spatiotemporal feature cube corresponds to the antenna array index N, the time dimension corresponds to the sampling time T, and the frequency dimension corresponds to the subcarrier index K. The length of the compact feature representation z is 128, which needs to be mapped to a three-dimensional space of dimension N×T×K. Since the product of N, T, and K is usually much greater than 128, direct reshaping is not feasible. Therefore, the dimension of z is first expanded to N×T×K through a fully connected extension layer. The weight matrix of this extension layer... The dimension is (NTK)×128, and the bias vector is... The extended vector is calculated with dimension NTK. Its length is NTK. Then, the one-dimensional vector e is reshaped into a three-dimensional tensor G in row-major order, with the reshaping rule being: for any antenna index... Time Index and frequency index The value of tensor G at this position is This three-dimensional tensor G is the final generated spatiotemporal feature cube. Each element of this cube... All of them carry high-level feature information after normalization, principal component extraction, tensor fusion and autoencoder compression, which fully encodes the space-time-frequency joint characteristics of low-altitude UAV signals in multi-antenna systems, providing structured feature inputs for subsequent classification or detection tasks.

[0040] Furthermore, step S3 includes the following steps: Step S31: Input the spatiotemporal frequency feature cube into the intent perception engine and call the preset behavior prototypes in the UAV behavior pattern library, including hovering, straight flight, circling, return to home and loss of contact behavior; In this embodiment, when the spatiotemporal feature cube G is input into the intent perception engine, a pre-trained UAV behavior pattern library is first loaded. This library contains five preset behavior prototypes, namely hovering behavior prototypes. Prototype of straight-line flight behavior , hovering behavior prototype Based on the return journey and the prototype of the missing person behavior Each behavior prototype is stored as a 3D template tensor, with dimensions consistent with the spatiotemporal feature cube G, i.e., N×T×K. The core operation of the intent-aware engine is 3D convolutional matching, which performs matching on each behavior prototype template tensor. Calculate the three-dimensional cross-correlation response plot between it and the input cube G. The formula for calculating the three-dimensional cross-correlation is as follows: ,in This represents the spatial displacement offset. The step size is set to 1, and the boundary padding method is zero-padding to ensure that the output response map has the same size as the input cube. For each behavior prototype, a response map with dimensions N×T×K is obtained. The location of the peak in the response graph indicates the region where the behavioral pattern is most likely to occur in the spatiotemporal space, and the amplitude of the peak reflects the confidence of the match. The intent-aware engine simultaneously calculates the response graphs of five behavioral prototypes and superimposes all the response graphs into a five-channel output tensor R with dimensions of 5×N×T×K, providing basic matching response information for subsequent behavioral feature signature extraction.

[0041] Step S32: Extract the corresponding signal feature signature for each behavior prototype, including phase jitter mode, spectral broadening trend and frame interval stability parameters, and generate a behavior prototype feature dictionary; In this embodiment, when extracting the corresponding signal feature signature for each behavioral prototype, three feature parameters are calculated for each of the five behavioral prototypes. The first parameter is the phase jitter mode, derived from the behavioral prototype template tensor. Extract along the time axis Given a phase change sequence, calculate the sum of the absolute values ​​of the differences in the sequence. ,in For time windows The first term represents the average phase value of all antennas and subcarriers within the array. The second term represents the spectral broadening trend, expressed as a function of the behavioral prototype template tensor along the frequency axis. Perform a Fourier transform to calculate the 3dB bandwidth of the spectrum. This refers to the width of the frequency range where the power spectral density drops to half of its peak value. The third term is the frame interval stability parameter, which calculates the rate of energy change between adjacent frames along the time axis using the behavior prototype template tensor. ,in For time windows The total energy within the frame interval stability parameter is defined as all Standard deviation The reciprocal of, that is , These three parameters are combined into a feature vector of length 3. The above extraction operation was performed on each of the five behavioral prototypes to obtain five feature vectors. , , , , These feature vectors are used to construct a behavior prototype feature dictionary, with the behavior type name as the key and the feature vector as the value. The dictionary contains five key-value pairs, each corresponding to a signal feature signature of a behavioral prototype.

[0042] Step S33: Perform sliding window matching between the feature points in the space-time-frequency feature cube and the behavioral prototype feature dictionary, and calculate the local matching degree score; In this embodiment, the feature points in the spatiotemporal feature cube G are compared with the behavioral prototype feature dictionary. When performing sliding window matching, the spatiotemporal feature cube G is first divided into several local windows along the time axis. The length W of each window is set to 16 time frames, and the sliding step S is set to 4 time frames. For each local window, the local feature cube inside it is extracted. The dimensions are N×W×K. Then, the same feature extraction operation as in step S32 is performed on this local feature cube to calculate its phase jitter mode. Spectrum broadening trend and frame interval stability parameters To obtain local feature vectors The matching score is calculated based on cosine similarity for a dictionary. Each behavioral prototype in ,calculate and Cosine similarity between The similarity score ranges from -1 to 1, with a value closer to 1 indicating a stronger match between the local features and the feature signature of the behavioral prototype. The matching scores of all five behavioral prototypes are combined into a score vector of length 5. After the sliding window traverses the entire timeline, a total of M score vectors are generated. These score vectors constitute a local matching degree score matrix of dimension M×5. Each row of the matrix corresponds to the matching result of a time window.

[0043] Step S34: Construct a temporal behavior chain based on the local matching degree score, identify the evolution path of behavior patterns in consecutive frames, and generate a behavior evolution trajectory map; In this embodiment, based on the local matching degree score matrix When constructing the temporal behavior chain, the maximum value index operation is first performed on each row of the matrix to find the behavior type with the highest score in each time window. This yields a sequence of behavior types of length M. The core of constructing a temporal behavior chain is identifying the evolution path of behavior patterns in consecutive frames. A state transition matrix is ​​defined. Its dimensions are 5×5. This represents the cumulative probability of transitioning from behavior type i to behavior type j. Initialization It is a zero matrix. Traversal behavior is a sequence of types. For two adjacent time windows t and t+1, let , ,Will Increase by 1. After the traversal is complete, Dividing each row by the sum of that row yields the normalized state transition probability matrix. The behavioral evolution trajectory diagram is constructed with time on the horizontal axis and behavioral type on the vertical axis; a point is plotted at each time window t. The color depth of a point is determined by its corresponding matching score. Decision. Connect all points in chronological order to form a broken line path; this path is the behavioral evolution trajectory diagram. The trajectory diagram clearly shows the changes in the drone's behavior over time, such as the transition from hovering to straight flight and then to circling. Each turning point in the path corresponds to the moment of switching behavior patterns.

[0044] Step S35: Based on the continuity and stability of the path in the behavior evolution trajectory map, output the final UAV identity fingerprint tag and behavior intention vector. The identity fingerprint tag includes the model, firmware version and encryption protocol type, and the behavior intention vector includes the current action category and the prediction probability of the next action, thereby generating a candidate set of behavior intentions.

[0045] In this embodiment, based on the behavioral evolution trajectory diagram When outputting the final result for the continuity and stability of the path, the continuity index of the path is calculated first. Defined as the proportion of adjacent time windows in a path that exhibit the same behavior type, i.e. Where I is the indicator function. Stability index Defined as the frequency of occurrence of the behavior type that accounts for the highest proportion in the entire sequence, i.e. .like Greater than the threshold 0.7 and If the value is greater than the threshold of 0.5, the behavior pattern is considered stable, and the behavior type with the highest proportion is directly used as the current action category. If the stability condition is not met, a weighted voting mechanism is adopted, using the matching score of each time window as the weight to calculate the weighted total score for each behavior type. The behavior type with the highest total score is taken as... The probability of predicting the next action is calculated based on the state transition probability matrix. ,Pick corresponding row vector The predicted probability distribution of the next action The generation of identity fingerprint tags relies on the association mapping between behavior types and a pre-set database. For identified fingerprints... Search for its corresponding drone model Firmware version and encryption protocol types The three combine to form a triplet. The final output is a behavioral intent vector. Based on the current action category and the probability distribution of the next action prediction It is pieced together, and its length is 1 + 5 = 6. and Package and generate a set of behavioral intent candidates. This set provides a complete description of the drone's identity attributes and behavioral intentions.

[0046] Furthermore, step S4 includes the following steps: Step S41: Activate the corresponding interference strategy module according to the behavioral intent vector. The module is pre-set with interference templates for different types of UAVs, including frequency band avoidance strategies, power ramp-up curves and modulation mode switching rules. In this embodiment, based on the behavioral intent vector When activating the corresponding interference strategy module, first from Extract the current action category This category is one of five behavior types. The interference strategy module has a pre-built interference template library, which is based on drone models. For indexing, the corresponding interference template is stored. Each interference template contains three core parameters. The first parameter is the frequency band avoidance strategy, defined as a tuple. ,in List of start and end frequencies for the protected frequency band. This is a list of the start and end frequencies of the interference band. The second item is the power ramp-up curve, defined as a piecewise linear function. ,in The starting power, For maximum power, The ramp-up time constant is denoted by t, and the interference duration is denoted by t. The third term is the modulation mode switching rule, defined as a state machine containing four modulation states: continuous wave modulation (CW), noise amplitude modulation (AM), noise frequency modulation (FM), and pulse modulation (PM). The state transition condition is a threshold value for the interference duration; for example, CW modulation occurs when t < 1 second, AM modulation occurs when 1 second ≤ t < 3 seconds, and so on. In Retrieve the corresponding template from the interference template library. Use this template with In The signal is sent to the activation controller, which then outputs an activation signal. This signal initiates the subsequent interference parameter adjustment process.

[0047] Step S42: Adjust the interference parameters based on the behavior intent vector, including reducing the interference power to avoid accidental damage when the behavior is identified as hovering, and increasing the interference to block the return path when the behavior is identified as returning. In this embodiment, based on behavioral intent vectors When adjusting interference parameters, first read the current action category. .like Identified as hovering behavior, the interference power adjustment module will adjust the power ramp-up curve. Reduced to 30% of the original value, that is At the same time, the climb time Extend to twice the original value, that is The purpose of this adjustment is to reduce the strength of interference signals and avoid accidental injury to nearby equipment or personnel that may be present while the device is hovering. If identified as a return-to-base behavior, the interference power adjustment module will... Increase to 150% of the original value, that is At the same time Shortened to 0.5 times the original value, that is This is to quickly reach maximum interference power and effectively block the drone's return path. If If the behavior is identified as disconnection, the interference power adjustment module will... Set to 200% of the original value, that is It also enables continuous wave modulation mode to continuously suppress the target communication link with maximum power. For straight-line flight and hovering behavior, the interference parameters maintain the default values ​​in the template. All adjusted parameters are encapsulated into a parameter structure. This structure is then passed to the subsequent waveform generation module.

[0048] Step S43: Combine the current electromagnetic environment background noise spectrum to dynamically calculate the optimal transmission power and radiation pattern of the interference signal, and generate a directional countermeasure signal waveform; In this embodiment, when dynamically calculating the optimal transmit power and radiation pattern of the interference signal based on the current electromagnetic environment background noise spectrum, the electromagnetic environment data in the current airspace is first collected in real time by the spectrum monitoring module. The sampling frequency is 100MHz and the sampling time is 10ms. The collected time-domain data is then subjected to a fast Fourier transform to obtain the background noise spectrum. The frequency resolution is 1kHz. The optimal transmit power is calculated based on the signal-to-interference-plus-noise ratio (SINR) constraint, aiming to minimize the transmit power while ensuring the jamming effect. The power of the jamming signal at the target UAV receiver is defined as... Background noise power is The signal-to-interference-plus-noise ratio threshold required for the interference effect is =10dB. Solve the equation ,in Let be the self-interference coefficient between the interfering signal and the target signal, set to 0.1. The minimum required interference power is then obtained. The beam pattern calculation is based on a beamforming algorithm for a uniform circular array with a radius of half a wavelength and eight elements. The desired beam pointing angle is... The location of the spatial peak in the three-dimensional cross-correlation response plot is determined. The steering vector is then calculated. Where d is the element spacing, For carrier wavelength. Beamforming weight vector The final generated directional countermeasure signal waveform is as follows: ,in This is the center frequency of the interference band.

[0049] Step S44: Perform pulse shaping and time slot alignment on the directional countermeasure signal waveform to ensure that it is synchronized with the target UAV's receiving window and generate a high-hit-rate jamming pulse sequence; In this embodiment, the waveform of the directional countermeasure signal is... When performing pulse shaping and time slot alignment, the receiving window period of the target UAV is first determined. This period is obtained by analyzing the frame structure in the time-frequency energy spectrum, measuring the time interval between two consecutive energy pulses, and taking the average value as the time. Pulse shaping uses a raised cosine filter with a roll-off factor α of 0.25 and a filter order L of 32. impulse response of raised cosine filter Perform convolution operation to obtain the shaped waveform. ,in The symbol period is used. The core of time slot alignment is adjusting the pulse transmission time to ensure it falls precisely within the target UAV's receiving window. The time offset between the current system clock and the target UAV's clock is calculated. This offset is obtained by performing preamble correlation detection on the received target signal. The shaped waveform is then time-shifted to obtain the aligned waveform. Subsequently, the continuous waveform is cut into a pulse sequence, with the duration of each pulse set to... The pulse interval is set to The cutting operation is as follows: For the nth pulse, extract... In the time interval The segment within is used as the waveform of the pulse. Arrange all pulse waveforms in chronological order to generate a high-hit-rate interference pulse sequence. Each pulse in the sequence is synchronized and aligned with the target receiving window.

[0050] Step S45: The high-hit-rate jamming pulse sequence is sent to the radio frequency power amplifier and radiated to the target airspace through the directional antenna array. At the same time, the transmission time, power, azimuth angle and frequency band information are recorded to generate a countermeasure operation log.

[0051] In this embodiment, a high-hit-rate interference pulse sequence is used. When the signal is fed into the RF power amplifier, the pulse waveform of the digital baseband is first... The signal is converted to an analog signal via a digital-to-analog converter (DAC) at a sampling rate of 1 GHz and a quantization bit depth of 14 bits. The analog signal then passes through an up-converter, up-converting it from the baseband frequency to the radio frequency carrier frequency. The upconversion uses quadrature modulation to convert the baseband... The two signals are multiplied by the cosine and sine components of the local oscillator respectively, and then added together to obtain the radio frequency signal. The radio frequency signal then enters the power amplifier, and the amplifier's gain... According to calculation and adjustments Configure the actual gain value. ,in This is the input power to the power amplifier. The amplified signal is transmitted to the directional antenna array via a feeder. The beam direction of the antenna array is controlled by a calculated weight vector w. The phase of each array element is adjusted by a phase shifter to achieve beam pointing. During the radiation process, the monitoring module records key parameters for each launch in real time: launch timestamp. Precision down to the microsecond level, transmission power The actual output power value and direction angle are obtained by sampling and detection through a coupler. Frequency band information is obtained from encoder readings of the antenna array. Record the start and end frequencies of the currently used interference band. All these records are written to the countermeasures operation log in chronological order. Each log entry is formatted as follows: This log provides complete and traceable data for subsequent evaluation of the effectiveness of countermeasures and operational audits.

[0052] Furthermore, step S43 includes the following steps: Step S431: Collect the spectral distribution of background noise in the current airspace, extract the signal-to-noise ratio and interference tolerance parameters of each frequency band, and generate an electromagnetic environment situation map; In this embodiment, when collecting the spectral distribution of background noise in the current airspace, a broadband spectrum monitoring receiver is first used to perform a frequency sweep measurement on the 20MHz to 6GHz band, with a frequency step set to 100kHz and a dwell time of 10 milliseconds for each frequency point. For each frequency point... A time-domain data segment with 1024 sampling points was collected, and its power spectral density was calculated. The PSD values ​​of all frequency points are stitched together to form a complete spectrum curve. The signal-to-noise ratio is calculated based on the target UAV's operating frequency band. Based on this, the frequency band is determined by the frequency band information in the Identity_tag generated in step S35. Within the frequency band, extract signal power Extracting noise power Signal-to-noise ratio The interference tolerance parameter is calculated based on the receiver's noise figure NF and the required demodulation threshold. Interference tolerance Repeat the above calculation for each frequency band to obtain a set of signal-to-noise ratio vectors corresponding to each frequency band. and interference tolerance vector These data are overlaid on a two-dimensional graph, with frequency on the horizontal axis and power density on the vertical axis. Color intensity represents the signal-to-noise ratio, generating an electromagnetic environment situation map. The situation map is stored in the form of an image matrix, where rows correspond to frequency indices, columns correspond to power density values, and the grayscale value of each pixel represents the signal-to-noise ratio level of that frequency.

[0053] Step S432: Based on the electromagnetic environment situation map and the receiver sensitivity curve of the target UAV, calculate the minimum effective power required for interference in each frequency band; In this embodiment, based on the electromagnetic environment situation map When calculating the minimum effective power required for interference in each frequency band based on the receiver sensitivity curve of the target UAV, the target UAV is first retrieved from the UAV model database. Corresponding receiver sensitivity curve This curve defines the minimum input power at which the receiver can correctly demodulate the signal at different frequencies. For each frequency band... ,from Extracting the background noise power of this frequency band Minimum effective power required for interference. The calculation formula is: ,in The required interference-to-signal ratio (ISR) threshold is set to 6 dB. The physical meaning of this formula is: when the interference signal reaches the input of the UAV receiver, its power must exceed the sum of the receiver sensitivity and background noise, plus the ISR margin, to effectively suppress the target signal. Performing the above calculation on all frequency bands yields a vector with a length equal to the total number of frequency bands. To ensure the interference effect, path loss factors are further considered, based on the predicted target distance. Calculate free space path loss ,in The carrier frequency is expressed in GHz. The unit is meters. The path loss is compensated to the minimum effective power to obtain the minimum transmit power required by the transmitter. This vector serves as the lower bound constraint for the subsequent power adjustment curve.

[0054] Step S433: Generate a dynamic power adjustment curve based on the power ramp-up rules preset by the minimum effective power and interference strategy module; In this embodiment, based on the minimum effective power When generating a dynamic power adjustment curve based on the preset power ramp-up rules of the interference strategy module, the preset power ramp-up rules are first extracted from the interference template. This rule is defined as a three-segment piecewise function. The first stage is a rapid ramp-up stage with a duration of T1 = 0.5 seconds, during which the power increases from... linear growth to The second stage is a stable adjustment stage, lasting for T2 = 1.0 second, with power at... The power fluctuates in a sinusoidal manner in the vicinity; the third stage is the steady-state maintenance stage, during which the power remains constant. The generation of the dynamic power adjustment curve requires combining preset rules with the calculated minimum effective power. Specifically, for the first stage, the actual transmit power... This means ensuring that the power output at any point during the climb does not fall below the minimum effective power. For the second stage, ,in The fluctuation range is set to 3dB. The modulation frequency is set to 2Hz. For the third stage, This generates a power curve that varies with time. The time resolution is 1 millisecond, and the total duration is consistent with the expected interference duration, set to 5 seconds. The curve is stored in array form with a length of 5000, and each element corresponds to the transmit power value at a specific time point.

[0055] Step S434: Combine the historical trajectory of the target UAV to predict the spatial position at the next moment, adjust the beam pointing of the directional antenna array, and generate a spatial pointing optimization vector; In this embodiment, when predicting the spatial location of the target UAV at the next moment by combining the target UAV's historical trajectory, the spatial location sequence of the target UAV within the past N time windows is first extracted from the behavioral evolution trajectory map. Where N is set to 10. A Kalman filter is used for state estimation and prediction, and the system's state vector is defined as... , including location and velocity components The state transition matrix F is defined for a constant velocity model, i.e. Where F is the state transition matrix, For process noise, the covariance matrix Q is set as Observation matrix Mapping the state vector to the observation space, i.e., the observation values. ,in For observation noise, the covariance matrix R is set as The prediction steps of Kalman filtering are as follows: The update steps are as follows: , , , where s pred To predict the state, p pred To predict the covariance, K is the Kalman gain, p est The updated covariance matrix, This is the updated state estimate, i.e., the state estimate at the current moment. After N iterations, the state estimate at the current moment is obtained. Then, a forward prediction is performed to obtain the state at the next time step. .from Extracting positional components This serves as a spatial position prediction for the next moment. Based on this predicted position, the azimuth angle relative to the center of the antenna array is calculated. and pitch angle Combining these two angles creates a spatial pointing optimization vector. .

[0056] Step S435: Perform joint constraint optimization on the dynamic power adjustment curve and the spatial pointing optimization vector to generate a set of disturbance parameters with joint space-power constraints; In this embodiment, the dynamic power adjustment curve is... With spatial pointing optimization vector When performing joint constraint optimization, the objective function and constraints of the optimization problem are first defined. The objective function is to maximize the equivalent isotropic radiated power of the interference signal at the target location. ,in For antenna array in direction The gain on the surface. The constraints include three aspects. The first is the power constraint, ,in It is determined by the adjusted maximum power. The second is the spatial pointing constraint, the beam pointing angle. and The deviation shall not exceed ±3 degrees, that is and The third constraint is time synchronization; power adjustment and beam pointing adjustment must be completed within the same time step, with a time step size of 10 milliseconds. The joint constraint optimization solution employs the alternating direction multiplier method, decomposing the problem into a power subproblem and a spatial pointing subproblem. In the power subproblem, the fixed beam pointing is... Solve Maximize EIRP. In the space-pointing subproblem, the fixed power is... Adjust the beam direction to make The maximum power curve is obtained after 10 alternating iterations. and optimal beam pointing angle Combine these two into a parameter set. This parameter set is the disturbance parameter set with space-power joint constraints.

[0057] Step S436: Input the set of interference parameters with space-power joint constraints into the waveform synthesizer, generate a time-domain baseband waveform through inverse Fourier transform, and then up-convert and modulate it to the target carrier frequency to generate a directional countermeasure signal waveform.

[0058] In this embodiment, the interference parameter set with space-power joint constraints is used. When input to the waveform synthesizer, it is first based on Optimal power curve and beam pointing angle Generate frequency domain waveforms. Spectrum template of the interference signal. Defined as in the target frequency band The frequency band consists of a complex sequence with an inner amplitude of 1 and a random phase, while the amplitude outside the band is set to zero. The length of the spectral template is set to 2048 points, corresponding to a frequency resolution of [missing information]. ,in The sampling rate is set to 100MHz. By applying power weighting, the weighted spectrum is obtained. Then, based on the beam pointing angle Calculate the steering vector of the antenna array Then, the spatial-frequency joint spectrum is obtained by outer productting it with the weighted spectrum. .right Performing the inverse Fourier transform yields the time-domain baseband waveform. The waveform has 2048 sampling points. After being converted into an analog signal by a digital-to-analog converter, the baseband waveform is sent to an up-converter module. The up-converter module combines the baseband signal with a carrier signal generated by a local oscillator. Multiply, where The target carrier frequency is used to obtain the radio frequency signal. This radio frequency signal is the final directional countermeasure signal waveform, and its power varies with time according to... Despite the changes, the beam direction remains consistently aligned with the predicted target location, achieving precise joint control of space and power.

[0059] Furthermore, step S435 includes the following steps: The dynamic power regulation curve is sliced ​​by time, and each slice corresponds to the power demand within a time window. In this embodiment, the dynamic power adjustment curve is... When slicing by time, first determine the total interference duration. Set to 5 seconds, time window length Set the time interval to 100 milliseconds and the overlap rate of adjacent windows to 0%. Calculate the total number of time slices accordingly. =50. For the first Each time slice has a time range of [number] time segments. .from Extract all power value sequences within the given time range from the array; the sequence length is [value missing]. ,in To optimize the time resolution of the power curve, we set it to 1 millisecond; therefore, each slice contains 100 power values. The arithmetic mean of these 100 power values ​​is then calculated. This represents the power demand for that slice. Simultaneously, the maximum power within that slice is calculated. and minimum value This serves as the constraint boundary for subsequent optimization. The calculation results of all slices are organized into a two-dimensional array. Its dimensions are ×3, each row contains three fields: This array comprehensively describes the temporal variation characteristics of power demand throughout the entire interference process. Each slice corresponds to an independent time window, and the power demand within the window is represented in the form of statistical parameters, providing discretized input data for subsequent combinatorial optimization.

[0060] Discretize the spatial pointing optimization vector according to angular resolution to generate a candidate beam pointing set; In this embodiment, the space is directed to the optimization vector. =[ , When discretizing by angular resolution, the azimuth angle is determined first. and pitch angle Dynamic range. The azimuth range is set to [-60°, 60°], and the elevation range is set to [-30°, 30°]. Angular resolution. Set to 1 degree, this means generating 121 discrete points in the azimuth direction and 61 discrete points in the elevation direction. Construct a two-dimensional grid covering the above-mentioned area, centered on [the grid point coordinates]. The calculation formula is: , where i ranges from -60 to 60; , where j ranges from -30 to 30. Each grid point corresponds to a candidate beam pointing direction. For each candidate direction... The gain of the antenna array in this direction is calculated in advance. Gain calculation is based on the pattern function of a uniform circular array, with the array radius set to half a wavelength and the number of array elements set to 8. All candidate directions and their corresponding gain values ​​are stored in a two-dimensional lookup table. Its dimensions are 121×61, and each element in the table is a tuple. This lookup table is the candidate beam pointing set, containing 7381 candidate directions, each with a corresponding antenna gain value, providing a complete spatial search space for subsequent combinatorial optimization.

[0061] The power demand and candidate beam pointing set within each time window are combined and optimized to find the power-pointing pairing scheme that maximizes interference efficiency. In this embodiment, when optimizing the power demand and candidate beam pointing set within each time window, the evaluation function for interference efficiency is first defined. Where P is the transmission power, For antenna gain, This represents path loss. The optimization objective is to maximize... It simultaneously satisfies two constraints. The first is a power constraint, whereby P must be within the specified time window. Within the interval. The second is the pointing constraint. Must be located in the predicted Within a ±3 degree range centered on the beam. The optimization solution employs a combination of exhaustive search and binary search. For each time window n, the search first proceeds from the candidate beam pointer to the set... Select a subset that satisfies the pointing constraint. The size of this subset is 7×7=49 candidate directions. Each candidate direction Calculate the corresponding optimal power .because Since it is a constant, the optimal power is obviously... Therefore, for each candidate direction, Select from among them The largest direction The optimal beam pointing for this time window corresponds to the optimal power. For all The optimization process is repeated over several time windows to obtain a set of optimized power-pointing pairing schemes. .

[0062] The optimized power-pointing pairing scheme is smoothly interpolated to generate a continuously changing time series of interference parameters; In this embodiment, the optimized power-pointing pairing scheme is used. When performing smooth interpolation, first obtain The start time of each time window At this time, the power sequence and angle sequence All values ​​are discrete, step-like values, and there may be abrupt changes between adjacent windows. The purpose of smooth interpolation is to generate a continuously changing time series and eliminate abrupt changes. Cubic spline interpolation is used to achieve this. As the independent variable, Construct a cubic spline function for the dependent variable. Constructing a cubic spline function requires solving a tridiagonal system of equations, with the boundary conditions set as natural boundaries, meaning the second derivative is zero at the endpoints. The interpolated power curve is shown below. The temporal resolution has been improved to 1 millisecond, compared to the original The resolution is consistent. Similarly, the azimuth sequence... and pitch angle sequence Perform cubic spline interpolation to obtain continuous angle curves. and During interpolation, the angle values ​​are seamlessly connected across the angle domain to avoid jumps caused by angles crossing the 360-degree boundary. The resulting continuously varying time series of disturbance parameters... The value of t ranges from 0 to... The sequence has a time resolution of 1 millisecond and a total length of 5000 time points. This sequence ensures a smooth transition in power and beam pointing over time.

[0063] The time series of interference parameters is converted into a sequence of instructions executable by the RF front end, generating a set of interference parameters with joint space-power constraints.

[0064] In this embodiment, the time series of interference parameters is... When converting to an executable instruction sequence for the RF front-end, the instruction format of the RF front-end is first defined. Each instruction contains four fields: timestamp. Power control word Phase shifter control word and switch control word The timestamp accuracy is 1 microsecond, based on the system startup time. The power control word is calculated based on... The mapping relationship with the digital-to-analog converter voltage is given by the mapping function. This linearly maps the power value to a 12-bit binary number ranging from 0 to 4095. The calculation of the phase shifter control word is based on... and The phase shifters are configured with a phase step of 5.625 degrees, corresponding to 6 bits of binary data. A total of 48 bits of control words are needed for the 8 array elements. The switch control word, consisting of 8 bits, controls the on / off state of each array element, with each bit corresponding to one element. The parameters at each time point are converted into the above instruction format to generate an instruction sequence. Since the time resolution is 1 millisecond, while the RF front-end's instruction execution cycle is 100 microseconds, the instruction sequence needs to be time-aligned. This involves copying the parameter values ​​within each 1-millisecond interval 10 times to match the instruction sequence's time resolution with the RF front-end's execution cycle. The final generated instruction sequence... The sequence consists of 50,000 instructions, each containing complete power, phase, and switching control information. This instruction sequence is then packaged into a set of disturbance parameters with joint space-power constraints. This parameter set can be directly sent to the RF front-end hardware for execution.

[0065] Furthermore, the present invention also provides a low-altitude unmanned aerial vehicle (UAV) signal analysis system based on multi-source features, including a processor, a memory, and a computer program stored in the memory and executable on the processor, for performing the low-altitude UAV signal analysis method based on multi-source features as described above.

[0066] The above description is merely a specific embodiment of the present invention, enabling those skilled in the art to understand or implement the invention. Various modifications to these embodiments will be readily apparent to those skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of the invention. Therefore, the present invention is not to be limited to the embodiments shown herein, but is to be accorded the widest scope consistent with the principles and novel features of the invention herein.

Claims

1. A method for analyzing low-altitude unmanned aerial vehicle (UAV) signals based on multi-source features, characterized in that, Includes the following steps: Step S1: Control the multi-source sensor array to perform full polarization scanning of the low-altitude target airspace, acquire the original radio frequency signal stream containing amplitude, phase and time-frequency characteristics, and generate multi-channel original signal tensors; Step S2: Perform nonlinear decoupling and channel estimation on the multi-channel original signal tensor to obtain the time-frequency energy spectrum through the multi-channel original signal tensor, extract the channel state information vector corresponding to each channel, and construct a space-time-frequency three-dimensional feature map based on the channel state information vector to generate a space-time-frequency feature cube; wherein, the extraction process of the channel state information vector specifically includes: Step S241: Extract the energy peak value corresponding to each subcarrier and its time window position from the time-frequency energy spectrum to generate a set of energy peak position pairs; Step S242: Construct an energy propagation path graph between subcarriers based on the set of energy peak positions. The nodes in the path graph are subcarrier indices, and the weight of the edges is the time difference of the energy peak of adjacent subcarriers, generating an energy propagation topology graph. Step S243: Perform graph neural network inference on the energy propagation topology graph, extract the neighbor aggregation features of each node, and generate the topology embedding vector; Step S244: Channel splicing of the topology embedding vector with the original time-frequency energy spectrum to generate an enhanced energy spectrum that fuses topology information; Step S245: Perform complex exponential fitting on the enhanced energy spectrum subcarrier by subcarrier, extract the channel fading factor and phase rotation factor corresponding to each subcarrier, and generate a channel state information vector; Step S3: Input the spatiotemporal feature cube into the intent perception engine, call the preset behavior prototype in the UAV behavior pattern library, extract the UAV identity fingerprint tag and behavior intent vector corresponding to the behavior prototype, and generate a behavior intent candidate set; Step S4: Activate the corresponding interference strategy module according to the candidate set of behavioral intentions, dynamically calculate the optimal transmission power and radiation pattern of the interference signal in combination with the background noise spectrum of the current electromagnetic environment, generate the directional countermeasure signal waveform, and perform pulse shaping and time slot alignment on the directional countermeasure signal waveform to generate a high hit rate interference pulse sequence. Step S5: The high-hit-rate jamming pulse sequence is sent to the radio frequency power amplifier and radiated to the target airspace through the directional antenna array. At the same time, the transmission time, power, azimuth angle and frequency band information are recorded to generate a countermeasure operation log.

2. The low-altitude UAV signal analysis method based on multi-source features according to claim 1, characterized in that, Step S2 includes the following steps: Step S21: Perform complex domain decomposition on the multi-channel original signal tensor, separate the real and imaginary components, and generate a complex signal matrix; Step S22: Calculate the instantaneous amplitude and instantaneous phase of each channel based on the complex signal matrix, construct the amplitude-phase joint distribution map, and generate the amplitude-phase feature map; Step S23: Perform time-frequency transformation on the amplitude-phase feature map, extract the time-varying spectral energy distribution, and generate a time-frequency energy spectrum; Step S24: Perform channel estimation based on the time-frequency energy spectrum, calculate the channel response vector of each subcarrier, and generate the channel state information vector; Step S25: Map the channel state information vector to a space-time-frequency three-dimensional grid to construct a space-time-frequency feature cube, where the spatial dimension corresponds to the antenna array index, the time dimension corresponds to the sampling time, and the frequency dimension corresponds to the subcarrier index.

3. The low-altitude UAV signal analysis method based on multi-source features according to claim 1, characterized in that, Step S243 includes the following steps: Each node in the energy propagation topology is initialized with the average energy value of that subcarrier, and an initial feature vector of the node is generated. Perform message passing operations on each node, aggregate the feature vectors of its first-order neighbor nodes and sum them with weights to generate neighbor aggregate features; The initial feature vector of a node is fused with the aggregated features of its neighbors using a gating system. The ratio between the two is dynamically adjusted through a learnable gating unit to generate a gating update feature vector. A nonlinear activation function is applied to the gated update feature vector, and residual connections are introduced to preserve the original information, thereby generating a deep graph feature vector; Stack the deep graph feature vectors of all nodes in subcarrier index order to generate a topology embedding vector.

4. The low-altitude UAV signal analysis method based on multi-source features according to claim 2, characterized in that, Step S25 includes the following steps: Step S251: Normalize the channel state information vector so that its mean is zero and its standard deviation is one, and generate a standardized channel vector; Step S252: Calculate the spatial correlation matrix based on the standardized channel vector, extract the principal component directions in the spatial dimension, and generate the spatial principal component vector; Step S253: Perform tensor product operation between the spatial principal component vector and the time and frequency dimensions to generate a high-dimensional feature tensor; Step S254: Dimensionality reduction and compression are performed on the high-dimensional feature tensor, and low-dimensional hidden layer representation is extracted using an autoencoder to generate compact spatiotemporal frequency feature representation; Step S255: Reshape the compact spatiotemporal feature representation into a spatiotemporal three-dimensional mesh structure to generate a spatiotemporal feature cube.

5. The low-altitude UAV signal analysis method based on multi-source features according to claim 1, characterized in that, Step S3 includes the following steps: Step S31: Input the spatiotemporal frequency feature cube into the intent perception engine and call the preset behavior prototypes in the UAV behavior pattern library, including hovering, straight flight, circling, return to home and loss of contact behavior; Step S32: Extract the corresponding signal feature signature for each behavior prototype, including phase jitter mode, spectral broadening trend and frame interval stability parameters, and generate a behavior prototype feature dictionary; Step S33: Perform sliding window matching between the feature points in the space-time-frequency feature cube and the behavioral prototype feature dictionary, and calculate the local matching degree score; Step S34: Construct a temporal behavior chain based on the local matching degree score, identify the evolution path of behavior patterns in consecutive frames, and generate a behavior evolution trajectory map; Step S35: Based on the continuity and stability of the path in the behavior evolution trajectory map, output the final UAV identity fingerprint tag and behavior intention vector. The identity fingerprint tag includes the model, firmware version and encryption protocol type, and the behavior intention vector includes the current action category and the prediction probability of the next action, thereby generating a candidate set of behavior intentions.

6. The low-altitude UAV signal analysis method based on multi-source features according to claim 1, characterized in that, Step S4 includes the following steps: Step S41: Activate the corresponding interference strategy module according to the behavioral intent vector. The module is pre-set with interference templates for different types of UAVs, including frequency band avoidance strategies, power ramp-up curves and modulation mode switching rules. Step S42: Adjust the interference parameters based on the behavior intent vector, including reducing the interference power to avoid accidental damage when the behavior is identified as hovering, and increasing the interference to block the return path when the behavior is identified as returning. Step S43: Combine the current electromagnetic environment background noise spectrum to dynamically calculate the optimal transmission power and radiation pattern of the interference signal, and generate a directional countermeasure signal waveform; Step S44: Perform pulse shaping and time slot alignment on the directional countermeasure signal waveform to ensure that it is synchronized with the target UAV's receiving window and generate a high-hit-rate jamming pulse sequence; Step S45: The high-hit-rate jamming pulse sequence is sent to the radio frequency power amplifier and radiated to the target airspace through the directional antenna array. At the same time, the transmission time, power, azimuth angle and frequency band information are recorded to generate a countermeasure operation log.

7. The low-altitude UAV signal analysis method based on multi-source features according to claim 6, characterized in that, Step S43 includes the following steps: Step S431: Collect the spectral distribution of background noise in the current airspace, extract the signal-to-noise ratio and interference tolerance parameters of each frequency band, and generate an electromagnetic environment situation map; Step S432: Based on the electromagnetic environment situation map and the receiver sensitivity curve of the target UAV, calculate the minimum effective power required for interference in each frequency band; Step S433: Generate a dynamic power adjustment curve based on the power ramp-up rules preset by the minimum effective power and interference strategy module; Step S434: Combine the historical trajectory of the target UAV to predict the spatial position at the next moment, adjust the beam pointing of the directional antenna array, and generate a spatial pointing optimization vector; Step S435: Perform joint constraint optimization on the dynamic power adjustment curve and the spatial pointing optimization vector to generate a set of disturbance parameters with joint space-power constraints; Step S436: Input the interference parameter set with space-power joint constraints into the waveform synthesizer, generate the time-domain baseband waveform through inverse Fourier transform, and then up-convert and modulate it to the target carrier frequency to generate the directional countermeasure signal waveform.

8. The low-altitude UAV signal analysis method based on multi-source features according to claim 7, characterized in that, Step S435 includes the following steps: The dynamic power regulation curve is sliced ​​by time, and each slice corresponds to the power demand within a time window. Discretize the spatial pointing optimization vector according to angular resolution to generate a candidate beam pointing set; The power demand and candidate beam pointing set within each time window are combined and optimized to find the power-pointing pairing scheme that maximizes interference efficiency. The optimized power-pointing pairing scheme is smoothly interpolated to generate a continuously changing time series of interference parameters; The time series of interference parameters is converted into a sequence of instructions executable by the RF front end, generating a set of interference parameters with joint space-power constraints.

9. A low-altitude unmanned aerial vehicle (UAV) signal analysis system based on multi-source features, characterized in that, The system includes a processor, a memory, and a computer program stored in the memory and executable on the processor, for performing the low-altitude unmanned aerial vehicle signal analysis method based on multi-source features as described in any one of claims 1-8.

Images