A method for extracting micro-motion parameters of a UAV at a low pulse repetition frequency

Abstract

The present application provides a kind of unmanned aerial vehicle micro-motion parameter extraction method under low pulse repetition frequency, comprising the following steps: calculating time-frequency spectrum entropy;Extracting rotational frequency;Blade length extraction.A kind of method for verifying the effectiveness of the above-mentioned unmanned aerial vehicle micro-motion parameter extraction method under low pulse repetition frequency is also provided.The present application can effectively estimate blade rotational frequency and blade length under low pulse repetition frequency conditions;Compared with VMD-OMP method, the relative error of blade rotational frequency and blade length estimation is smaller;Compared with traditional OMP algorithm, the relative error of blade rotational frequency and blade length estimation is smaller, and the operation time is shorter.The present application can realize micro-motion target feature parameter extraction and identification under low pulse repetition frequency, and improve the overall working performance of radar.

Description

TECHNICAL FIELD

[0001] The present application relates to signal and information processing technology, and particularly relates to a multi-rotor unmanned aerial vehicle micro-motion parameter extraction method under a low pulse repetition frequency. BACKGROUND

[0002] Unmanned aerial vehicles have many advantages such as low cost, vertical take-off and landing, and hovering at any time, and have been widely used in military and civilian fields. In the scene of detecting unmanned aerial vehicles by radar, in addition to the translational flight, the blades of the unmanned aerial vehicle also rotate around the rotation shaft during flight. These fine movements will periodically modulate the backscattering echo of the target, that is, sidebands are generated on both sides of the Doppler spectrum of the main body of the unmanned aerial vehicle. This additional Doppler frequency modulation is called micro-Doppler effect. The micro-Doppler effect provides effective and rich information for the motion state monitoring, imaging and classification recognition of the unmanned aerial vehicle. In recent years, the analysis and extraction of the micro-motion characteristics of the unmanned aerial vehicle have become one of the research hotspots at home and abroad, and have been continuously and widely concerned in the field of radar automatic target recognition (ATR).

[0003] Based on the radar incidence and reflection mechanism, the rotating blades of unmanned aerial vehicles (UAVs) or helicopters will produce multiple sinusoidal frequency modulation signal components with a 'flickering effect'. When the radar has a high pulse repetition frequency (HPRF), it can obtain a time-frequency spectrogram containing multiple sinusoidal frequency modulation signal components and blade flickering characteristics. The flickering effect period of the UAV echo signal and the instantaneous frequency of the sinusoidal frequency modulation signal are directly related to the radius and rotation frequency of the UAV blades. In recent years, methods for extracting the micro-motion characteristics of UAVs based on the flickering effect and micro-motion instantaneous frequency have been proposed. The most commonly used method is to extract the micro-motion parameters of UAVs from the time-frequency spectrogram of the echo. Professor Li Gang's team at Tsinghua University converts the time-frequency spectrogram of the UAV echo signal into a rhythm velocity graph and estimates the flickering period of the UAV blades by extracting the frequency of the rhythm velocity graph. Professor Wan Xianrong's team at Wuhan University uses an external radar source to estimate the micro-motion parameters of the helicopter rotor by using time-frequency analysis and the orthogonal matching pursuit (OMP) algorithm. Harmanny et al. proposed using robust singular value decomposition (SVD) to decompose the left and right singular vectors from the time-frequency spectrogram of the UAV echo and extract the flickering period and spectral width. For sinusoidal frequency modulation signal feature extraction, Hou et al. used fractional Fourier transform to remove the echo modulation caused by the acceleration and speed of the UAV, and used fast synchronous squeezing transform to estimate the instantaneous frequency of the sinusoidal frequency modulation signal. Song Chen et al. from the Innovation Institute of Aerospace Science and Technology of the Chinese Academy of Sciences projected the UAV signal from the time-frequency domain to the time-frequency rotation domain and used the particle swarm algorithm to extract the time-frequency concentration index of the time-frequency rotation domain to estimate the instantaneous frequency of the sinusoidal frequency modulation signal. Qin et al. used iradon to extract the instantaneous frequency information of the sinusoidal frequency modulation from the UAV signal detected by the terahertz radar. However, when the maximum Doppler shift caused by the rotation of the UAV blades exceeds half of the pulse repetition frequency (PRF), the micro-Doppler signal will exhibit spectral aliasing, wrapping, and other phenomena in the time-frequency domain, making it difficult to effectively extract the micro-Doppler features using the above methods.

[0004] When the radar has a low pulse repetition frequency (LPRF), i.e. the PRF is less than four times the blade rotation speed, the data with length of several rotation periods can be used to observe the jet engine modulation (JEM) and helicopter rotation modulation (HERM) features of the unmanned aerial vehicle (UAV) from the frequency domain and time-frequency domain, respectively. By extracting the spectral peak / spectral line interval of the JEM / HERM, the corresponding micro-motion parameters can be estimated. Huang et al. found the fundamental frequency of the HERM feature by using a longer short time fourier transform (STFT) window, and estimated the spectral line interval of the HERM by using a logarithmic harmonic summation algorithm. Fioranelli et al. classified UAVs with different loads by extracting the center of gravity and the spectral width of the first harmonic of the HERM feature. Peter et al. converted the HERM feature extraction into a pitch detection problem, and proposed a multi-frequency detector to obtain the spectral line interval of the HERM. However, the above methods can only extract the rotation frequency of the UAV, and cannot estimate the radius of the UAV blade. SUMMARY

