Lift-off platform unmanned aerial vehicle signal detection and accurate positioning method

Abstract

The invention relates to a signal detection and accurate positioning method for a lift-off platform unmanned aerial vehicle, and the method comprises the steps: deploying an unmanned aerial vehicle-mounted platform, capturing a global electromagnetic signal in real time, outputting multichannel original signal data, carrying out the time-frequency domain joint analysis, extracting a frequency hopping period, a modulation mode and pseudo-random code features, and recognizing a target unmanned aerial vehicle signal. A control instruction link signal and an image transmission link signal are separated; the arrival time difference and the arrival phase difference of the target unmanned aerial vehicle signal are converted into a dynamic space coordinate system, a three-dimensional positioning equation set is constructed, and real-time space position coordinates of the target unmanned aerial vehicle are obtained through calculation; and extracting transmitting end characteristic parameters of control instruction link signals, constructing a control signal propagation path model in combination with the real-time position of the target unmanned aerial vehicle, and positioning geographic coordinates of an unmanned aerial vehicle control source by reversely tracking a signal transmitting direction and an attenuation law. The terrain shielding influence is effectively avoided, the target identification accuracy and positioning precision are improved, and the unmanned aerial vehicle control source traceability is realized.

Description

Technical Field

[0001] This invention belongs to the field of unmanned aerial vehicle (UAV) signal detection technology, and particularly relates to a method for signal detection and precise positioning of UAVs on airborne platforms. Background Technology

[0002] Existing UAV signal detection and positioning technologies mainly include radio frequency signal-based detection and positioning methods, vision-based detection and positioning methods, and acoustic-based detection and positioning methods. The existing UAV signal detection and positioning technologies are as follows: Radar detection technology: Radar detects drones by emitting electromagnetic waves and receiving the echoes reflected back from targets. When the electromagnetic waves emitted by the radar encounter a drone, some of the waves are reflected back to the radar. The radar determines the drone's distance, speed, and azimuth based on parameters such as the time delay and frequency changes of the echo. However, in urban environments, reflections and scattering from buildings and other objects generate significant clutter interference, affecting the accurate identification and location of drone signals. Furthermore, radar equipment is typically large, expensive, and complex to deploy and maintain. Acoustic wave recognition technology utilizes the acoustic characteristics generated by drones during flight for detection and location. During flight, the rotation of the drone's propellers and the operation of its motors produce sound waves with specific frequencies and characteristics. By deploying multiple microphone arrays, acoustic signals from the environment are collected, and then these signals are analyzed and processed to extract acoustic features associated with the drone, thereby determining its presence and location. However, acoustic wave recognition technology is susceptible to environmental noise. In noisy environments, ambient noise may mask the acoustic signals generated by the drone, leading to a decrease in the accuracy of detection and location. Furthermore, sound waves are affected by distance and obstacles during propagation, limiting their propagation distance and making it less effective at detecting drones at long distances. Radio reception and detection technology: This technology uses equipment such as spectrum analyzers to monitor and analyze radio signals in the air, identifying characteristics such as the frequency, modulation method, and signal strength of drone communication signals, thereby determining the presence and location of the drone. However, in complex electromagnetic environments, such as city centers or near communication base stations, there is a large amount of interference from other radio signals, which can easily lead to misjudgment and missed detection of drone signals. In addition, some advanced drones may use encrypted communication technology, increasing the difficulty of signal identification and decryption.

[0003] Infrared detection technology: During flight, drones generate heat from their engines, motors, and other components, emitting infrared radiation. Infrared detectors capture these infrared radiation signals to detect the presence and location of drones. However, infrared detection technology has a relatively short detection range and is greatly affected by weather conditions. In adverse weather conditions such as rain and fog, the propagation of infrared signals is severely attenuated, leading to a decrease in detection effectiveness. Furthermore, its ability to identify drones is limited, making it difficult to accurately distinguish between different drone models. Therefore, this invention aims to achieve full-domain three-dimensional capture of UAV signals by constructing a three-dimensional detection network, thereby overcoming the problem of signal interference in complex environments and improving the reliability of signal detection. Summary of the Invention

[0004] The purpose of this invention is to provide a method for signal detection and precise positioning of unmanned aerial vehicles (UAVs) on airborne platforms, in order to solve the technical problems existing in the prior art.

[0005] To solve the above-mentioned technical problems, the technical solution adopted by the present invention is as follows: A method for signal detection and precise positioning of an aerial platform UAV includes the following steps: S1: Deploy multiple UAV platforms equipped with detection equipment, including multi-band array antennas with multiple subarrays, to capture electromagnetic signals across the entire domain in real time and output multi-channel raw signal data; S2: Perform joint time-frequency domain analysis on the multi-channel raw signals, extract the frequency hopping period, modulation mode and pseudo-random code features of the signals, match them with the preset UAV signal feature library to identify the target UAV signals, and separate the control command link and image transmission link signals of the target signals; S3: Establish a dynamic spatial coordinate system, transform the time difference of arrival (TDoA) and phase difference of arrival (PDoA) of the identified target UAV signal between multiple subarrays to this coordinate system, construct a three-dimensional positioning equation set, and solve it to obtain the real-time spatial position coordinates of the target UAV. S4: For the control command link signal, extract its transmitter characteristic parameters, combine them with the real-time spatial coordinates of the target UAV, establish a control signal propagation path model, and locate the geographical coordinates of the UAV control source by reverse tracing the signal transmission direction and attenuation law.

[0006] Preferably, the specific process of step S1 is as follows: S11: Deployment of the launch platform: Deploy the launch equipment equipped with a multi-band adaptive array antenna to the preset altitude, complete the platform attitude calibration through the ground control system, and complete the initial spatial coordinate anchoring based on GNSS positioning data; S12: Array Antenna Activation and Frequency Band Configuration: Activates the multi-band adaptive array antenna, which contains 3 or more subarrays, and achieves seamless coverage of the 1.2-6GHz frequency band through the frequency band overlap design between the subarrays; S13: Beamforming coordination: Each subarray independently generates an adjustable beam through phased array technology. The ground control unit dynamically adjusts the beam pointing angle and gain distribution between subarrays based on real-time electromagnetic environment monitoring data. S14: Signal Acquisition and Multi-channel Output: The array antenna identifies the direction of the incoming signal through spatial spectrum estimation, drives the corresponding subarray beam to focus on the direction of the incoming signal to enhance the receiving gain, and outputs it to the signal processing unit in real time according to the subarray channel number to form a multi-channel raw dataset containing timestamp, signal strength and phase information.

[0007] Preferably, the specific process of performing time-frequency domain joint analysis on the multi-channel original signal in step S2 to extract the frequency hopping period, modulation scheme, and pseudo-random code features is as follows: S20: Multi-channel signal preprocessing: Synchronous calibration of the original digital signals of the multi-channel output of each subarray, and filtering out environmental noise and co-frequency interference by wavelet threshold denoising. S21: Joint Time-Frequency Domain Analysis: Perform a short-time Fourier transform on the preprocessed signal to generate a time-frequency distribution map, identify the energy distribution characteristics of the signal in the time-frequency plane, and locate the frequency jump time and corresponding carrier frequency value of the frequency hopping signal. For the steady-state signal segment in the time-frequency distribution map, transient features are extracted by continuous wavelet transform, and instantaneous frequency and amplitude are calculated by combining Hilbert transform to determine the frequency hopping period; S22: Modulation mode identification: For non-frequency hopping signal segments or single-frequency dwell segments of frequency hopping signals, extract the amplitude histogram, phase change rate, and spectral correlation features of the signal; The modulation scheme is identified by matching features using a support vector machine classifier. S23: Pseudo-random code feature extraction: For signals identified as spread spectrum modulation, the pseudo-random code period is determined by sliding correlation calculation; The code sequence features are extracted based on the signal phase reversal points within the code period, and compared with common UAV pseudo-random codes in a preset feature library to determine the code type and code length.

[0008] Preferably, the specific process of matching the extracted features with a preset UAV signal feature library in step S2 to identify the target UAV signal is as follows: S24: Standardization Processing and Feature Library Call: Standardize the extracted target signal features; call the preset UAV signal feature library, which contains a set of baseline features for known UAV models, with each parameter associated with a corresponding weight coefficient; S25: Multi-dimensional feature matching calculation: For frequency hopping period: calculate the deviation rate between the target signal frequency hopping period and the reference period of each model in the feature library. Deviation rate = |target value - reference value| / reference value. When the deviation rate is ≤5%, it is judged as a preliminary match. Regarding modulation methods: A type matching method is used. If the modulation method of the target signal is exactly the same as a certain model in the feature library, the score is 100% of the corresponding weight score of that parameter. If it is a compatible type, the score is 60% of the corresponding weight score. S26: Comprehensive matching degree evaluation and target recognition: The multi-dimensional matching scores above are weighted and summed to calculate the comprehensive matching degree between the target signal and each UAV model in the feature library, and the UAV model is matched according to the comprehensive matching degree.