[0005] In view of the deficiencies in the prior art, the present application provides a method for extracting micro-motion parameters of an unmanned aerial vehicle (UAV) under low pulse repetition frequency (LPRF), which specifically comprises the following steps:

[0006] Step 1: Calculate the time-frequency spectrum entropy

[0007] In order to select a signal with obvious HERM feature, the UAV echo signal is first divided into M segments with a time length of t m , i.e. s(t) = [s d1 (t m ),..., s dm (t m ),..., s dM (t m )] T , where s(t) is the UAV echo signal, s dm (t m ) is the mth segment signal, and m ∈ [1, M]; the short time fourier transform (STFT) is performed on the segmented signal s dm (t m ) to obtain its corresponding time-frequency representation STFT(t i , f j ), where i represents the horizontal coordinate and j represents the vertical coordinate, and the lengths of i and j are both t m ;

[0008] The time-frequency spectrum entropy is used to describe the energy distribution difference of different segmented signals in the time-frequency distribution

[0009]

[0010] H m represents the time-frequency spectrum entropy of the mth segmented signal, and the time-frequency plane is equally divided into i x j time-frequency blocks with equal area; the energy X ij of the (i, j)th time-frequency block corresponds to the amplitude of the time-frequency representation |STFT(t i , f j )| ; P ij represents the probability of the energy of the (i, j)th time-frequency block relative to the energy of the entire time-frequency plane;

[0011] Second step: Extracting the rotation frequency

[0012] The VMD is used to extract the center frequency point of each spectral peak, and the spectral line interval of the micro-Doppler frequency of the unmanned aerial vehicle blade is calculated; the decomposition problem of the multi-peak signal is constructed as:

[0013]

[0014] In the formula: is the partial derivative with respect to time t; s k is the time-domain signal of the kth spectral peak; ω k is the center frequency of the kth spectral peak, which is a variable with respect to time t; K is the number of peaks; is the segmented signal with the maximum time-frequency spectrum entropy;

[0015] The center frequency ω k of the multi-peak signal is solved using the alternating multiplier method, and according to the relationship between the HERM spectral line interval and the rotation frequency, the estimated value of the rotation frequency of the unmanned aerial vehicle rotor blade of N blades is represented as:

[0016]

[0017] In the formula, mean(·) represents the mean value operation, and diff(·) represents the first-order difference operation, is the estimated value of ω k , and N is the number of blades;

[0018] Third step: Blade length extraction

[0019] On the basis of extracting the estimated value of the rotation frequency , the OMP method of atomic scaling is proposed to estimate the blade length; the dictionary established based on the unmanned aerial vehicle blade length estimation where the ith atom under the matching condition of the xth scaling is:

[0020]

[0021] where η x is the scaling factor of the xth time, η1=1, i is the serial number of the atom in the dictionary; α, β, R0are the azimuth angle, the pitch angle and the distance of the blade rotation center relative to the radar at the initial time, respectively; is the blade rotation frequency; f c is the radar transmitting carrier frequency; c is the speed of light; is the xth blade length estimate value; is the phase of the nth blade at time t; λ is the wavelength; all the above are prior information; the initial phase of the blade rotation is ignored to reduce the dimension of the dictionary;

[0022] The OMP algorithm is improved as follows:

[0023] (1) Signal decomposition: when the carrier frequency is f c and the scaling factor is η1, the sparse decomposition is performed on the segmented signal with the maximum entropy , and the atom obtained after g iterations is energy-normalized as The energy estimate maximum value of the ith atom under the first scaling matching condition is selected , and the corresponding length unit position is At this time , the relationship between it and the blade length estimate value is:

[0024]

[0025] In the formula, l0is the initial value of the iteration search; L is the total length of the search; g is the counter for the iteration count of i; is the energy estimate maximum value corresponding to the g value; G is the maximum iteration number;

[0026] (2) Scaling matching: under the scaling factor η x condition, the length unit position is scaled in size, and the theoretical position of the length unit after scaling matching is , which is expressed as:

[0027]

[0028] In the formula, x is the scaling number, which is a positive integer greater than 1;

[0029] At this time, , the energy value corresponding to it is expressed as:

[0030]

[0031] wherein, δ is a threshold value; is a scaling coefficient η x the energy estimation maximum value of the i-th atom under the condition, is a scaling coefficient η x the theoretical position of the length unit after scaling matching under the condition;

[0032] (3) Energy accumulation; the energy accumulation maximum value is selected by accumulating the energy of the corresponding position and and The step is represented as:

[0033]

[0034] The estimation value obtained from the formula (6), (9) as the final estimation result of the blade length

[0035] In one specific embodiment of the present application, the scaling coefficients η1=1, η2=0.75, η3=0.5, η4=0.25, and the threshold δ=10.

[0036] In another specific embodiment of the present application, the initial value l0=0.2m, the equivalent length L=2.8m, and the maximum iteration number G=2000.

[0037] An effectiveness verification method of the unmanned aerial vehicle micro-motion parameter extraction method under the low pulse repetition frequency is also provided, and specifically includes the following steps:

[0038] First step: data preprocessing

[0039] ①Determination of state; the motion state of the unmanned aerial vehicle is determined by analyzing the time and distance information collected by the radar. When the unmanned aerial vehicle is in a hovering state, the distance does not change with time. When the unmanned aerial vehicle is in a uniform motion state, the distance and time are in a linear relationship. When the unmanned aerial vehicle is in an accelerated motion state, the distance and time are in a nonlinear relationship.

[0040] ②Signal synthesis; the echo data collected by the radar is synthesized into a complex signal by formula (10):

[0041] s(t) = I(t) + jQ(t) (9)

[0042] wherein, I(t) + jQ(t) is simply written as I and Q, and the I / Q signal is a same-phase quadrature signal; the phase of Q is 90° different from that of I.

[0043] ​③ Slice selection; take 1000 (i.e. 1s into 1000ms) as a slice unit, slice segment the entire data segment according to the total sampling time R, divide the entire data segment into T / 1000 time segments, calculate the time-frequency spectrum entropy of all segmented signals according to formula (2), select the data slice with the maximum time-frequency spectrum entropy for further processing;

[0044]

[0045] In the formula: H m represents the time-frequency spectrum entropy of the mth signal, and the time-frequency plane is equally divided into i x j time-frequency blocks with equal area; the energy X ij of the (i, j)th time-frequency block corresponds to the amplitude |STFT(t i , f j ); ij P k represents the probability of the energy of the (i, j)th time-frequency block relative to the energy of the entire time-frequency plane;

[0046] Second step: extract the rotation frequency

[0047] The VMD algorithm is used to process the data with the maximum time-frequency spectrum entropy, the spectral interval of the HERM line of the unmanned aerial vehicle echo signal is estimated by formula (3), and the rotation frequency of the unmanned aerial vehicle blade is calculated by formula (4);

[0048] and

[0049] In the formula: is the segmented signal with the maximum time-frequency spectrum entropy; is the partial derivative with respect to time t; s k is the time-domain signal of the kth spectral peak; ω k is the center frequency of the kth spectral peak, which is a variable with respect to time t; K is the number of peaks;

[0050] According to the relationship between the HERM spectral line interval and the rotation frequency, the estimated value of the rotation frequency of the unmanned aerial vehicle rotor blade of N blades can be represented as:

[0051]

[0052] In the formula, mean(·) represents the mean value operation, diff(·) represents the first-order difference operation, is the estimated value of ω k , and N is the number of blades;

[0053] Third step: extract the blade length

[0054] On the basis of the acquisition of the blade rotation frequency, the atomic dilatation OMP method is used to concentrate the micro-Doppler feature energy of the echo data according to the formula 5 to formula 9, so as to estimate the blade length;

[0055] Dictionary established based on the unmanned aerial vehicle blade length estimation Wherein the i-th atom under the x-th dilatation matching condition is:

[0056]

[0057] Wherein, η x is the dilatation coefficient of the x-th time, η1=1, i is the serial number of the atom in the dictionary; α, β, R0 are the azimuth angle, the pitch angle and the distance of the blade rotation center relative to the radar at the initial moment respectively; is the blade rotation frequency; f c is the radar transmitting carrier frequency; c is the speed of light; is the x-th blade length estimation value; is the phase of the n-th blade at t moment; λ is the wavelength; the above are all prior information; in order to reduce the dimension of the dictionary, the initial phase of the blade rotation is ignored;

[0058] When the carrier frequency is f c and the dilatation coefficient is η1, the sparse decomposition is carried out on the segmented signal with the maximum entropy After g times of iterative operation, the atom is obtained, and the energy normalization is The energy estimation maximum value of the i-th atom under the first dilatation matching condition is selected The corresponding length unit position At this time The relationship between the length unit position and the blade length estimation value is:

[0059]

[0060] In the formula, l0 is the initial value of the iteration search; L is the total length of the search; g is a counter for the iteration count of i; is the energy estimation maximum value corresponding to the g value; G is the maximum iteration number;

[0061] The length unit position is dilated in scale under the condition of the dilatation coefficient η x The theoretical position of the length unit after dilatation matching is expressed as:

[0062]

[0063] In the formula, x is the dilatation number sequence, which is a positive integer greater than 1.

[0064] At this time, The corresponding energy value is expressed as:

[0065]

[0066] In the formula, δ is a threshold value; is a scaling coefficient η x The energy estimation maximum value of the i-th atom under the condition, is a scaling coefficient η x The length unit theoretical position after scaling and matching under the condition;

[0067] The above obtained and The energy of the corresponding position and are accumulated, and the maximum energy accumulation value is selected, which is expressed as:

[0068]

[0069] The estimation value obtained from the formula (6), (9) is the final estimation result of the blade length

[0070] In one specific embodiment of the present application, the scaling coefficients η1=1, η2=0.75, η3=0.5, η4=0.25, and the threshold δ=10.

[0071] In another specific embodiment of the present application, the initial value l0=0.2m, the equivalent length L=2.8m, and the maximum iteration number G=2000.

[0072] The present application has the following advantages:

[0073] 1. Under the condition of low pulse repetition frequency, the blade rotation frequency and the blade length can be effectively estimated;

[0074] 2. Compared with the VMD-OMP method, the relative error of the estimation of the blade rotation frequency and the blade length is smaller;

[0075] 3. Compared with the traditional OMP algorithm, the relative error of the estimation of the blade rotation frequency and the blade length is smaller, and the operation time is shorter. BRIEF DESCRIPTION OF DRAWINGS

[0076] Figure 1 A schematic flow chart of a rotor unmanned aerial vehicle micro-motion parameter extraction method under the condition of LPRF according to the present application is shown;

[0077] Figure 2 A rotating schematic diagram of an unmanned aerial vehicle blade is shown;

[0078] Figure 3 Results of time-frequency spectrum entropy estimation are shown, wherein Figure 3 (a) time-domain echo signal is shown, Figure 3 (b) normalized time-frequency spectrum entropy is shown, Figure 3 (c) time-frequency representation of paragraph 20 is shown;