[0009] Preferably, the specific process of separating the control command link and image transmission link signal of the target signal in step S2 is as follows: S27: Preset feature library for distinguishing between control command link and image transmission link; S28: Signal framing and parameter extraction: Perform frame synchronization processing on the identified target UAV signal, segment the continuous signal into independent data frames; extract the key parameters of each frame signal; S29: Link type classification and separation: The decision tree classification model is used to determine the link type of the framed signal.

[0010] Preferably, the specific process of establishing a dynamic spatial coordinate system based on the real-time GNSS positioning data of the launch platform and the attitude parameters output by the inertial navigation system in step S3 is as follows: S31: Reference coordinate anchoring: Taking the center of mass of the launch platform as the origin of the coordinate system O, the geodetic coordinates (L, B, H) of the platform are obtained in real time through the GNSS module. L is longitude, B is latitude, and H is altitude. The coordinates are converted into absolute coordinates (X0, Y0, Z0) in the WGS-84 geocentric rectangular coordinate system and used as the spatial anchor point of the dynamic coordinate system. Set the initial orientation of the coordinate system: with the origin O as the starting point, the vector pointing from the Earth's center to point O is the Z-axis. The geodetic coordinates are converted into planar coordinates through Gaussian projection to determine the initial east-west axis (X-axis) and north-north axis (Y-axis), forming the initial ENU coordinate system framework; S32: Attitude parameter fusion and dynamic coordinate system correction: The attitude parameters output by the inertial navigation system (INS) are acquired in real time, including pitch angle θ, roll angle φ, and heading angle ψ; a rotation matrix is ​​constructed based on the attitude angles: according to the Euler angle transformation formula, the initial ENU coordinate system is rotated around the Z-axis, Y-axis, and X-axis in sequence to generate a dynamic rotation matrix; The initial ENU coordinate system is corrected in real time by using a rotation matrix: the axis vectors (X0, Y0, Z0) of the initial ENU coordinate system are multiplied by the rotation matrix to obtain the corrected instantaneous coordinate axis directions (X', Y', Z'). S33: Time synchronization and coordinate update mechanism: The PPS signal of GNSS is used as the time reference to synchronize the attitude parameters output by INS with the GNSS positioning data using timestamps. Perform a coordinate system update at specified intervals: combine the latest GNSS absolute coordinates (X0', Y0', Z0') and INS attitude angles (ψ', θ', φ'), recalculate the rotation matrix, correct the position of the coordinate system origin and the axis orientation, and compensate for platform drift.

[0011] Preferably, the specific process in S3 of transforming the arrival time difference and arrival phase difference of the target UAV signal among the multiple subarrays to this coordinate system, constructing a three-dimensional positioning equation set, and solving it to obtain the real-time spatial position coordinates of the target UAV is as follows: S34: Multi-subarray spatial coordinate calibration: In the established dynamic spatial coordinate system, the relative coordinates of each subarray antenna are pre-calibrated: let the coordinates of subarray i be (xᵢ, yᵢ, zᵢ), i=1, 2, ..., n, n≥3. The physical distance between subarrays is obtained by laser ranging and converted into three-dimensional coordinate values ​​in the coordinate system to form the subarray position matrix. S35: Spatiotemporal synchronization conversion between TDoA and PDoA: TDoA Conversion: For target signals received by multiple subarrays, using subarray 1 as a reference, calculate the arrival time difference Δtᵢ1 between the remaining subarrays i and subarray 1. Combining this with the electromagnetic wave propagation speed c, convert it into a distance difference Δdᵢ1 = Δtᵢ1 × c, resulting in the relationship: [(X-xᵢ)² + (Y-yᵢ)² + (Z-zᵢ)²] 1 / 2 -[(X-x1)²+(Y-y1)²+(Z-z1)²] 1 / 2 =Δdᵢ1, where (X, Y, Z) are the coordinates of the target in the dynamic coordinate system to be solved; PDoA Conversion: For signals with the same carrier frequency f, calculate the phase difference Δφᵢ1 between subarray i and subarray 1, and convert it into the distance difference Δd'ᵢ1=(Δφᵢ1×λ) / (2π), where λ is the signal wavelength and λ=c / f. Based on the baseline length bᵢ1 between subarrays, the supplementary relationship is obtained: [(X-xᵢ)²+(Y-yᵢ)²+(Z-zᵢ)²]1 / 2 -[(X-x1)²+(Y-y1)²+(Z-z1)²] 1 / 2 =Δd'ᵢ1+kᵢ1×λ, where kᵢ1 is an integer, and the range of values ​​for kᵢ1 is constrained by the TDoA result.

[0012] S36: Construction and Solution of Three-Dimensional Positioning Equations Select three or more subarray combinations and construct a nonlinear equation system containing three unknowns (X, Y, Z) based on the distance difference relationship between TDoA and PDoA conversions mentioned above. The nonlinear equation system is linearized using the Taylor series expansion method: taking the initial position estimate of the target as the starting point of the iteration, the Jacobian matrix is ​​calculated and the equation system is transformed into the linear form Ax=b; The linear equations are solved by the least squares method to obtain the target coordinate correction. The correction is iteratively updated until the correction is ≤0.1m. The converged (X, Y, Z) is output as the real-time spatial coordinates of the target in the dynamic coordinate system.

[0013] S37: Coordinate system transformation and absolute position output: Transform the target coordinates (X, Y, Z) in the dynamic coordinate system to the WGS-84 geocentric rectangular coordinate system: using the rotation matrix R in step S3. -1 The attitude is calculated using the inverse of the rotation matrix to eliminate the influence of the platform attitude on the coordinates and obtain the absolute offset of the target relative to the origin of the launch platform. By superimposing the GNSS absolute coordinates (X0, Y0, Z0) of the launch platform, the WGS-84 geodetic coordinates (L, B, H) of the target UAV are obtained, completing the closed loop of the solution from signal time difference / phase difference to geographic coordinates.

[0014] Preferably, the specific process of step S4 is as follows: S41: Extraction of characteristic parameters from the control command link transmitter: Carrier frequency drift rate extraction: Carrier frequency tracking is performed on the control command link signal. The carrier frequency is continuously collected within a specified time period, and the carrier frequency drift rate over time is calculated through linear fitting. Power fluctuation period extraction: The received signal power is sampled in real time, and the power spectral density is analyzed by Fourier transform to identify the frequency value corresponding to the main peak. Its reciprocal is the power fluctuation period, and the fluctuation amplitude range is recorded at the same time. Encoding format parsing: Decode the control command frame data, extract the frame header structure, check method, data frame length and bit synchronization method, and generate encoding feature code; S42: Construction of the control signal propagation path model: Basic path model: Using the calculated real-time UAV position coordinates (X1, Y1, Z1) as relay points, assuming the control signal travels directly from the control source (X0, Y0, Z0) to the UAV, and combining the relative position of the UAV and the launch platform (X... p Y p Z p Establish a two-stage propagation model: control source → UAV → launch platform; Attenuation correction: Based on the received power P of the control command signal r Combined with the free space decay formula P r =P t -20lg(4πd / λ), where Pt is the transmit power, d is the propagation distance, and λ is the signal wavelength. An environmental attenuation coefficient is introduced, and terrain occlusion information scanned by the lidar of the launch platform is used for real-time correction. The occlusion area coefficient increases by 5-20dB, and the propagation distance d from the control source to the UAV is calculated. 01 ; S43: Control source direction reverse tracing: Angle of arrival calculation: Based on the PDoA measurement value of the control command signal from the multi-subarray antenna of the launch platform, combined with dynamic spatial coordinate system transformation, the angle of arrival, azimuth α1 and elevation β1 of the signal from the UAV to the launch platform are obtained; Reverse derivation: Based on the UAV's position (X1, Y1, Z1) and motion state (speed, heading), calculate the incident direction when the control signal reaches the UAV. Assuming the UAV is an omnidirectional receiver, the azimuth angle α0 of the control source relative to the UAV is α1 + 180° ± Δα, where Δα is the UAV attitude correction amount; the pitch angle β0 is calculated by the difference between the signal incident altitude angle and the UAV's altitude. Trajectory intersection localization: Using the UAV's position as the vertex, construct rays along the reverse directions α0 and β0, and combine this with the inversely calculated propagation distance d. 01 The candidate coordinates (X0', Y0', Z0') of the control source at that moment are obtained. S44: Multi-time data fusion and coordinate calculation: The coordinates of candidate control sources are continuously collected at 5-10 time intervals of 0.5-2 seconds, synchronized with the UAV positioning update cycle; Kalman filtering is used to fuse candidate coordinates at multiple time points: with the control source moving speed ≤10m / s as the state constraint, coordinate drift is dynamically corrected to obtain the filtered control source coordinates (X0, Y0, Z0). Elevation constraints are applied using ground digital maps, with Z0 typically being ±5m above ground level. The final output is the WGS-84 geodetic coordinates of the control source.