[0079] Figure 4 Results of variational modal decomposition (VMD) processing are shown, wherein Figure 4 (a) JEM spectral peak feature is shown, Figure 4 (b) VMD extracted spectral peak is shown;

[0080] Figure 5 Results of micro-motion parameter extraction of different methods are shown, wherein Figure 5 (a) OMP algorithm micro-motion parameter extraction result is shown, Figure 5 (b) VMD-OMP method micro-motion parameter extraction result is shown;

[0081] Figure 6 Results of VMD-atom expansion OMP method processing are shown, wherein Figure 6 (a) OMP decomposition result when η1=1 is shown, Figure 6 (b) OMP decomposition result when η2=0.75 is shown, Figure 6 (c) OMP decomposition result when η3=0.5 is shown, Figure 6 (d) OMP decomposition result when η4=0.25 is shown, Figure 6 (e) energy accumulation result is shown;

[0082] Figure 7 A certain six-rotor test unmanned aerial vehicle is shown;

[0083] Figure 8 Data preprocessing diagram is shown, wherein Figure 8 (a) unmanned aerial vehicle slow-time range image is shown, Figure 8 (b) signal time-domain diagram is shown, Figure 8 (c) signal time-frequency spectrum entropy is shown, Figure 8 (d) time-frequency spectrum entropy maximum slice signal time-frequency diagram is shown;

[0084] Figure 9 Results of VMD processing of data slices are shown, wherein Figure 9 (a) data slice JEM spectral peak feature is shown, Figure 9 (b) VMD extracted spectral peak is shown;

[0085] Figure 10 Results of VMD-atom expansion OMP based processing of data slices are shown, whereinFigure 10 (a) shows the OMP decomposition result when η1=1, Figure 10 (b) shows the OMP decomposition result when η2=0.75, Figure 10 (c) shows the OMP decomposition result when η3=0.5, Figure 10 (d) shows the OMP decomposition result when η4=0.25, Figure 10 (e) shows the energy accumulation result. DETAILED DESCRIPTION

[0086] The application will be further described below in conjunction with the embodiments and the accompanying drawings.

[0087] The application proposes a UAV micro-motion parameter extraction method under LPRF detection conditions, comprising the following steps:

[0088] Step 1: Calculate the time-frequency spectrum entropy

[0089] In order to select the signal with obvious HERM characteristics, first divide the UAV echo signal into M segments with a time length of t m , i.e. s(t) = [s d1 (t m ),..., s dm (t m ),..., s dM (t m )] T , wherein s(t) is the UAV echo signal, s dm (t m ) is the mth segment signal, m ∈ [1, M]; perform short-time Fourier transform STFT on the segmented signal s dm (t m ) to obtain its corresponding time-frequency representation STFF(t i , f j ) (for convenience, i represents the horizontal coordinate and j represents the vertical coordinate, and the lengths of i and j are both t m );

[0090] The time-frequency spectrum entropy is used to describe the energy distribution difference of different segmented signals in the time-frequency distribution

[0091]

[0092] In the formula: H m represents the time-frequency spectrum entropy of the mth segment signal, the time-frequency plane is equally divided into i × j time-frequency blocks with equal areas; the energy X ij of the (i, j)th time-frequency block corresponding to the time-frequency representation is the amplitude |STFT(t i , f j )|; P ijrepresents the probability of the energy of the (i, j)th time-frequency block relative to the energy of the entire time-frequency plane.

[0093] Second step: Extracting the rotation frequency

[0094] The present application uses VMD to extract the center frequency of each spectral peak, so as to calculate the spectral line interval of the micro-Doppler frequency of the unmanned aerial vehicle blade. According to the definition of the VMD algorithm, the decomposition problem of the multi-peak signal can be constructed as:

[0095]

[0096] In the formula: is the time-frequency spectrum entropy maximum segmented signal; is the partial derivative with respect to time t; s k is the time domain signal of the kth spectral peak; ω k is the center frequency of the kth spectral peak, which is a variable with respect to time t; K is the number of peaks. The problem described in formula (3) can be solved using the Alternating Direction Method of Multipliers (ADMM) to obtain the center frequency ω k of the multi-peak signal (LI Y, XIA W J, DONG S Q. Time-based multi-component irregular FM micro-Doppler signals decomposition via STVMD [J]. IET Radar, Sonar and Navigation, 2020, 14(10): 1502-1511.)

[0097] According to the relationship between the HERM spectral line interval and the rotation frequency, the rotation frequency estimate of the unmanned aerial vehicle rotor blade of N blades can be represented as:

[0098]

[0099] In the formula, mean(·) represents the mean operation, diff(·) represents the first-order difference operation, is the estimate of ω k , and N is the number of blades.

[0100] Third step: blade length extraction

[0101] Based on the extraction of the rotation frequency estimate , the present application proposes an atomic scaling OMP method to estimate the blade length. The dictionary established based on the unmanned aerial vehicle blade length estimation is where the ith atom under the xth scaling matching condition is:

[0102]

[0103] where η x is the scaling factor of the xth iteration, η1=1, i is the serial number of the atom in the dictionary (used for the scaling of the blade length). α, β, R0are the azimuth angle, the elevation angle and the distance of the blade rotation center relative to the radar at the initial moment, respectively; is the blade rotation frequency; f c is the radar transmitting carrier frequency; c is the speed of light; is the xth blade length estimate; is the phase of the nth blade at time t; λ is the wavelength; all the above are prior information. Compared with the blade length and the rotation frequency of the UAV, the initial phase parameter of the UAV blade has a lower impact factor. In order to reduce the dimension of the dictionary, the initial phase of the blade rotation is ignored.

[0104] In addition to the blade echo in the UAV echo signal, there are also main body echoes and noise signals. These signals will produce false peaks in the sparse domain, which will have a great impact on the accurate estimation of the blade length. To solve this problem, the scaling factor η x is used to adjust the atom in the dictionary D. In order to keep the whole unchanged, when the carrier frequency f c is reduced (increased) by η x times, the blade length estimate is increased (reduced) by η x times, while the false peak position is difficult to change according to the linear rule. By scaling the atom, the blade length estimate is selected, so as to reduce the negative impact of the interference signal on the blade length estimate . Based on this principle, the OMP algorithm is improved, and the algorithm process is as follows:

[0105] (1) Signal decomposition: when the carrier frequency is f c and the scaling factor is η1, the segmented signal with the maximum time-frequency spectrum entropy is sparsely decomposed, and the atom obtained after g iterations is normalized as The energy estimate maximum value of the ith atom under the first scaling matching condition is selected , and the corresponding length unit position is selected (at this time ). The relationship between the length unit position and the blade length estimate is:

[0106]

[0107] In the formula, l0 is the initial value for the iterative search; L is the total length of the search; g is a counter used for counting the iterations of i; Maximum value for energy estimation The corresponding g value; G is the maximum number of iterations.

[0108] (2) Scaling matching: When the scaling coefficient is η x Under the conditions of length element position The theoretical position of the length element after scaling and matching. Represented as:

[0109]

[0110] In the formula, x is the index of the scaling order, which is a positive integer greater than 1.

[0111] at this time, The corresponding energy values ​​are represented as follows:

[0112]

[0113] In the formula, δ is the threshold value, used to filter the energy extreme values ​​that meet the conditions; The scaling factor η x The maximum energy estimate of the i-th atom under the given conditions (obtained from Equation 8) The scaling factor η x The theoretical position of the length element after scaling and matching under the given conditions (obtained by Equation 7).

[0114] (3) Energy accumulation. The energy obtained above... and Energy at the corresponding position and Accumulate ( and According to Formula 8), select the maximum accumulated energy value. This step can be expressed as:

[0115]

[0116] The estimated values ​​obtained from equations (6) and (9) As the final estimate of blade length

[0117] Example: A method for estimating the rotational frequency and blade length of a UAV at low pulse repetition frequencies.

[0118] The parameters of the continuous wave radar used in the example are shown in Table 1, and the rotation diagram of the UAV blades is shown in the figure. Figure 2 As shown.

[0119] Table 1. Operating parameters of the radar in the simulation experiment

[0120]

[0121] Additive white Gaussian noise is added to the simulated unmanned aerial vehicle echo signal, and the signal-to-noise ratio after pulse compression is 10 dB. The time-domain representation of the echo signal is shown in Figure 3 (a). The data is segmented with 300 sampling points, Figure 3 (b) shows the normalized time-frequency spectrum entropy of the segmented signal. As shown in Figure 3 (c), the time-frequency representation of the segmented signal with the maximum time-frequency spectrum entropy contains relatively obvious HERM characteristics.

[0122] As shown in Figure 4 (a), the segmented signal with the maximum time-frequency spectrum entropy shows obvious JEM characteristics in the frequency domain, Figure 4 (b) is the result of the VMD algorithm processing. The VMD algorithm separates the peak signal, and the interval of the spectrum is calculated by formula (4) to be 49.8767 Hz, so the rotating frequency of the rotor blade is 24.9383 Hz, and the true rotating speed of the blade is 25 r / s, and the estimation error is 0.2468%. The simulation experiment shows that by using VMD to extract the center frequency point of each spectrum peak, the spectrum line interval of HERM can be calculated, and the rotating frequency of the rotor blade can be further extracted.

[0123] The method is further compared with OMP and VMD-OMP in performance. Figure 5 (a) is a frequency spectrum diagram using the OMP algorithm to estimate the rotating frequency and blade length simultaneously, wherein the estimated value of the rotating frequency is 18 Hz, the estimated value of the blade length is 0.671 m, and the time is 1.98 s. Figure 5 (b) is the energy distribution of the reconstructed signal using OMP to extract the blade length based on the VMD extracted rotating frequency, wherein the estimated value of the rotating frequency is 24.95 Hz, the estimated value of the blade length is 0.437 m, and the time is 0.11 s.

[0124] Figure 6 (a), (b), (c) and (d) are echo signal results using OMP decomposition under the conditions that the carrier frequency is f c = 10 GHz, the initial value l0 = 0.4, the equivalent length L = 0.4, the length unit (i.e. the maximum number of iterations) is 1000, the threshold value δ = 10 (the simulation experiment finds that when δ = 10, the algorithm can run more effectively), and the scaling coefficients are η1 = 1, η2 = 0.75, η3 = 0.5 and η4 = 0.25, respectively.

[0125] Figure 6 (a) is the length unit position and the energy estimation is maximum, and the corresponding blade length estimation value is and

[0126] Figure 6 (b) the estimated position of the scaled matched length unit is chosen as and and and The energy estimate corresponding to and The theoretical position of the length unit is given by equation (7) as and Since the estimated position of the length unit and are between and respectively, and the corresponding energies reach a maximum, we have and The energy estimate corresponding to and

[0127] Figure 6 (c) the estimated position of the scaled matched length unit is chosen as and and or The energy estimate corresponding to and The theoretical position of the length unit is given by equation (7) as and Since the estimated position of the length unit and are between and respectively, we have or The energy estimate corresponding to and