[0015] Preferably, the process also includes the following: Multiple UAV platforms equipped with detection equipment are deployed to form a dynamic collaborative network by combining a preset frequency domain / time domain / code domain joint communication protocol with dynamic spectrum sensing technology, thereby constructing a three-dimensional detection network that can avoid terrain obstruction such as buildings and mountains. The three-dimensional detection network is used to capture various signals of the target UAV in a three-dimensional manner to obtain multi-target signals; the time difference arrival data and frequency difference arrival data of each UAV platform are collected simultaneously, as well as the real-time positioning data generated by the inertial navigation unit on each platform, to form a multi-source dataset. After preprocessing the multi-source dataset, the target UAV is located using a distributed positioning algorithm. The multi-target signals are imported into a signal processing unit based on an FPGA hardware acceleration architecture. First, the signals are quickly sorted through hardware acceleration. Then, the frequency hopping / spread spectrum signal features are extracted to form real-time multi-target signal processing data that is coordinated with the positioning results, supporting subsequent source tracing and positioning. Based on the real-time processing data of the multi-target signals, the control source of the target UAV is identified by matching signal fingerprint features. At the same time, the arrival angle data of the multi-UAV in the three-dimensional detection network is called and the geospatial data of the geographic information system is integrated to finally realize the source tracing and positioning of the control source.

[0016] The beneficial effects of this invention include: 1. By deploying multi-band array antennas with three or more subarrays and combining them with the collaborative networking of multiple UAV platforms, a seamless coverage network for the commonly used frequency bands of UAVs in the 1.2-6GHz range can be built, solving the frequency blind spot problem of traditional single-band equipment and capturing various types of signals such as control commands and image transmissions from different UAV models.

[0017] 2. By extracting frequency hopping period, modulation scheme, and pseudo-random code features through joint time-frequency domain analysis, and combining this with multi-dimensional matching, accurate identification of UAV models can be achieved. Compared with traditional single-feature identification methods, the accuracy of model identification is improved, and the error in classifying unknown signals is reduced.

[0018] 3. A dynamic spatial coordinate system is established through GNSS and INS fusion, updating the coordinates every 50ms to compensate for platform drift. The coordinate system axis pointing and platform attitude synchronization deviation is ≤0.5°. Compared with a fixed coordinate system, this eliminates positioning reference errors caused by platform tilt and jitter, providing a stable spatial framework for TDoA / PDoA conversion. Strong positioning equation constraints and small solution errors: Combining the dual constraints of TDoA and PDoA, a three-dimensional nonlinear equation system is constructed. Through Taylor series expansion and least squares iteration, the target UAV positioning error is small, improving accuracy compared to single TDoA positioning. Even if GNSS briefly loses lock, positioning continuity is maintained through INS track extrapolation, reducing positioning interruption time to less than 0.5s. Distributed positioning and computation enables parallel processing of multi-target signals. A single platform can complete positioning updates for 10 targets per second, improving processing efficiency by 3 times compared to centralized computation, adapting to scenarios where multiple UAVs simultaneously intrude.

[0019] 4. A two-segment propagation model is established using the UAV's location as a relay point. Combining the free-space attenuation formula and lidar terrain correction, the propagation distance from the control source to the UAV is calculated in reverse. Compared to the model without terrain correction, the accuracy of distance calculation is improved. Through PDoA measurement and reverse direction derivation, a ray is constructed based on the propagation distance to obtain the candidate coordinates of the control source. Then, through Kalman filtering and multi-time fusion, the source tracing error is reduced. Compared to the traditional signal direction of arrival + ground investigation method, the tracing efficiency is effectively improved. Attached Figure Description

[0020] Figure 1 This is a flowchart illustrating the method for signal detection and precise positioning of an unmanned aerial vehicle (UAV) on an airborne platform according to the present invention.

[0021] Figure 2 This is a schematic diagram of the architecture of the multi-band array antenna of the present invention. Detailed Implementation

[0022] The following is in conjunction with the appendix Figures 1-2 The present invention will be further described in detail below: Example 1 See appendix Figure 1 As shown, a method for signal detection and precise positioning of an aerial platform UAV includes the following steps: S1: Deploy multiple UAV platforms equipped with detection equipment. The detection equipment includes multi-band array antennas with three or more subarrays. Each subarray covers the commonly used UAV frequency bands from 1.2GHz to 6GHz, and the frequency band overlap between adjacent subarrays is no less than 15%. The spatial spectrum estimation capability is formed through beamforming collaboration between subarrays, and the full-domain electromagnetic signals are captured in real time and multi-channel raw signal data is output. S2: Perform joint time-frequency domain analysis on the multi-channel raw signals, extract the frequency hopping period, modulation mode and pseudo-random code features of the signals, match the extracted features with the preset UAV signal feature library to identify the target UAV signals, and separate the control command link and image transmission link signals of the target signals to provide a layered signal source for positioning and tracing. S3: Based on the real-time GNSS positioning data of the launch platform and the attitude parameters output by the inertial navigation system, a dynamic spatial coordinate system is established. The time difference of arrival (TDoA) and phase difference of arrival (PDoA) of the identified target UAV signal between multiple subarrays are transformed into this coordinate system. A three-dimensional positioning equation system is constructed and the real-time spatial position coordinates of the target UAV are obtained by solving it. S4: For the control command link signal, extract its transmitter characteristic parameters, including carrier frequency drift rate, power fluctuation period and encoding format. Combined with the real-time spatial coordinates of the target UAV, establish a control signal propagation path model. By tracing the signal transmission direction and attenuation law in reverse, locate the geographical coordinates of the UAV control source.

[0023] In another embodiment of this example, UAV signal detection and positioning are achieved through the following process: Multiple UAV platforms equipped with detection equipment are deployed, and dynamic collaborative networking is achieved through a pre-defined frequency domain / time domain / code domain joint communication protocol combined with dynamic spectrum sensing technology. This constructs a three-dimensional detection network that can avoid obstructions from buildings, mountains, and other terrain features. The detection equipment includes multi-band receiving antennas, spectrum sensing modules, and communication units. The UAVs are initially deployed in a layered high- and low-altitude + regional grid layout. Low-altitude UAVs can cover an altitude of 10-100 meters, and high-altitude UAVs can cover an altitude of 100-500 meters, ensuring that signal blind spots caused by terrain obstruction are avoided.

[0024] The three-dimensional detection network is used to capture various signals of the target UAV in a three-dimensional manner to obtain multi-target signals; the time difference arrival (TDOA) data and frequency difference arrival (FDOA) data of each UAV platform are collected simultaneously, as well as the real-time positioning data generated by the high-precision inertial navigation unit (IMU+RTK positioning module) on each platform, to form a multi-source dataset. After preprocessing the multi-source dataset, the target UAV is located using a distributed positioning algorithm: based on the preprocessed multi-source data, the UAV platform acts as a local node, forming a local observation cluster with 3-5 neighboring UAV nodes to complete the initial positioning in parallel. Construct a system of positioning equations: TDOA constraint: Let the target coordinates be ( x,y,z ),node i With nodes j The coordinates are ( x i,y i ,z i ), ( x j ,y j ,z j ),node i With nodes j The coordinates are provided by a high-precision inertial navigation unit, then the signal arrival time difference Δ t ij satisfy: ; in, c It is the speed of light.

[0025] FDOA constraints: combined with nodes i With nodes j Speed ​​of movement ( v ix ,v iy ,v iz ), ( v jx , v jy , v jz (), provided by a high-precision inertial navigation unit, frequency offset difference Δ f ij satisfy: in, f 0 The signal center frequency.

[0026] The above system of equations is solved by iterative least squares method to obtain the initial position of the target in the local coordinate system. x loc ,y loc ,z loc The residual is calculated and used as the weight basis for subsequent global fusion. When the residual is ≤0.5m, the weight coefficient is ≥0.7.

[0027] Each local node uploads the preliminary positioning results and residuals to the network central node. For example, in a high-altitude reference UAV, the final positioning result is output through a global fusion algorithm. The weight allocation center node assigns weights to each initial position based on the residuals of each local result. w k : The weighted average method is used to fuse local results, and the formula is as follows: ; The multi-target signals are imported into a signal processing unit based on an FPGA hardware acceleration architecture. First, the signals are quickly sorted through hardware acceleration. Then, the frequency hopping / spread spectrum signal features are extracted to form real-time multi-target signal processing data that is coordinated with the positioning results, supporting subsequent source tracing and positioning. Based on the real-time processing data of the multi-target signals, the control source of the target UAV is identified by signal fingerprint feature matching (including radio frequency fingerprint + modulation ripple feature). At the same time, the multi-UAV angle of arrival (AOA) data of the three-dimensional stereo detection network is called up and the geospatial data of the geographic information system (GIS) is integrated to finally realize the source tracing and positioning of the control source, forming a complete technical closed loop of signal acquisition-target positioning-control source tracing.

[0028] The radio frequency fingerprint includes the signal's carrier frequency offset, phase noise, and power spectral density, while the modulation ripple features include the signal's modulation depth, inter-symbol interference mode, and symbol timing deviation. The target UAV's control source is identified by matching its "radio frequency fingerprint + modulation ripple features" with those of UAV control sources in a pre-defined fingerprint database.

[0029] The calibrated AOA data for each detection drone is combined with the global coordinates of the detection drone itself provided by the high-precision inertial navigation unit. X i ,Y i ,Z i ), generate a spatial direction line from the UAV to the control source, the expression is: ; in θ i For the first i The global azimuth angle of the drone. ϕ i X represents the global elevation angle, and X, Y, and Z represent the global coordinates of any point on the direction line.

[0030] The preliminary candidate location of the control source is calculated using the spatial straight line intersection method: every two AOA direction lines will form an intersection area in space, and the intersection of all intersection areas is calculated to obtain the preliminary candidate location area; Within the initial candidate region, for each candidate point ( X c ,Y c ,Z c ), calculate its distance to each drone. d c , Based on theoretical received powerP r,calc With the actual received power of the drone P r,meas Compare and retain | P r,calc - P r,meas Candidate points with a depth of less than 2dB are selected based on GIS data. Points that are physically inaccessible are removed from the candidate area, including inside buildings, on lakes, and on steep mountains. Candidate points in open areas such as plazas, rooftops, and roadsides are retained.

[0031] For the selected candidate points, the optimal positioning coordinates are calculated using the weighted least squares method. X final ,Y final ,Z final (Among them, since the control sources are mostly operated by ground personnel,) Z final Ground elevation is usually used.

[0032] Example 2 Based on Example 1, see Figure 2 The specific process of step S1 is as follows: S11: Deployment of the launch platform: Deploy the launch equipment equipped with a multi-band adaptive array antenna to a preset altitude of 50-300 meters, complete the platform attitude calibration through the ground control system to ensure that the array antenna is in a stable state with a horizontal reference plane error of ≤0.5°, and complete the initial spatial coordinate anchoring based on GNSS positioning data; S12: Array Antenna Activation and Frequency Band Configuration: Activate the multi-band adaptive array antenna, which contains 3 or more subarrays, such as subarray A covering 1.2-2.4GHz, subarray B covering 2.0-3.8GHz, and subarray C covering 3.4-6GHz. Through the frequency band overlap design between subarrays, the frequency band overlap between adjacent subarrays is ≥15%. For example, subarrays A and B overlap in 2.0-2.4GHz, and subarrays B and C overlap in 3.4-3.8GHz, achieving seamless coverage of the 1.2-6GHz frequency band. S13: Beamforming Synergy: Each subarray independently generates an adjustable beam through phased array technology. The ground control unit dynamically adjusts the beam pointing angle and gain distribution between subarrays based on real-time electromagnetic environment monitoring data, so that the beams of adjacent subarrays form complementary spatial coverage in overlapping frequency bands. A global spatial spectrum estimation model is constructed through coherent synthesis of signals between subarrays, realizing spatial scanning of electromagnetic signals within the range of 360° azimuth angle and -30° to +30° elevation angle. S14: Signal Acquisition and Multi-channel Output: When the target signal enters the scanning range, the array antenna identifies the direction of the incoming signal through spatial spectrum estimation, drives the corresponding subarray beam to focus on the direction of the incoming signal (beamwidth ≤ 3°) to enhance the receiving gain, and simultaneously activates the multi-channel AD sampling module to convert the original RF signals captured by each subarray into digital signals, and outputs them to the signal processing unit in real time according to the subarray channel number, such as channel 1 corresponding to subarray A and channel 2 corresponding to subarray B, forming a multi-channel raw dataset containing timestamps, signal strength and phase information.

[0033] The specific process of performing time-frequency domain joint analysis on the multi-channel original signal in step S2 to extract the frequency hopping period, modulation scheme, and pseudo-random code features is as follows: S20: Multi-channel signal preprocessing: Synchronous calibration of the original digital signals of the multi-channel output of each subarray, correction of time deviation based on timestamp, and filtering out environmental noise and co-frequency interference by wavelet threshold denoising to retain effective signal components; S21: Joint Time-Frequency Domain Analysis: Perform a short-time Fourier transform on the preprocessed signal to generate a time-frequency distribution map, identify the energy distribution characteristics of the signal in the time-frequency plane, and locate the frequency jump time and corresponding carrier frequency value of the frequency hopping signal. For the steady-state signal segment in the time-frequency distribution map, transient features such as the rising / falling edge slope of the signal during frequency hopping are extracted by continuous wavelet transform. The instantaneous frequency and amplitude are calculated by combining Hilbert transform to determine the frequency hopping period. S22: Modulation mode identification: For non-frequency hopping signal segments or single-frequency dwell segments of frequency hopping signals, extract the amplitude histogram, phase change rate, and spectral correlation features of the signal; The modulation scheme is identified by matching features using a support vector machine (SVM) classifier, including AM, FM, FSK, PSK, and QAM.

[0034] S23: Pseudo-random code feature extraction: For signals identified as spread spectrum modulation, the pseudo-random code period is determined by sliding correlation calculation; The code sequence features are extracted based on the signal phase reversal points within the code period, and compared with common UAV pseudo-random codes in a preset feature library to determine the code type and code length.

[0035] The specific process of matching the extracted features with a preset UAV signal feature library in step S2 to identify the target UAV signal is as follows: S24: Standardization and Feature Library Retrieval: The target signal features (frequency hopping period, modulation scheme, and pseudo-random code features) extracted in step 2 are standardized and converted into a numerical format consistent with the baseline parameters of the feature library. For example, the frequency hopping period is converted into a millisecond-level value, and the modulation scheme is mapped to an enumeration type code in the feature library. A preset UAV signal feature library is retrieved. This library contains a set of baseline features for known UAV models, including parameters such as frequency hopping period range, modulation scheme type, pseudo-random code type / code length, and carrier frequency fluctuation threshold. Each parameter is associated with a corresponding weight coefficient, where the weight of the modulation scheme and pseudo-random code features accounts for ≥60%. S25: Multi-dimensional feature matching calculation: For frequency hopping period: calculate the deviation rate between the target signal frequency hopping period and the reference period of each model in the feature library. Deviation rate = |target value - reference value| / reference value. When the deviation rate is ≤5%, it is judged as a preliminary match. Regarding modulation methods: a type matching method is used. If the target signal modulation method is exactly the same as a certain model in the feature library, the score is 100% of the weighted score of that parameter. If it is a compatible type, such as the target being 4FSK and belonging to the same FSK family as 2FSK in the feature library, the score is 60% of the weighted score. For pseudo-random code characteristics: the correlation coefficient is calculated by code sequence correlation. A correlation coefficient ≥ 0.85 is considered a high match. Combined with the code length deviation, the higher the correlation coefficient and the smaller the code length deviation, the higher the score.

[0036] S26: Comprehensive Matching Evaluation and Target Recognition: The multi-dimensional matching scores are weighted and summed to calculate the comprehensive matching degree (out of 100) between the target signal and each UAV model in the feature library; a matching threshold (≥70 points) is set. When the comprehensive matching degree of a certain model exceeds the threshold, the target signal is determined to be the signal of the corresponding UAV model; if multiple models have matching degrees exceeding the threshold, the model with the highest score is selected as the recognition result, and the second highest score model is marked as a candidate; if the matching degree of all models is below the threshold, the signal is classified as an unknown UAV signal, and the feature library update mechanism is triggered.

[0037] The specific process of separating the control command link and image transmission link signals of the target signal in step S2 is as follows: S27: Preset distinguishing feature library between control command link and image transmission link: Control command link usually has the characteristics of low data rate (≤1Mbps), short frame length (≤256 bytes), high real-time performance (transmission interval ≤100ms), and FSK / PSK modulation; Image transmission link usually has the characteristics of high data rate (≥2Mbps), long frame length (≥1024 bytes), periodic burst transmission (frame interval 10-50ms), and QAM / OFDM modulation, and the image transmission signal often contains video coding features; S28: Signal Framing and Parameter Extraction: By detecting preset synchronization headers or frame period characteristics, the identified target UAV signal is processed for frame synchronization, dividing the continuous signal into independent data frames; key parameters of each frame are extracted, including data rate, frame length, modulation method, and transmission interval (time difference between adjacent frames). S29: Link Type Classification and Separation: A decision tree classification model is used to determine the link type of the framed signal. If the signal frame meets the following conditions: data rate ≤ 1Mbps, frame length ≤ 256 bytes, transmission interval ≤ 100ms, and modulation method is FSK / PSK, then it is determined to be a control command link. If the signal frame meets the following conditions: data rate ≥ 2Mbps, frame length ≥ 1024 bytes, transmission interval 10-50ms, and modulation method is QAM / OFDM, then it is determined to be an image transmission link; Physical layer separation is performed on the two types of link signals. Independent signal buffers are allocated to store frame data of the control command link and the image transmission link respectively, while the timestamp association between the two types of links is recorded.

[0038] Link consistency verification: Verify whether the two types of link signals after separation come from the same UAV: ​​By comparing the carrier frequency drift trend (drift difference ≤ 5kHz), transmit power fluctuation correlation (correlation coefficient ≥ 0.7) and spatial direction of arrival angle deviation (deviation ≤ 1°) of the two, if all meet the consistency conditions, the separation result is confirmed to be valid; if not, the framing and classification process is re-executed until the consistency requirements are met.