[0128] Figure 10 (d) the estimated position of the scaled matched length unit is outside the chosen range, in which case The energy estimate corresponding to The theoretical position of the length unit is given by equation (7) as and Since the estimated position of the length unit and are not between or respectively, we have The energy estimate corresponding to ​​​​

[0129] Figure 10 (e) Comparison of accumulated energy values. The position of the length unit when the energy estimate is at its maximum. and Comparison of the cumulative energy values ​​after three atomic scaling OMP processes, and the position of the length unit. At that time, the accumulated energy value is larger, so the final estimated value of its corresponding blade length is... As the final estimate of blade length

[0130] Table 2 compares the performance of three methods—OMP, VMD-OMP, and VMD-atomic scaling OMP—in extracting leaf length and leaf rotation frequency through 100 Monte Carlo experiments.

[0131] Table 2 Performance Comparison of Three Micro-Motion Parameter Extraction Methods

[0132]

[0133] Table 2 shows the computational performance of the three methods for blade length, rotation frequency, and time. The computer hardware parameters were: Intel(R) Core(TM)2 Duo CPU 2.4GHz, 32GB RAM, and Windows 7 64-bit operating system. Simulation results show that the absolute errors in blade length estimation for the OMP algorithm, VMD-OMP method, and the proposed method are 40%, 16%, and 5.8%, respectively; the absolute errors in blade rotation frequency estimation for the VMD algorithm and OMP algorithm are 0.36% and 16.32%, respectively. The computation time for the VMD-OMP method and the proposed method is comparable in magnitude and both outperform the OMP algorithm. Based on the comparison of the three sets of simulation experiments, the proposed method demonstrates the best overall performance in extracting micro-motion parameters under low pulse repetition frequency radar conditions.

[0134] Furthermore, the effectiveness of the proposed method was verified using experimental data, including the following steps:

[0135] Step 1: Data preprocessing.

[0136] ① Determine the state. The motion state (hovering, constant speed, acceleration) of the UAV is determined by analyzing the time and distance information collected by radar. When hovering, the distance of the UAV does not change with time. When moving at a constant speed, the distance of the UAV has a linear relationship with time. When accelerating, the distance of the UAV has a non-linear relationship with time.

[0137] ② Signal synthesis. The echo data collected by the radar are combined into a complex signal using formula (10) (this technique is well known to those skilled in the art).

[0138] ③ Slice selection. Take 1000 (i.e. divide 1s into 1000ms) as a slice unit, slice the whole data segment according to the total sampling time R, divide the whole data segment into T / 1000 time segments, calculate the time-frequency spectrum entropy (obtained from formula 2) of all segmented signals, and select the data slice with the maximum time-frequency spectrum entropy for further processing.

[0139] Second step: Extract the rotation frequency.

[0140] The VMD algorithm is used to process the data with the maximum time-frequency spectrum entropy, estimate the spectral interval of the HERM line of the UAV echo signal (obtained from formula 3), and calculate the rotation frequency of the UAV blade (obtained from formula 4).

[0141] Third step: Extract the blade length.

[0142] Based on the obtained blade rotation frequency, the atomic scaling OMP method is used to concentrate the micro-Doppler feature energy of the echo data (according to formulas 5 to 9), so as to estimate the blade length.

[0143] This method can realize the extraction and identification of the feature parameters of the micro-motion target under low pulse repetition frequency, and improve the overall working performance of the radar. Now the measured data of a six-rotor test UAV collected by an X-band radar of a certain institute is used for verification. The test UAV is shown in Figure 7 .

[0144] The X-band single-base narrowband radar and the six-rotor UAV parameters are shown in Tables 3 and 4.

[0145] Table 3 X-band radar operating parameters

[0146]

[0147] Table 4 Test six-rotor UAV parameters

[0148]

[0149] First step, data preprocessing.

[0150] The measured data is the I / Q two-way sampling data. In order to extract the information contained in the data to the greatest extent, the following method is adopted for data preprocessing.

[0151] (1) Determine the state. Determine the motion state of the UAV through time and distance information analysis. Figure 8 (a) is the slow-time range image of the UAV, which shows that the distance between the UAV and the radar is between 7476m and 7486m within the observation time, and the UAV is basically in a relatively static hovering state.

[0152] (2) Signal synthesis. The I / Q two-way data is combined into a complex signal by formula (10), and the time-domain representation of the complex signal is as follows Figure 8 (b) shows.

[0153] s(t) = I(t) + jQ(t) (9)

[0154] Where I(t) + jQ(t) is abbreviated as I and Q respectively, and the I / Q signal is the in-phase quadrature signal; I is in-phase; Q is quadrature, which is 90° out of phase with I (this technology is well known to those skilled in the art).

[0155] (3) Slice selection. The time-frequency spectrum entropy of all segmented signals is calculated for slice selection. Figure 8 (c) is the slice segmentation processing of the entire data segment with 1000 sampling points as the slice unit, the time-frequency spectrum entropy of all segmented signals is compared, and the data slice with the maximum time-frequency spectrum entropy (39th segment) is selected for subsequent processing; as shown in Figure 8 (d) shows that the data slice signal with the maximum time-frequency spectrum entropy contains obvious HERM characteristics in the time-frequency representation.

[0156] Second step, rotation frequency extraction.

[0157] Based on the above data preprocessing analysis, the VMD algorithm is used to process the data with the maximum time-frequency spectrum entropy. Figure 9 (a) is the data slice signal with the maximum time-frequency spectrum entropy, which shows obvious JEM characteristics in the frequency domain. Figure 9 (b) is the spectral peak of the data slice signal extracted by the VMD algorithm, and the echo signal can be decomposed into 12 modes, as shown in Table 5.