[0039] Example 3 Based on Example 1 or Example 2, the specific process of establishing a dynamic spatial coordinate system in step S3 based on the real-time GNSS positioning data of the launch platform and the attitude parameters output by the inertial navigation system is as follows: S31: Reference coordinate anchoring: With the center of mass of the launch platform as the origin O of the coordinate system, the geodetic coordinates (L, B, H) of the platform are obtained in real time through the GNSS module and converted into absolute coordinates (X0, Y0, Z0) in the WGS-84 geocentric rectangular coordinate system, which serve as the spatial anchor point of the dynamic coordinate system. Set the initial orientation of the coordinate system: Starting from the origin O, the vector pointing from the Earth's center to point O is the Z-axis (the celestial axis, pointing in the direction of the Earth's outer normal). The geodetic coordinates are converted into planar coordinates through Gaussian projection, and the initial east-west axis (X-axis) and north-north axis (Y-axis) are determined, forming the initial ENU (East-North-Celestial) coordinate system framework.

[0040] S32: Attitude parameter fusion and dynamic coordinate system correction: Real-time acquisition of attitude parameters output by INS, including pitch angle (θ, rotation angle around the X-axis), roll angle (φ, rotation angle around the Y-axis), and yaw angle (ψ, rotation angle around the Z-axis), with a sampling frequency ≥100Hz; Constructing a rotation matrix based on attitude angles: According to the Euler angle transformation formula, the initial ENU coordinate system is rotated around the Z-axis, Y-axis, and X-axis in sequence to generate a dynamic rotation matrix R(ψ, θ, φ); The initial coordinate system is corrected in real time by using a rotation matrix: the axis vectors (X0, Y0, Z0) of the initial ENU coordinate system are multiplied by the rotation matrix R to obtain the corrected instantaneous coordinate axis directions (X', Y', Z'), ensuring that the coordinate system always keeps in line with the real-time attitude of the launch platform. For example, when the platform tilts, the coordinate system tilts synchronously with it.

[0041] S33: Time synchronization and coordinate update mechanism: The PPS (pulse per second) signal of GNSS is used as the time reference to synchronize the attitude parameters output by INS with the GNSS positioning data by timestamp, and the synchronization error is controlled within 1ms. A coordinate system update is performed every 50ms: combining the latest GNSS absolute coordinates (X0', Y0', Z0') and INS attitude angles (ψ', θ', φ'), the rotation matrix is ​​recalculated and the coordinate system origin position and axis orientation are corrected to compensate for platform drift; When the GNSS signal is briefly lost (≤3s), the coordinate system origin position is temporarily updated using INS track estimation data, and the position offset is predicted by forward differential to maintain the continuity of the coordinate system.

[0042] In S3, the time difference of arrival (TDoA) and phase difference of arrival (PDoA) of the target UAV signal between multiple subarrays are transformed into this coordinate system, and a three-dimensional positioning equation system is constructed to solve for the real-time spatial position coordinates of the target UAV. The specific process is as follows: S34: Multi-subarray spatial coordinate calibration: In the established dynamic spatial coordinate system (origin is the centroid O of the launch platform, and the axes are X', Y', Z'), the relative coordinates of each subarray antenna are pre-calibrated: let the coordinates of subarray i be (xᵢ, yᵢ, zᵢ), i=1, 2, ..., n, n≥3. The physical distance between subarrays is obtained by laser ranging and converted into three-dimensional coordinate values ​​in the coordinate system to form the subarray position matrix. Record the phase center offset of each subarray (≤0.05m) for error compensation in subsequent PDoA calculations.

[0043] S35: Spatiotemporal synchronization conversion between TDoA and PDoA: TDoA Conversion: For the target signal received by multiple subarrays, taking subarray 1 as the reference, calculate the arrival time difference Δtᵢ1 between the remaining subarrays i and subarray 1 (by extracting the signal leading edge time difference through a cross-correlation algorithm), and combine this with the electromagnetic wave propagation speed c (3×10⁻⁶). 8 (m / s), converted to distance difference Δdᵢ1=Δtᵢ1×c, yielding the relationship: [(X-xᵢ)²+(Y-yᵢ)²+(Z-zᵢ)²] 1 / 2 -[(X-x1)²+(Y-y1)²+(Z-z1)²] 1 / 2 =Δdᵢ1, where (X, Y, Z) are the coordinates of the target in the dynamic coordinate system to be solved; PDoA Conversion: For signals with the same carrier frequency f, calculate the phase difference Δφᵢ1 between subarray i and subarray 1, and convert it to the distance difference Δd'ᵢ1=(Δφᵢ1×λ) / (2π), where λ is the signal wavelength and λ=c / f. Considering the baseline length between subarrays bᵢ1 (≤λ / 2 to avoid phase ambiguity), we obtain the supplementary relation: [(X-xᵢ)²+(Y-yᵢ)²+(Z-zᵢ)²] 1 / 2 -[(X-x1)²+(Y-y1)²+(Z-z1)²] 1 / 2 =Δd'ᵢ1+kᵢ1×λ, where kᵢ1 is an integer, and the range of values ​​for kᵢ1 is constrained by the TDoA result.

[0044] S36: Construction and Solution of Three-Dimensional Positioning Equations Select three or more subarray combinations and construct a nonlinear equation system containing three unknowns (X, Y, Z) based on the distance difference relationship between TDoA and PDoA conversions mentioned above. The nonlinear equation system is linearized using the Taylor series expansion method: taking the initial position estimate of the target (a rough positioning result based on the signal direction angle) as the starting point of the iteration, the Jacobian matrix is ​​calculated, and the equation system is transformed into the linear form Ax=b; The linear equations are solved by the least squares method to obtain the target coordinate correction. The correction is iteratively updated until the correction is ≤0.1m. The converged (X, Y, Z) is output as the real-time spatial coordinates of the target in the dynamic coordinate system.

[0045] S37: Coordinate system transformation and absolute position output: Transform the target coordinates (X, Y, Z) in the dynamic coordinate system to the WGS-84 geocentric rectangular coordinate system: using the rotation matrix R in step 3. -1 The attitude is calculated using the inverse of the rotation matrix to eliminate the influence of the platform attitude on the coordinates and obtain the absolute offset of the target relative to the origin of the launch platform. By superimposing the GNSS absolute coordinates (X0, Y0, Z0) of the launch platform, the WGS-84 geodetic coordinates of the target UAV are obtained, completing the closed loop of the solution from signal time difference / phase difference to geographic coordinates.

[0046] The specific process of step S4 is as follows: S41: Extraction of characteristic parameters from the control command link transmitter: Carrier frequency drift rate extraction: Carrier frequency tracking is performed on the control command link signal using phase-locked loop technology with a tracking bandwidth ≤1kHz. Carrier frequency values ​​are continuously collected within 100ms with a sampling interval of 1ms. The carrier frequency drift rate over time is calculated through linear fitting, with the unit being Hz / s. Power fluctuation period extraction: The received signal power is sampled in real time at a sampling rate of 100Hz. The power spectral density is analyzed by Fourier transform to identify the frequency value corresponding to the main peak. Its reciprocal is the power fluctuation period. At the same time, the fluctuation amplitude range is recorded. Encoding format parsing: Decode the control command frame data, extract the frame header structure (synchronization code sequence), check method (such as CRC16 / 32), data frame length and bit synchronization method (such as Manchester encoding, differential encoding), and generate encoding feature code (used to distinguish different control sources of the same model).

[0047] S42: Construction of the control signal propagation path model: Basic path model: Using the real-time UAV position coordinates (X1, Y1, Z1) calculated in step S3 as the relay point, assuming the control signal travels directly from the control source (X0, Y0, Z0) to the UAV (line-of-sight propagation), and combining the relative position of the UAV and the launch platform (X... p Y p Z p Establish a two-stage propagation model: control source → UAV → launch platform; Attenuation correction: Based on the received power P of the control command signal r Combined with the free space decay formula P r =P t -20lg(4πd / λ), where Pt is the transmit power, d is the propagation distance, and λ is the signal wavelength. An environmental attenuation coefficient is introduced, and terrain occlusion information scanned by the lidar of the launch platform is used for real-time correction. The occlusion area coefficient increases by 5-20dB, and the propagation distance d from the control source to the UAV is calculated. 01 ; Feature constraint supplement: Carrier frequency drift rate and power fluctuation period are used as model constraints. The carrier frequency drift rate of the same control source changes by ≤10% within a 1km range. The power fluctuation period is strongly correlated with the vibration characteristics of the control source antenna. For example, the fluctuation period of a handheld remote control is mostly 0.5-2s. This is used to eliminate false paths.

[0048] S43: Control source direction reverse tracing: Angle of arrival calculation: Based on the PDoA measurement value of the control command signal from the multi-subarray antenna of the launch platform, combined with dynamic spatial coordinate system transformation, the angle of arrival, azimuth α1 and elevation β1 of the signal from the UAV to the launch platform are obtained; Reverse derivation: Based on the UAV's position (X1, Y1, Z1) and motion state (speed, heading), calculate the incident direction when the control signal reaches the UAV. Assuming the UAV is an omnidirectional receiver, the azimuth angle of the control source relative to the UAV is α0 = α1 + 180° ± Δα, where Δα is the UAV attitude correction amount, ≤ 2°, and the pitch angle β0 is calculated by the difference between the signal incident altitude angle and the UAV's altitude. Trajectory intersection localization: Using the UAV's position as the vertex, construct rays along the reverse directions α0 and β0, and combine this with the propagation distance d calculated in step 2. 01 The candidate coordinates (X0', Y0', Z0') of the control source at that moment are obtained.

[0049] S44: Multi-time data fusion and coordinate calculation: The coordinates of candidate control sources are continuously collected at 5-10 time intervals of 0.5-2 seconds, synchronized with the UAV positioning update cycle; Kalman filtering is used to fuse candidate coordinates at multiple time points: with the control source moving speed ≤10m / s as the state constraint, coordinate drift is dynamically corrected to obtain the filtered control source coordinates (X0, Y0, Z0). Elevation constraints are applied using ground digital maps, with Z0 typically being ±5m above ground level. The final output is the WGS-84 geodetic coordinates of the control source.

[0050] Example 4 Based on Example 1, Example 2, or Example 3, a consumer-grade drone, model M1 of a certain brand, is found flying illegally in a certain area. Its signal characteristics are preset as follows: Control command link: FSK modulation, carrier frequency 1.8GHz, pseudo-random code is M sequence (code length 32 bits), data rate 0.6Mbps, frame length 192 bytes; Image transmission link: QAM modulation, carrier frequency 2.7GHz, frequency hopping period 600μs, data rate 8Mbps, frame length 2048 bytes, including H.265 video encoding synchronization code.

[0051] Taking this UAV as an example, the specific implementation process of the UAV signal detection and precise positioning method of the present invention using the aerial platform is illustrated as follows: The system deploys three reconnaissance drones (numbered U1, U2, and U3) to form an airborne platform. Its core configuration includes: Multi-band array antenna: 3 subarrays: Subarray A: 1.2-2.4GHz; Subarray B: 2.0-3.8GHz; Subarray C: 3.4-6GHz, with 20% overlap between adjacent subarray bands; GNSS module: Supports WGS-84 coordinate system; INS module: Sampling frequency 100Hz; AD sampling module: sampling rate 250MSps, quantization bit depth 16bit.

[0052] Detection platform deployment height: 180m; The average elevation of the area (to be used later to control the elevation of the source): 52.3m; Electromagnetic environment: No strong interference sources, terrain is mainly flat and open, with some low buildings in some areas, and the LiDAR scan confirms that the obstruction coefficient correction is 5dB.

[0053] First, platform deployment and calibration were performed: three reconnaissance drones were deployed to an altitude of 180m, and attitude calibration was completed through the ground control system to ensure that the horizontal reference plane error of the array antenna was ≤0.3°. Initial anchoring coordinates were then obtained based on GNSS. U1: (L1=116.412358°E, B1=39.920145°N, H1=180m) U2: (L2=116.413892°E, B2=39.921567°N, H2=180m) U3: (L3=116.411203°E, B3=39.922018°N, H3=180m) Antenna activation and frequency band configuration: Activate the multi-band array antenna, and the three sub-arrays work together to achieve seamless coverage of 1.2-6GHz. Sub-arrays A and B overlap in the frequency band of 2.0-2.4GHz, and sub-arrays B and C overlap in the frequency band of 3.4-3.8GHz.

[0054] Beamforming coordination: The ground control unit monitors the electromagnetic environment and dynamically adjusts the beam pointing angle of the subarray. The antennas of the three UAVs work together to form a full-domain scanning range with an azimuth angle of 360° and an elevation angle of -30° to +30°. Coherent synthesis improves the accuracy of spatial spectrum estimation.

[0055] Signal Acquisition and Output: After the target UAV signal enters the scanning range, the array antenna identifies the direction of arrival through spatial spectrum estimation and drives the corresponding subarray beam to focus, with a beamwidth of 2.5°. The AD sampling module is activated, converting the acquired 1.8GHz and 2.7GHz signals into digital signals and outputting them to the signal processing unit by channel to form a multi-channel raw dataset containing timestamps, signal strength, and phase information.

[0056] Then, signal feature extraction and target recognition are performed: Signal preprocessing: Synchronous calibration of the multi-channel raw signals is performed, and the deviation is corrected based on the timestamp; wavelet threshold denoising is performed using a soft-hard threshold trade-off function, with the threshold set to T=0.02V for the low-frequency band (1.8GHz) and 1.3T=0.026V for the high-frequency band (2.7GHz) to filter out environmental noise.

[0057] Time-frequency domain analysis and feature extraction: Short-time Fourier transform: Hanning window is used, with window length = carrier frequency period × 8, to generate time-frequency distribution map, locate the frequency hopping time and carrier frequency value of 2.7GHz signal, extract transient features through continuous wavelet transform, and calculate the frequency hopping period as 598μs by combining Hilbert transform; Modulation mode identification: The amplitude histogram and phase change rate features of the 1.8GHz signal were extracted and identified as FSK modulation by an SVM classifier; the 2.7GHz signal was identified as QAM modulation. Pseudo-random code extraction: Sliding correlation calculation is performed on the 1.8GHz FSK signal. The correlation window length is equal to the code length range × 1.5. The pseudo-random code period is determined to be 1ms. The extracted code sequence features are compared with the preset library and determined to be an M sequence with a code length of 32 bits.

[0058] Target recognition: Standardization process: Convert the frequency hopping period (598μs), modulation scheme (FSK / QAM), and pseudo-random code features (M sequence 32bit) into a feature library compatible format, call the preset feature library, with modulation scheme weighting 40%, pseudo-random code weighting 30%, and frequency hopping period weighting 30%. Multi-dimensional matching: Frequency hopping cycle deviation rate = |598-600| / 600≈0.33%, ≤5%, 30 points; Modulation mode perfectly matched, 70 points; Pseudo-random code correlation coefficient 0.93 (≥0.85), 30 points; Comprehensive matching degree = 30+70+30=130 points (full score 150 points, far exceeding the threshold of 80 points), determined to be a certain brand of M1 drone.

[0059] Link separation: Preset feature library: Control command link (≤1Mbps, ≤256 bytes, FSK / PSK), Image transmission link (≥2Mbps, ≥1024 bytes, QAM / OFDM); Framing and parameter extraction: The target signal is framed. The 1.8GHz signal has a data rate of 0.6Mbps, a frame length of 192 bytes, and a transmission interval of 60ms; the 2.7GHz signal has a data rate of 8Mbps, a frame length of 2048 bytes, and a transmission interval of 35ms. H.265 synchronization code is detected. The decision tree model determined that 1.8GHz was the control command link and 2.7GHz was the image transmission link. The carrier frequency drift difference between the two was ≤2kHz, the power fluctuation correlation coefficient was 0.85, and the direction of arrival angle deviation was ≤0.5°, confirming that the separation was effective.

[0060] Then, perform precise positioning of the target drone: Establishment of a dynamic spatial coordinate system: Baseline anchoring: With the centroid of U1 as the origin O, its GNSS coordinates are converted to WGS-84 geocentric rectangular coordinates (X0=12938546.32m, Y0=4875623.15m, Z0=4145892.78m), and an initial ENU coordinate system (X East, Y North, Z Sky) is established. Attitude correction: Collect the INS attitude parameters of U1, pitch angle θ=0.4°, roll angle φ=0.2°, yaw angle ψ=35°, construct the rotation matrix R(ψ, θ, φ) in the order of Z→Y→X, and correct the initial axis vector to obtain the instantaneous coordinate axis (X', Y', Z'); Time synchronization and update: The INS and GNSS data are synchronized using the GNSS PPS signal, and the coordinate system is updated every 50ms to compensate for platform drift.

[0061] Subarray coordinate calibration: In the dynamic coordinate system, calibrate the relative coordinates of the three subarrays of the U1 array antenna: subarray A (0.6m, 0m, 0m), subarray B (0m, 0.6m, 0m), and subarray C (-0.6m, 0m, 0m), with a phase center offset ≤ 0.03m.

[0062] TDoA to PDoA conversion: TDoA Conversion: Using subarray A as the reference, calculate the arrival time difference Δt between subarrays B and C and A. B1 =12ns、Δt C1 =15ns, converted to distance difference Δd B1 =3.6m, Δd C1 =4.5m, thus obtaining the following relationship: [(X-0)²+(Y-0.6)²+(Z-0)²] 1 / 2 -[(X-0.6)²+(Y-0)²+(Z-0)²] 1 / 2 =3.6; [(X+0.6)²+(Y-0)²+(Z-0)²] 1 / 2 -[(X-0.6)²+(Y-0)²+(Z-0)²] 1 / 2 =4.5; PDoA Conversion: For a 1.8GHz signal with wavelength λ=0.167m, calculate the phase difference Δφ between subarrays B and C and subarray A. B1=2π×3.2、Δφ C1 =2π×4.1 (after unwinding), converted to distance difference Δd' B1 =0.534m、Δd' C1 =0.684m, combined with the TDoA constraint k value, k B1 =21,k C1 =26, thus obtaining the supplementary relation.

[0063] Solving the positioning equations: A nonlinear equation system is constructed by selecting the relationship between subarrays AB and AC. Using the initial position (X0=500m, Y0=300m, Z0=100m) roughly estimated by the signal direction angle as the starting point for iteration, the Jacobian matrix is ​​calculated to linearize the equation system, which is then solved iteratively using the least squares method. Convergence occurs when the coordinate correction is ≤0.08m, yielding the target coordinates in the dynamic coordinate system (X=523.4m, Y=315.6m, Z=98.2m).

[0064] Absolute coordinate output: The platform attitude effect is eliminated by using the inverse matrix R⁻¹ of the rotation matrix to obtain the absolute offset of the target relative to the origin U1; the GNSS geocentric rectangular coordinates of U1 are then superimposed to convert it to WGS-84 geodetic coordinates. Longitude: Lt=116.417892°E Latitude: Bt=39.923456°N Altitude: Zt = 112.5m After verification, the positioning error is ≤3m.

[0065] Source tracing control: Emitter Feature Extraction: Carrier frequency drift rate: The 1.8GHz control signal was tracked using a phase-locked loop with a bandwidth of 0.8kHz. The carrier frequency value was collected within 100ms with a sampling interval of 1ms. The drift rate was obtained by linear fitting and was 18Hz / s. Power fluctuation period: The received power was sampled at a sampling rate of 100Hz. Fourier transform analysis showed that the main peak frequency of the power spectral density was 0.9Hz, corresponding to a fluctuation period of 1.11s and a fluctuation amplitude of ±1.8dB. Encoding format: Decode the control frame, extract the synchronization code sequence, CRC16 check, frame length 192 bytes, Manchester encoding, and generate the encoding feature code.

[0066] Propagation path model construction: based on the target UAV coordinates (L... t B t Z t Using [a specific point] as the relay point, a two-segment model is established: control source → UAV → U1. Based on the received power Pr = -58dBm, combined with the free space attenuation formula, and introducing a 5dB environmental attenuation coefficient corrected by the lidar, the propagation distance d from the control source to the UAV is calculated.01 =950m.

[0067] Control source direction tracing: Angle of arrival calculation: Based on PDoA measurements, the azimuth angle α1 = 75° and the pitch angle β1 = 6° of the signal from the UAV to U1 are obtained; Reverse derivation: The azimuth angle of the control source relative to the UAV is α0 = 75° + 180° ± 1° = 254° ~ 256°, take 255°, and the pitch angle is β0 = 3° (calculated by altitude difference). Candidate coordinate calculation: Using the UAV's position as the vertex, construct rays along α0 and β0, and combine them with d 01 =950m, and the candidate coordinates were calculated as (Lc1=116.408921°E, Bc1=39.918765°N, Hc1=89.3m).

[0068] Multi-time data fusion and coordinate output: Candidate coordinates were collected at 1-second intervals for 8 time points, fused using Kalman filtering, and constrained to a control source movement speed of ≤10m / s to dynamically correct drift. Combined with a ground digital map (the area's ground elevation is 52.3m), the elevation was corrected by a ground height +2m constraint, and the final output was the control source's WGS-84 geodetic coordinates. Longitude: Lc=116.409123°E Latitude: Bc = 39.918978°N Altitude: Hc = 54.3m The tracing error is ≤7m, enabling complete positioning of the target UAV and its control source.

[0069] The final result is: The target UAV WGS-84's geodetic coordinates are: 116.417892°E, 39.923456°N, 112.5m; Control source WGS-84 geodetic coordinates: 116.409123°E, 39.918978°N, 54.3m.

Claims

1. An aerial platform UAV signal detection and accurate positioning method, characterized in that, The method comprises the following steps: S1: deploying multiple unmanned aerial platforms carrying detection devices, the detection devices comprising multiple subarray multi-band array antennas, capturing global electromagnetic signals in real time and outputting multi-channel raw signal data; S2: performing joint time-frequency domain analysis on the multi-channel raw signal, extracting the frequency hopping period, modulation mode and pseudo-random code features of the signal, matching with the preset unmanned aerial vehicle signal feature library to identify the target unmanned aerial vehicle signal, and separating the control instruction link and the image transmission link signal of the target signal; S3: establishing a dynamic space coordinate system, converting the time difference of arrival TDoA and the phase difference of arrival PDoA of the identified target unmanned aerial vehicle signal between multiple subarrays to the coordinate system, constructing a three-dimensional positioning equation set, and solving to obtain the real-time spatial position coordinates of the target unmanned aerial vehicle; S4: for the control instruction link signal, extracting the transmission end feature parameters, combining the real-time spatial position coordinates of the target unmanned aerial vehicle, establishing a control signal propagation path model, and positioning the geographical coordinates of the unmanned aerial vehicle control source by reversely tracking the signal transmission direction and attenuation law.

2. The method of claim 1, wherein the method further comprises: The specific process of step S1 is as follows: S11: platform deployment: the unmanned aerial platform carrying the multi-band adaptive array antenna is launched to the preset height, the platform attitude is calibrated through the ground control system, and the initial spatial coordinate anchoring is completed based on the GNSS positioning data; S12: array antenna start and frequency band configuration: the multi-band adaptive array antenna is activated, which contains 3 or more subarrays, and the seamless coverage of 1.2-6GHz frequency band is realized through the frequency band overlap design between subarrays; S13: beamforming cooperation: each subarray independently generates a controllable beam through phased array technology, and the ground control unit dynamically adjusts the beam pointing angle and gain distribution between subarrays based on real-time electromagnetic environment monitoring data; S14: signal capture and multi-channel output: the array antenna identifies the signal coming direction through spatial spectrum estimation, drives the corresponding subarray beam to focus on the signal direction to enhance the receiving gain, and outputs to the signal processing unit in real time according to the subarray channel number, forming a multi-channel raw data set containing time stamp, signal strength and phase information.

3. The method of claim 1, wherein the method further comprises: The specific process of step S2 is as follows: S20: multi-channel signal preprocessing: the multi-channel raw digital signal output by each subarray is synchronized and calibrated, and the environmental noise and same frequency interference are filtered out through wavelet threshold denoising; S21: joint time-frequency domain analysis: Performing short-time Fourier transform on the preprocessed signal to generate a time-frequency distribution map, identifying the energy distribution characteristics of the signal on the time-frequency plane, and locating the frequency hopping time and corresponding carrier frequency value of the frequency hopping signal; For the steady-state signal segment in the time-frequency distribution map, the transient feature is extracted through continuous wavelet transform, the instantaneous frequency and amplitude are calculated combined with Hilbert transform, and the frequency hopping period is determined; S22: modulation mode identification: For non-frequency hopping signal segments or single frequency resident segments of frequency hopping signals, the amplitude histogram, phase change rate and spectrum correlation features of the signal are extracted; The features are matched through a support vector machine classifier to identify the modulation mode; S23: Pseudo-random code feature extraction: For the signals identified as spread spectrum modulation, the pseudo-random code period is determined by sliding correlation calculation; Based on the signal phase flip point within the code period, the code sequence features are extracted and compared with the common drone pseudo-random code in the preset feature library to determine the code type and length.

4. The method of claim 3, wherein the method further comprises: The specific process of matching the extracted features with the preset drone signal feature library in step S2 to identify the target drone signal is as follows: S24: Standardization processing and feature library calling: standardize the extracted target signal features; call the preset drone signal feature library, which contains a reference feature set of known drone models, and each parameter is associated with a corresponding weight coefficient; S25: Multi-dimensional feature matching calculation: For frequency hopping period: calculate the deviation rate of the target signal frequency hopping period and the reference period of each type in the feature library, deviation rate = |target value-reference value| / reference value, when the deviation rate ≤5%, it is determined as preliminary matching; For modulation mode: use type matching method, if the target signal modulation mode is completely consistent with a certain type in the feature library, the score is 100% of the corresponding weight value, if it is a compatible type, the score is 60% of the corresponding weight value; S26: Matching degree comprehensive evaluation and target recognition: weighted sum of multi-dimensional matching scores, calculate the comprehensive matching degree of target signal and each drone model in the feature library, and match the drone model according to the comprehensive matching degree.

5. The method of claim 4, wherein the method further comprises: The specific process of separating the control instruction link and the image transmission link signal of the target signal in step S2 is as follows: S27: Pre-set control instruction link and image transmission link distinguishing feature library; S28: Signal framing and parameter extraction: frame synchronization processing is performed on the identified target drone signal, and the continuous signal is segmented into independent data frames; the key parameters of each frame signal are extracted; S29: Link type classification and separation: use decision tree classification model to determine the link type of the framed signal.

6. The method of claim 1, wherein the method further comprises: The specific process of establishing a dynamic space coordinate system based on the GNSS real-time positioning data and the attitude parameters output by the inertial navigation system in step S3 is as follows: S31: Reference coordinate anchoring: taking the centroid of the airborne platform as the origin O of the coordinate system, the platform's geodetic coordinates are obtained in real time through the GNSS module, and converted into absolute coordinates (X0, Y0, Z0) in the WGS-84 geocentric rectangular coordinate system, which are used as the spatial anchor point of the dynamic coordinate system; Set the initial pointing of the coordinate system: taking the origin O as the starting point, the vector from the center of the earth to the O point is the Z axis, the geodetic coordinates are converted into plane coordinates through Gauss projection, the initial east axis X and north axis Y are determined, and the initial ENU coordinate system framework is formed; S32: Attitude parameter fusion and coordinate system dynamic correction: Real-time collection of attitude parameters output by the inertial navigation system INS, including pitch angle θ, roll angle φ, and heading angle ψ; rotation matrix construction based on attitude angle: according to the Euler angle transformation formula, the initial ENU coordinate system is rotated around Z axis, Y axis, and X axis in turn to generate a dynamic rotation matrix; Real-time correction of the initial ENU coordinate system by rotation matrix: multiply the axis vector (X0, Y0, Z0) of the initial ENU coordinate system by the rotation matrix to obtain the corrected instantaneous coordinate axis direction (X', Y', Z'); S33: Time synchronization and coordinate updating mechanism: Use the PPS signal of GNSS as the time reference to synchronize the time stamps of the attitude parameters output by INS and the GNSS positioning data; Update the coordinate system every specified time: combine the latest GNSS absolute coordinates (X0', Y0', Z0') and INS attitude angles (ψ', θ', φ') to recalculate the rotation matrix and correct the coordinate system origin position and axis direction, compensating for platform drift.

7. The method of claim 6, wherein the method further comprises: The specific process of converting the time difference of arrival and phase difference of arrival of the target UAV signal among multiple sub-arrays to the coordinate system in S3 to construct a three-dimensional positioning equation set and solve for the real-time spatial position coordinates of the target UAV is as follows: S34: Multi-subarray spatial coordinate calibration: In the established dynamic spatial coordinate system, the relative coordinates of each subarray antenna are calibrated in advance: let the coordinates of subarray i be (xᵢ, yᵢ, zᵢ), i = 1, 2,..., n, n ≥ 3, the physical distance between subarrays is obtained by laser ranging and converted to three-dimensional coordinate values in the coordinate system to form a subarray position matrix; S35: Spatio-temporal synchronization conversion of TDoA and PDoA: TDoA conversion: the target signal received by multiple sub-arrays is converted into a distance difference Δdᵢ1=Δtᵢ1×c, combined with the electromagnetic wave propagation speed c, to obtain the relationship: [(X-xᵢ)²+(Y-yᵢ)²+(Z-zᵢ)²] 1 / 2 -[(X-x1)²+(Y-y1)²+(Z-z1)²] 1 / 2 =Δdᵢ1, where (X, Y, Z) is the coordinate to be solved of the target in the dynamic coordinate system; PDoA conversion: for the signal of the same carrier frequency f, the phase difference Δφᵢ1 of the subarray i and the subarray 1 is calculated, which is converted into the distance difference Δd'ᵢ1=(Δφᵢ1×λ) / (2π), λ is the signal wavelength, λ=c / f, based on the baseline length bᵢ1 between the subarrays, a supplementary relationship is obtained: [(X-xᵢ)²+(Y-yᵢ)²+(Z-zᵢ)²] 1 / 2 -[(X-x1)²+(Y-y1)²+(Z-z1)²] 1 / 2 =Δd'ᵢ1+kᵢ1×λ, kᵢ1 is an integer; S36: Construction and solution of three-dimensional positioning equation set: Select three or more subarrays, and based on the distance difference relationship formula converted by TDoA and PDoA, construct a nonlinear equation set containing three unknowns (X, Y, Z); Linearize the nonlinear equation set using Taylor series expansion: take the initial position estimate of the target as the iteration starting point, calculate the Jacobian matrix, and convert the equation set to linear form Ax = b; Solve the linear equation set by least squares method to obtain the target coordinate correction, and iteratively update until the correction is ≤0.1m, and output the converged (X, Y, Z) as the real-time spatial coordinates of the target in the dynamic coordinate system; S37: Coordinate system conversion and absolute position output: Convert the target coordinates (X, Y, Z) in the dynamic coordinate system to the WGS-84 geocentric rectangular coordinate system: through the rotation matrix R -1 Perform attitude solution to eliminate the influence of platform attitude on coordinates and obtain the absolute offset of the target relative to the origin of the airborne platform. Superimpose the GNSS absolute coordinates (X0, Y0, Z0) of the ascending platform to obtain the WGS-84 geodetic coordinates of the target UAV.

8. The method of claim 7, wherein the method further comprises: The specific process of step S4 is as follows: S41: Feature parameter extraction of control instruction link transmission end: Carrier frequency drift rate extraction: track the carrier frequency of the control instruction link signal, continuously collect the carrier frequency within a specified time, and calculate the drift rate of the carrier frequency over time by linear fitting; Power fluctuation period extraction: real-time sample the signal receiving power, analyze the power spectral density by Fourier transform, identify the frequency value corresponding to the main peak, and its reciprocal is the power fluctuation period, and record the fluctuation amplitude range; Encoding format analysis: decode the control instruction frame data, extract the frame header structure, check method, data frame length, and bit synchronization method, and generate an encoding feature code; S42: Control signal propagation path model construction: The basic path model: taking the calculated real-time position coordinates (X1, Y1, Z1) of the unmanned aerial vehicle as a relay point, assuming that the control signal is directly from the control source (X0, Y0, Z0) to the unmanned aerial vehicle, combining the relative position (X p , Y p , Z p ) of the unmanned aerial vehicle and the take-off platform, a two-section propagation model is established: control source→unmanned aerial vehicle→take-off platform; Attenuation law correction: based on the received power P r of the control instruction signal, combined with the free space attenuation formula P r =P t -20lg(4πd / λ), P t is the transmit power, d is the propagation distance, λ is the signal wavelength, the environmental attenuation coefficient is introduced, and the terrain shielding information scanned by the airborne platform laser radar is used for real-time correction. The shielding area coefficient increases by 5-20dB, and the propagation distance d 01 from the control source to the unmanned aerial vehicle is inversely calculated. S43: Control source direction reverse tracking: Arrival angle calculation: based on the PDoA measurement value of the control instruction signal by the multi-subarray antenna of the airborne platform, combined with dynamic spatial coordinate conversion, the arrival angle of the signal from the UAV to the airborne platform is obtained, azimuth angle α1, pitch angle β1; Reverse direction derivation: according to the position (X1, Y1, Z1) and motion state of the UAV, the incident direction of the control signal when it reaches the UAV is calculated, assuming that the UAV is omnidirectional reception, the azimuth angle of the control source relative to the UAV α0=α1+180°±Δα, Δα is the attitude correction amount of the UAV; the pitch angle β0 is calculated by the signal incident elevation angle and the altitude difference of the UAV; Trajectory intersection positioning: taking the position of the UAV as the vertex, constructing a ray along the reverse direction α0, β0, combining the back-propagated distance d 01 , to obtain the candidate coordinates (X0', Y0', Z0') of the control source at this moment; S44: multi-time data fusion and coordinate solution: Collect 5-10 time control source candidate coordinates synchronously with UAV positioning update period; Use Kalman filter to fuse multi-time candidate coordinates: take the control source moving speed ≤10m / s as the state constraint, dynamically correct the coordinate drift, and get the filtered control source coordinates (X0, Y0, Z0); Combined with the altitude constraint of the digital map, Z0 is the ground height ±5m, and finally output the WGS-84 geodetic coordinates of the control source.

9. The method of claim 1, wherein the method further comprises: Also includes the following process: Deploy multiple UAV platforms carrying detection equipment, complete dynamic cooperative networking through the preset frequency domain / time domain / code domain joint communication protocol combined with dynamic spectrum sensing technology, and build a three-dimensional detection network that can avoid building and mountain terrain shielding; Through the three-dimensional detection network, the target UAV's various signals are captured to obtain multi-target signals; simultaneously collect the time difference arrival data, frequency difference arrival data of each UAV platform, and the real-time positioning data generated by the inertial navigation unit of each platform to form a multi-source data set; After preprocessing the multi-source data set, the target UAV is located by distributed positioning solution; Import multi-target signals into signal processing units based on FPGA hardware acceleration architecture, first complete signal rapid sorting through hardware acceleration, then extract frequency hopping / spreading signal features, form multi-target signal real-time processing data coordinated with positioning results, and support subsequent traceability positioning; Based on multi-target signal real-time processing data, identify the control source of the target UAV through signal fingerprint feature matching, simultaneously call the multi-UAV arrival angle data of the three-dimensional detection network, and fuse geographic spatial data of the geographic information system, finally realize the traceability positioning of the control source.

Images