[0158] Table 5 Center frequency f of each IMF k

[0159]

[0160] The average estimation of the frequency interval of the center frequency of the corresponding mode in Table 5 is 81.6654 Hz, so the estimated value of the blade rotation frequency is 40.8327 Hz, and the blade speed is 41.7 r / s, which is relatively close, and the absolute error of the estimation is about 2.08%.

[0161] Third step, blade length extraction.

[0162] The blade length extraction method proposed in this paper is based on the VMD algorithm to obtain the blade rotation frequency, and uses the atomic scaling OMP method to concentrate the energy of the echo data micro-Doppler characteristics, so as to estimate the blade length.

[0163] Figure 10(a), (b), (c) and (d) are the results of using OMP algorithm to decompose the data slice signal under the condition that the carrier frequency is f c = 9.5 GHz, the initial value l0= 0.2, the equivalent length L= 2.8, the length unit (i.e. the maximum iteration number) is 2000, the threshold value δ= 10, the scaling coefficients are η1= 1, η2= 0.75, η3= 0.5 and η4= 0.25, respectively.

[0164] Figure 10 (a) is the energy estimation value at the length unit position of 100, 164 and 225, and the corresponding blade length estimation values are 0.3387, 0.4283 and 0.5138, respectively; Figure 10 (b) is the length unit estimation position after scaling and matching, which is 135, 219 and 295; Figure 10 (c) is the length unit estimation position after scaling and matching, which is 192, 320 and 440; Figure 10 (d) is the length unit estimation position after scaling and matching, which is 393, 661 and 899. Figure 10 (e) is the energy accumulation value after 3 times of atomic scaling OMP processing at the length unit position of 100, 164 and 225, and the corresponding blade length estimation value 0.3387 at the length unit position of 100 is finally taken as the final blade length estimation result.

[0165] From Figure 10 It can be seen from (e) that the blade length estimation value is 0.3387 m, while the blade length is 0.355 m, and the estimation error is 4.59%. Therefore, the method proposed in this paper can extract the blade length parameter of the six-rotor unmanned aerial vehicle in the actual scene.

Claims

1. A method for extracting micro-motion parameters of a UAV at a low pulse repetition frequency, characterized in that, Specifically comprising the following steps: First step: calculate time-frequency spectrum entropy To select signals with distinct HERM characteristics, the UAV echo signal is first divided into M segments of duration t. m The segment, i.e., s(t) = [s d1 (t m ), ..., s dm (t m ), ..., s dM (t m )] T Where s(t) is the UAV echo signal, s dm (t m Let be the m-th segment of the signal, where m ∈ [1, M]; for the segmented signal s dm (t m Performing a short-time Fourier transform (STFT) yields its corresponding time-frequency representation as STFT(t). i f j For ease of representation, i represents the x-coordinate and j represents the y-coordinate, with both i and j having a length of t. m ; Adopt time-frequency spectrum entropy to describe the energy distribution difference of different segmented signals in time-frequency distribution wherein: H m represents the spectral entropy of the mth segment signal, the time-frequency plane is equally divided into i x j time-frequency blocks with equal area; the energy X ij of the (i, j)th time-frequency block corresponds to the amplitude |STFT(t i , f j )| of the time-frequency representation; P ij represents the probability of the energy of the (i, j)th time-frequency block relative to the energy of the entire time-frequency plane. Second step: extract rotation frequency Use VMD to extract the center frequency of each spectrum peak, calculate the spectrum line interval of the micro-Doppler frequency of the UAV blade; the decomposition problem of multi-peak signal is constructed as: and where: is the partial derivative with respect to time t; s k is the time-domain signal of the kth spectral peak; ω k is the center frequency of the kth spectral peak, which is a function of time t; K is the number of peaks; is the segmented signal with the maximum time-frequency spectral entropy; The center frequency ω of the multi-peak signal k Using the alternating multiplier method, according to the relationship between the HERM spectral line interval and the rotation frequency, the estimated value of the rotation frequency of the unmanned aerial vehicle rotor blade of N blades is represented as: where mean( ) represents the averaging operation, and diff( ) represents the first order difference operation, is an estimate of ω k , and N is the number of blades. Third step: blade length extraction On the basis of extracting the rotation frequency estimate value The atomic expansion and contraction OMP method is proposed to estimate the blade length; the dictionary established based on the unmanned aerial vehicle blade length estimation The i-th atom under the x-th expansion and contraction matching condition is: where η x is the scaling factor of the xth iteration, η1= 1, i is the order of the atom in the dictionary; a, β, R0are the azimuth angle, elevation angle and distance of the blade rotation center relative to the radar at the initial moment, respectively; is the blade rotation frequency; f c is the radar transmitting carrier frequency; c is the speed of light; is the xth iteration of the blade length estimate; is the phase of the nth blade at time t; λ is the wavelength; the above are all prior information; in order to reduce the dimension of the dictionary, the initial phase of blade rotation is ignored; The OMP algorithm is improved as follows: (1) Signal decomposition: at a carrier frequency of f c When the scaling factor is η1, the segmented signal with the largest time-frequency entropy. The atoms obtained after sparse decomposition and g iterations are Its energy is normalized to Select the maximum energy estimate of the i-th atom under the first scaling matching condition. Corresponding length unit position at this time Its comparison with the estimated blade length The relationship is: where lo is the initial value of the iterative search; L is the total length of the search; g is a counter for the iteration count of i; maxima of the energy estimate corresponding g value; G is the maximum number of iterations; (2) Scale matching: under the condition of scale coefficient η x respectively, the length unit position is scaled, and the length unit theoretical position after scale matching is represented as: In the formula, x is the number of scaling times, which is a positive integer greater than 1; At this time, The corresponding energy value is expressed as: In the formula, δ is the threshold value; The scaling factor η x The maximum energy estimate of the i-th atom under the given conditions. The scaling factor η x The theoretical position of the length element after scaling and matching under certain conditions; (3) Energy accumulation; the above obtained with energy of the corresponding position with accumulation, select the maximum value of energy accumulation, this step is represented as: The resulting estimates from equations (6), (9) As a final estimate of the blade length 2. The method of claim 1, wherein, Scaling coefficients η1=1, η2=0.75, η3=0.5, η4=0.25, threshold δ=10.

3. The method of claim 1, wherein, Initial value l0=0.2m, equivalent length L=2.8m, maximum iteration number G=2000.

4. The method of claim 1, wherein the effectiveness of the method of extracting the UAV micro-motion parameters at a low pulse repetition frequency is verified. Specifically comprising the following steps: First step: data preprocessing ①Determine the state; analyze and determine the motion state of the UAV through the time and distance information collected by the radar. When hovering, the distance of the UAV does not change with time, when moving at a constant speed, the distance of the UAV is linearly related to time, and when accelerating, the distance of the UAV is nonlinearly related to time; ②Signal synthesis; combine the echo data collected by the radar into a complex signal through formula (10): s(t) = I(t) + jQ(t) (9) Where I(t) + jQ(t) is simply written as I and Q, and the I / Q signal is in-phase quadrature signal; the phase of Q is 90° different from that of I; ③Slice selection; take 1000 (i.e. divide 1s into 1000ms) as a slice unit, slice and segment the entire data segment according to the total sampling time T, divide the entire data segment into T / 1000 time slices, calculate the time-frequency spectrum entropy of all segmented signals according to formula (2), and select the data slice with the maximum time-frequency spectrum entropy for further processing; In the formula: H m Let X represent the time-frequency entropy of the m-th signal segment. This time-frequency plane is divided into i×j equal time-frequency blocks; the energy X corresponding to the (i,j)-th time-frequency block. ij The amplitude |STFT(t) is a time-frequency representation. i f j )|;P ij This represents the probability of the energy of the (i, j)th time-frequency block relative to the energy of the entire time-frequency plane; Second step: extract rotation frequency Use VMD algorithm to process the data with the maximum time-frequency spectrum entropy, estimate the spectrum interval of the HERM line of the UAV echo signal from formula (3), and calculate the rotation frequency of the UAV blade from formula (4); wherein: is the time-frequency spectrum entropy maximum segmented signal; is the partial derivative with respect to time t; s k is the time-domain signal of the kth spectral peak; ω k is the center frequency of the kth spectral peak, which is a variable with respect to time t; K is the number of peaks; According to the relationship between HERM spectrum interval and rotation frequency, the estimated value of the rotation frequency of the N-blade UAV rotor blade can be expressed as: where mean( ) represents the averaging operation, and diff( ) represents the first order difference operation, is an estimate of ω k , and N is the number of blades. Third step: extract blade length Based on the obtained blade rotation frequency, use the atomic scaling OMP method to concentrate the echo data micro-Doppler feature energy according to formulas 5 to 9, thereby estimating the blade length; Dictionary based on estimation of length of drone blade where the i-th atom in the x-th rescaling match condition is: where η x is the scaling factor of the xth iteration, η1= 1, i is the order of the atom in the dictionary; a, β, R0are the azimuth angle, elevation angle and distance of the blade rotation center relative to the radar at the initial moment, respectively; is the blade rotation frequency; f c is the radar transmitting carrier frequency; c is the speed of light; is the xth iteration of the blade length estimate; is the phase of the nth blade at time t; λ is the wavelength; all the above are prior information; the initial phase of the blade rotation is ignored to reduce the dimension of the dictionary; When the carrier frequency is f c , the expansion and contraction coefficient is η1, the sparse decomposition is performed on the segmented signal with the maximum time-frequency spectrum entropy , and the atom obtained after g times of iterative operation has an energy normalization of The energy estimation maximum value of the i-th atom under the first expansion and contraction matching condition is selected , and the corresponding length unit position is determined. At this time , the relationship between the atom and the blade length estimation value is: where lo is the initial value of the iterative search; L is the total length of the search; g is a counter for the iteration count of i; the energy estimate maximum the corresponding g value; G is the maximum number of iterations; Under the condition of scale coefficient η x The length unit position is scaled under the condition of scale coefficient η The scaled length unit position In the formula, x is the number of scaling times, which is a positive integer greater than 1; At this time, The corresponding energy value is expressed as: In the formula, δ is a threshold value; is a scaling coefficient η x Energy estimation maximum value of the i-th atom under the condition that the i-th atom is in the position of the i-th length unit, is a scaling coefficient η x Theoretical position of the length unit after scaling and matching under the condition that the i-th atom is in the position of the i-th length unit. The above obtained With The energy of the corresponding position With The energy accumulation is added, and the maximum energy accumulation value is selected. This step is represented as: The resulting estimates from equations (6), (9) As a final estimate of the blade length 5. The method of claim 4, wherein, Scaling coefficients η1=1, η2=0.75, η3=0.5, η4=0.25, threshold δ=10.

6. The method of claim 4, wherein the low pulse repetition frequency UAV fine motion parameter extraction method is characterized by, Initial value l0=0.2m, equivalent length L=2.8m, maximum iteration number G=2000.

Images