CN112287796B - Radiation source identification method based on VMD-Teager energy operator - Google Patents

Radiation source identification method based on VMD-Teager energy operator Download PDF

Info

Publication number
CN112287796B
CN112287796B CN202011144744.2A CN202011144744A CN112287796B CN 112287796 B CN112287796 B CN 112287796B CN 202011144744 A CN202011144744 A CN 202011144744A CN 112287796 B CN112287796 B CN 112287796B
Authority
CN
China
Prior art keywords
signal
radiation source
component
components
phase component
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202011144744.2A
Other languages
Chinese (zh)
Other versions
CN112287796A (en
Inventor
金燕华
李秋雪
颜松涛
阎啸
王童樾
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
University of Electronic Science and Technology of China
Original Assignee
University of Electronic Science and Technology of China
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by University of Electronic Science and Technology of China filed Critical University of Electronic Science and Technology of China
Priority to CN202011144744.2A priority Critical patent/CN112287796B/en
Publication of CN112287796A publication Critical patent/CN112287796A/en
Application granted granted Critical
Publication of CN112287796B publication Critical patent/CN112287796B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/08—Feature extraction
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/214—Generating training patterns; Bootstrap methods, e.g. bagging or boosting
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/24—Classification techniques
    • G06F18/241—Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches
    • G06F18/2411—Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches based on the proximity to a decision surface, e.g. support vector machines
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/02—Preprocessing
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/12—Classification; Matching

Abstract

The invention discloses a radiation source identification method based on VMD-Teager energy operator, which comprises the steps of respectively collecting signal samples of each radiation source, respectively adopting an improved VMD algorithm to carry out signal decomposition on an in-phase component and an orthogonal component of each signal sample, selecting a low-frequency modal component to reconstruct to obtain a reconstructed in-phase component and a reconstructed orthogonal component, respectively calculating the energy of the reconstructed in-phase component and the reconstructed orthogonal component based on the Teager energy operator, then calculating the modal component corresponding to the in-phase component and the orthogonal component of each signal sample and the box dimension of the reconstructed in-phase component and the reconstructed orthogonal component, constructing the energy and the box dimension to obtain a characteristic vector of the signal sample, taking the characteristic vector as input, taking the serial number of the corresponding radiation source as output, training a classification model to obtain a radiation source identification model, and obtaining the characteristic vector of the signal when the radiation source is required to be identified, and inputting the radiation source identification model to obtain an identification result. The invention can effectively improve the identification rate of the radiation source in the environment with low signal-to-noise ratio.

Description

Radiation source identification method based on VMD-Teager energy operator
Technical Field
The invention belongs to the technical field of radiation source identification, and particularly relates to a radiation source identification method based on a VMD-Teager energy operator.
Background
The communication radiation source individual identification technology has important significance for improving the safety of a wireless communication system and improving the military communication reconnaissance capability. By analyzing the characteristics of the communication signal of the radio device, i.e. the radio frequency fingerprint, the individual sources of the radiation sources of the signal can be distinguished. Radiation source individual identification has been widely used in radio, internet of things security, and the like.
At present, in the aspect of radio frequency fingerprint feature extraction, many different algorithms have been developed, which can be summarized as two main types extracted from time domain and frequency domain: (1) wavelet features extracted from the frequency domain, Hilbert transform, integral bispectrum, frequency domain features of spurious components, and the like. (2) And extracting mean value characteristics, natural measurement and the like from the time domain to be used as the radio frequency fingerprint characteristics of the radiation source. The radio frequency fingerprint feature extraction algorithms can effectively identify the individual radiation source in an ideal environment, but the radio frequency fingerprint feature extraction algorithms are not researched specially aiming at the low signal to noise ratio environment at present, and the identification effect is not ideal in the low signal to noise ratio environment by using the existing method. In the document "Klein R W, sample M A, Mendenhall M J. application of Wavelet-Based RF matching to enhanced Wireless Network Security [ J ]. Journal of Communications and Networks,2012,11(6): 544) 555", Wavelet domain fingerprint features Based on dual-tree Wavelet transform were extracted from the non-transient preamble response of the signal, but at a signal-to-noise ratio of 5dB, the recognition rate was below 60%. A communication radiation source individual identification method based on empirical mode decomposition [ J ] is reported in Chinese institute of electronic science, 2013(04):393 and 397. "utilizes an Empirical Mode Decomposition (EMD) method to extract frequency domain characteristics of a steady-state signal spurious component, and the identification rate is about 78% in an environment with a signal-to-noise ratio of 5 dB. The document "Lin Y, Zhu X, Zheng Z, et al, the identification method of wireless device based on dimensional reduction and machine learning [ J ]. Journal of Supercomputing, 2017" using Hilbert feature, proposes a method based on a combination of dimensionality reduction and machine learning, with an average recognition rate of 82% in the range of 0-20 dB. The document "Ding L, Wang S, Wang F, et al. specific Emitter Identification of visual Networks [ J ]. IEEE Communications Letters,2018,22(12): 2591-. The document "Ali A M, Uzundaukan E, Kara A. Association of Features and classes for Bluetooth RF converting [ J ]. IEEE Access,2019,7: 50524-. The documents "Wang X, Zhang Y, Zhang H, et al, identification and authentication for Wireless transmission security based on RF-DNA finger print [ J ]. EURASIP Journal on Wireless Communications and Networking,2019 (1)" propose to add mean features to a radio frequency fingerprint, to identify by means of mean, bias and kurtosis statistics and by means of multiple discriminant analysis, the recognition rate being around 80% at a signal-to-noise ratio of 5 dB.
Among the algorithms involved in radio frequency fingerprint feature extraction, EMD is an algorithm that decomposes a multi-component signal into several individual single components, i.e., eigenmode functions and a residual signal, but the method involves envelopes, mode mixing, end effects and unexplained negative frequencies caused by the Hilbert transform. The Local Mean Decomposition (LMD) method improved according to EMD also has problems of distortion components, modal mixing, and long decomposition time. In response to the problem, in 2014, dragomirtski and Zosso propose Variational Mode Decomposition (VMD), and since the algorithm is combined with wiener filtering, the algorithm has good noise robustness and has been widely applied to detection signal variation and feature extraction. And the VMD algorithm can effectively alleviate the defects of modal aliasing and noise interference and the like generated by EMD and LMD decomposition, and is widely applied to mechanical fault diagnosis, detection signal change and feature extraction at present.
The document "Aghnaiya A, Ali AM, Kara A. variable Mode composition Based Radio Frequency fingerprint identification of Bluetooth Devices [ J ]. IEEE Access,2019, PP (99): 1-1." applies the VMD algorithm to the Radio Frequency fingerprint identification, decomposes the Bluetooth transient signal into a series of band-limited modes using the VMD, and then reconstructs the transient signal according to these modes. And extracting high-order statistical features from the complex form of the reconstructed transient, and then identifying the BT device by using a linear support vector machine classifier. Meanwhile, based on a VMD algorithm, high-order modal features are directly extracted from the reconstructed transient signals, the average recognition rate is 91% in 0-5dB, but the average recognition rate is only 70.1% in-5-0 dB. It can be found that although the VMD algorithm has good noise robustness, the VMD is still under study in the field of radiation source identification under low signal-to-noise ratio environment.
Disclosure of Invention
The invention aims to overcome the defects of the prior art, provides a radiation source identification method based on a VMD-Teager energy operator, and improves the identification rate of a radiation source in a low signal-to-noise ratio environment.
In order to achieve the purpose, the radiation source identification method based on the VMD-Teager energy operator comprises the following steps:
s1: for M radiation sources to be identified, respectively acquiring a plurality of signal samples for each radiation source, and recording the d signal sample as xd(t)＝Id(t)+jQd(t)，Id(t)、Qd(t) respectively represent signal samples xd(t) in-phase and quadrature components, D ═ 1,2, …, D, t, time;
s2: for the in-phase component I of each signal sampled(t), orthogonal component Qd(t) respectively adopting improved VMD algorithm to carry out signal decomposition to obtain N modal components
Figure GDA0003501120990000038
Figure GDA0003501120990000038
1,2, …, N, and the specific steps of signal decomposition by using the modified VMD algorithm include:
s2.1: the number n of initialization modes is 2;
s2.2: recording the signal to be decomposed as x (t), and performing VMD decomposition on the signal to be decomposed x (t) by adopting a VMD algorithm to obtain n modal components imfi′(t)，i′＝1,2,…,n；
S2.3: calculating each signal x (t) to be decomposed and each modal component imfi′(t) coefficient of correlation Ri′；
S2.4: judging n correlation coefficients Ri′Whether the mode number is decreased along with the increase of the mode number i', if so, the step S2.5 is carried out, otherwise, the step S2.6 is carried out;
s2.5: judging whether N is less than N, wherein N represents the preset number of modes needing to be decomposed, if so, entering a step S2.6, otherwise, entering a step S2.7;
s2.6: making n equal to n +1, and returning to the step S2.2;
s2.7: judging whether N is larger than N, if so, entering step S2.8, otherwise, entering step S2.10;
s2.8: calculating correlation coefficients R corresponding to the n modal componentsi′Selecting 2 correlation coefficients with the minimum difference value from the difference values of the adjacent 2 correlation coefficients, combining the corresponding 2 modal components, and summing the 2 correlation coefficients to obtain the correlation coefficient corresponding to the modal component;
s2.9: making n equal to n-1, and returning to step S2.7;
s2.10: taking the current N modal components as a signal decomposition result;
s3: for the in-phase component I of each signal sampled(t), orthogonal component Qd(t) N modal components obtained by decomposition
Figure GDA0003501120990000031
Selecting the first N' modal components to be combined to obtain a reconstructed in-phase component
Figure GDA0003501120990000032
Reconstructing the orthogonal components
Figure GDA0003501120990000033
The value of N' is set according to actual needs;
s4: reconstructing the in-phase components using Teager energy operator pairs, respectively
Figure GDA0003501120990000034
Reconstructing the orthogonal components
Figure GDA0003501120990000035
Demodulating to obtain corresponding instantaneous amplitude
Figure GDA0003501120990000036
Then respectively calculating and reconstructing in-phase components
Figure GDA0003501120990000037
Reconstructing the orthogonal components
Figure GDA0003501120990000041
Energy of
Figure GDA0003501120990000042
Figure GDA0003501120990000043
Figure GDA0003501120990000044
S5: for each signal sample, its in-phase component I is calculated separatelyd(t), orthogonal component QdN modal components of (t)
Figure GDA0003501120990000045
Respectively, of
Figure GDA0003501120990000046
And corresponding reconstructed in-phase component
Figure GDA0003501120990000047
Reconstructing the orthogonal components
Figure GDA0003501120990000048
Respectively, of
Figure GDA0003501120990000049
S6: constructing and obtaining a feature vector by adopting the energy and box dimension of each signal sample
Figure GDA00035011209900000410
Taking the characteristic vector as input, and taking the radiation source serial number corresponding to the signal sample as a label to construct and obtain D training samples;
s7: training the M classification model by using the training sample constructed in the step S6 to obtain a radiation source identification model;
s8: when the radiation source needs to be identified, the radiation source signal x '(t) ═ I' (t) + jQ '(t) is acquired, I' (t) and Q '(t) respectively represent the in-phase component and the quadrature component of the signal x' (t), and the same method in step S2 is adopted to carry out signal decomposition to obtain N modal components of the in-phase component I '(t) and the quadrature component Q' (t)
Figure GDA00035011209900000422
Figure GDA00035011209900000411
Reconstructing by the same method in the step S3 to obtain a reconstructed in-phase component
Figure GDA00035011209900000412
Reconstructing the orthogonal components
Figure GDA00035011209900000413
The reconstructed in-phase component is calculated by the same method in step S4
Figure GDA00035011209900000414
Reconstructing the orthogonal components
Figure GDA00035011209900000415
Energy E ofI′、EQ′N modal components of the in-phase component I ' (t) and the quadrature component Q ' (t) of the signal x ' (t) are calculated by the same method as in step S5
Figure GDA00035011209900000416
Figure GDA00035011209900000417
Of box dimension Yi I′、Yi Q′And reconstructing the in-phase component
Figure GDA00035011209900000418
Reconstructing the orthogonal components
Figure GDA00035011209900000419
Of (2) box dimension
Figure GDA00035011209900000420
Constructing to obtain a feature vector
Figure GDA00035011209900000421
And inputting the data into the radiation source recognition model trained in the step S7 to obtain a radiation source recognition result.
The invention relates to a radiation source identification method based on VMD-Teager energy operator, which comprises the steps of respectively collecting signal samples of each radiation source, respectively adopting an improved VMD algorithm to carry out signal decomposition on an in-phase component and an orthogonal component of each signal sample, selecting a low-frequency modal component to carry out signal reconstruction to obtain a reconstructed in-phase component and a reconstructed orthogonal component, respectively calculating the energy of the reconstructed in-phase component and the reconstructed orthogonal component based on the Teager energy operator, then calculating the box dimensions of the modal component corresponding to the in-phase component and the orthogonal component of each signal sample, and the box dimensions of the reconstructed in-phase component and the reconstructed orthogonal component, constructing the energy and the box dimensions to obtain the characteristic vector of the signal sample, taking the characteristic vector as input, taking the radiation source serial number corresponding to the signal sample as output, training a classification model to obtain a radiation source identification model, and when the radiation source identification is required, and obtaining the characteristic vector of the signal by adopting the same method, and inputting the radiation source identification model to obtain an identification result.
The invention has the following beneficial effects:
1) the VMD algorithm is improved, and the aliasing phenomenon of each mode after the suppression decomposition is effectively enhanced;
2) the radio frequency fingerprint feature extraction is carried out by combining the Teager energy operator and the improved VMD algorithm, and the radiation source identification rate under the environment with low signal-to-noise ratio is effectively improved.
Drawings
FIG. 1 is a flow chart of an embodiment of the radiation source identification method based on VMD-Teager energy operators according to the present invention;
FIG. 2 is a flow chart of signal decomposition using the modified VMD algorithm of the present invention;
FIG. 3 is a decomposition result of the decomposition performed by the conventional VMD algorithm in the present embodiment;
FIG. 4 is a decomposition result of the improved VMD algorithm for decomposition in the present embodiment;
FIG. 5 is a schematic diagram of a feature vector set of each intercom in this embodiment;
fig. 6 is a statistical chart of the hybrid identification rate of the present invention under different signal-to-noise ratios in this embodiment. .
Detailed Description
The following description of the embodiments of the present invention is provided in order to better understand the present invention for those skilled in the art with reference to the accompanying drawings. It is to be expressly noted that in the following description, a detailed description of known functions and designs will be omitted when it may obscure the subject matter of the present invention.
FIG. 1 is a flow chart of an embodiment of the radiation source identification method based on VMD-Teager energy operators according to the present invention. As shown in FIG. 1, the radiation source identification method based on VMD-Teager energy operator of the present invention comprises the following specific steps:
s101: collecting radiation source signal samples:
for M radiation sources to be identified, respectively acquiring a plurality of signal samples for each radiation source, and recording the d signal sample as xd(t)＝Id(t)+jQd(t)，Id(t)、Qd(t) respectively represent signal samples xdThe in-phase and quadrature components of (t), D1, 2, …, D, t representing time.
S102: signal sample decomposition using the modified VMD algorithm:
then for the in-phase component I of each signal sampled(t), orthogonal component Qd(t) respectively adopting improved VMD algorithm to carry out signal decomposition to obtain N modal components
Figure GDA0003501120990000051
i＝1,2,…,N。
In practical application, the acquired signals are decomposed directly by using the VMD algorithm, so that the phenomenon of insufficient or excessive decomposition exists, and occasionally, the phenomenon of modal aliasing also occurs. Therefore, the invention improves the VMD algorithm by utilizing the correlation coefficient of the mode and the source signal after the VMD decomposition and the center frequency of the mode, and adopts the improved VMD algorithm to decompose the signal. Fig. 2 is a flow chart of signal decomposition using the improved VMD algorithm of the present invention. As shown in fig. 2, the specific steps of signal decomposition by using the improved VMD algorithm in the present invention include:
s201: the number of initialization modes n is 2.
S202: VMD decomposition of the signal:
recording the signal to be decomposed as x (t), and performing VMD decomposition on the signal to be decomposed x (t) by adopting a VMD algorithm to obtain n modal components imfi′(t), i ═ 1,2, …, n. Then there are:
Figure GDA0003501120990000061
the specific principle and steps of the VMD algorithm can be found in the literature "Lv, Zhongliang, Tang, Baoping, Zhou, Yi, and Zhou," A Novel Method for Mechanical factory floor Based on spatial Mode Decomposition and Multikernel Support Vector machine "," Shock and Vibration (2016) (2016: 1-11.Web. ", and the specific parameters of the VMD algorithm during operation can be set according to actual needs.
S203: calculating a correlation coefficient:
calculating each signal x (t) to be decomposed and each modal component imfi′(t) coefficient of correlation Ri′。
S204: judging n correlation coefficients Ri′Whether to decrement with the increase of the modality number i', if so, step S205 is entered, otherwise, step S206 is entered.
S205: and judging whether N is less than N, wherein N represents the preset number of modes needing to be decomposed, if so, entering step S206, otherwise, entering step S207.
S206: let n be n +1, return to step S202.
S207: and judging whether N is larger than N, if so, entering step S208, and otherwise, entering step S210.
S208: and (3) modal component reconstruction:
calculating correlation coefficients R corresponding to the n modal componentsi′Selecting 2 correlation coefficients with the minimum difference value from the difference values of the adjacent 2 correlation coefficients, combining the corresponding 2 modal components, and summing the 2 correlation coefficients to obtain the correlation coefficient corresponding to the modal component. It is apparent that after the combined reconstruction, the n modal components are converted to n-1 modal components.
S209: let n be n-1, return to step S207.
S210: obtaining a signal decomposition result:
and taking the current N modal components as a signal decomposition result.
S103: signal reconstruction:
in the invention, the low-frequency modal component is selected for signal reconstruction, and since the modal components obtained by signal decomposition are arranged in ascending order according to the central frequency, the in-phase component I of each signal sampled(t), orthogonal component Qd(t) N modal components obtained by decomposition
Figure GDA0003501120990000071
Selecting the first N' modal components to be combined to obtain a reconstructed in-phase component
Figure GDA0003501120990000072
Reconstructing the orthogonal components
Figure GDA0003501120990000073
The value of N' is set according to actual needs. Reconstructing in-phase components
Figure GDA0003501120990000074
Reconstructing the orthogonal components
Figure GDA0003501120990000075
Are respectively:
Figure GDA0003501120990000076
Figure GDA0003501120990000077
s104: extracting energy characteristics of the reconstructed signal by using a Teager energy operator:
the Teager energy operator is a nonlinear energy operator, has better time resolution for the transient change of a discrete time signal, can effectively calculate the transient amplitude and the transient frequency of the signal, and can analyze and track the energy of the narrow-band signal through a simple mathematical theory.
Reconstructing the in-phase components using Teager energy operator pairs, respectively
Figure GDA0003501120990000078
Reconstructing the orthogonal components
Figure GDA0003501120990000079
Demodulating to obtain corresponding instantaneous amplitude
Figure GDA00035011209900000710
Instantaneous amplitude
Figure GDA00035011209900000711
The calculation formula of (2) is as follows:
Figure GDA00035011209900000712
wherein the content of the first and second substances,
Figure GDA00035011209900000713
respectively representing reconstructed in-phase components
Figure GDA00035011209900000714
First and second derivatives of (a).
Instantaneous amplitude
Figure GDA00035011209900000715
The calculation formula of (2) is as follows:
Figure GDA0003501120990000081
wherein the content of the first and second substances,
Figure GDA0003501120990000082
respectively representing reconstructed orthogonal components
Figure GDA0003501120990000083
First and second derivatives of (a).
Then respectively calculating and reconstructing in-phase components
Figure GDA0003501120990000084
Reconstructing the orthogonal components
Figure GDA0003501120990000085
Energy of
Figure GDA0003501120990000086
Figure GDA0003501120990000087
Figure GDA0003501120990000088
S105: extracting box-dimension features of the reconstructed signal:
for each signal sample, its in-phase component I is calculated separatelyd(t), orthogonal component QdN modal components of (t)
Figure GDA0003501120990000089
Respectively, of
Figure GDA00035011209900000810
Andcorresponding reconstructed in-phase component
Figure GDA00035011209900000811
Reconstructing the orthogonal components
Figure GDA00035011209900000812
Respectively, of
Figure GDA00035011209900000813
The fractal box dimension of the signal is a nonlinear characteristic characterization mode of the signal. In modal component
Figure GDA00035011209900000814
For example, its box dimension
Figure GDA00035011209900000815
The calculation method of (2) is briefly described as follows:
let the maximum side length of the square (box) be Imax：
Figure GDA00035011209900000816
Wherein w represents a modal component
Figure GDA00035011209900000817
The signal length of (2).
To pair
Figure GDA00035011209900000818
Performing preprocessing to shift the whole signal sequence upwards so that the minimum value of the sequence is shifted on the x axis, and the minimum value of the signal is zero at the moment, so as to obtain a new signal sequence Z (t):
Figure GDA00035011209900000819
resampling the new signal sequence Z (t) by using an interpolation function to make the total number of points of the sequence Imax+1, the resulting interpolation sequence is X (t). Setting the specific ratio to stretch the resampled sequence x (t) so that the maximum value is equal to the sequence length, namely:
Figure GDA00035011209900000820
when the sequence is of length Imax+1, maximum value of Imax. Dividing the length of the square grid along the x-axis direction of a coordinate axis to form R strips:
Figure GDA0003501120990000091
wherein L represents the side length of a square (box).
Respectively calculating the number of squares occupied by each sequence in each strip, and setting the segment of the sequence X (t) in the r-th strip as Xr(t) the number of squares contained therein is yrThen, there are:
Figure GDA0003501120990000092
wherein the content of the first and second substances,
Figure GDA0003501120990000093
which means that the rounding is made up,
Figure GDA0003501120990000094
indicating a rounding down.
The total number of squares, i.e. box dimensions, is then:
Figure GDA0003501120990000095
and similarly, the box dimensions of other signals can be calculated.
S106: constructing a training sample:
constructing and obtaining a feature vector by adopting the energy and box dimension of each signal sample
Figure GDA0003501120990000096
It can be seen that the feature vector is a 2+2(N +1) -dimensional vector. And D training samples are constructed by taking the characteristic vector as input and the radiation source serial numbers corresponding to the signal samples as labels.
S107: training a classification model:
and (4) training the M classification model by adopting the training sample constructed in the step (S106) to obtain a radiation source identification model.
S108: and (3) identifying a radiation source:
when the radiation source needs to be identified, the radiation source signal x '(t) ═ I' (t) + jQ '(t) is acquired, I' (t) and Q '(t) respectively represent the in-phase component and the quadrature component of the signal x' (t), and the same method in step S102 is adopted to carry out signal decomposition to obtain N modal components of the in-phase component I '(t) and the quadrature component Q' (t)
Figure GDA0003501120990000097
Figure GDA0003501120990000098
Reconstructing by the same method in step S103 to obtain a reconstructed in-phase component
Figure GDA0003501120990000099
Reconstructing the orthogonal components
Figure GDA00035011209900000910
The reconstructed in-phase component is calculated by the same method in the step S104
Figure GDA00035011209900000911
Reconstructing the orthogonal components
Figure GDA00035011209900000912
Energy E ofI′、EQ′N modal components of the in-phase component I ' (t) and the quadrature component Q ' (t) of the signal x ' (t) are calculated by the same method in step S105
Figure GDA00035011209900000913
Of box dimension Yi I′、Yi Q′And reconstructing the in-phase component
Figure GDA00035011209900000914
Reconstructing the orthogonal components
Figure GDA00035011209900000915
Of (2) box dimension
Figure GDA00035011209900000916
Constructing to obtain a feature vector
Figure GDA00035011209900000917
And inputting the data into the radiation source recognition model trained in the step S107 to obtain a radiation source recognition result.
Examples
In order to better illustrate the technical effects of the invention, the invention is experimentally verified by using a specific example. In the embodiment, interphones produced by five factories are selected as radiation sources of signals and are respectively marked as A (A1/A2), B (B1/B2), C (C1/C2), D (D1/D2) and E (E1/E2). And setting the sampling rate to be 1MS/s to collect the signal samples of the five pairs of interphones. Under the noiseless environment, the signal sample collection is respectively carried out under the conditions of 1 meter (signal-to-noise ratio 10dB) and 5 meters (signal-to-noise ratio 5dB) from the interphone. During each acquisition, 400 signal samples are respectively acquired for each device, the number of sampling points of each sample is 1024, and 10 × 400 samples in total form a pure signal set.
The effectiveness of the improved VMD algorithm proposed by the present invention is first verified. The signals of a certain existing interphone are decomposed by adopting a traditional VMD algorithm and the improved VMD algorithm provided by the invention respectively. Fig. 3 is a decomposition result of the decomposition performed by the conventional VMD algorithm in the present embodiment. Fig. 4 is a decomposition result of the decomposition performed by the modified VMD algorithm in the present embodiment. Comparing fig. 3 and fig. 4, it can be known that the improved VMD algorithm improves the phenomenon of modal aliasing, which provides great convenience for extracting the radio frequency fingerprint feature information of the signal.
The method of the invention is adopted to carry out the radio frequency fingerprint feature extraction on a pure signal set with the acquisition distance of 1 meter to form a feature vector set. Fig. 5 is a schematic diagram of a feature vector set of each intercom in this embodiment. It can be seen from fig. 5 that the differences of the radio frequency fingerprint characteristics of the interphones are obvious, which lays a good foundation for the subsequent identification.
In this embodiment, an M classification model is constructed based on an SVM (support vector machine), an obtained feature vector set is divided into a training set and a test set, the M classification model is trained by using the training set to obtain a radiation source recognition model, the accuracy of the radiation source recognition model obtained by using the test set is tested, and the mixed recognition rate is calculated to be 93.4%. Then, the radio frequency fingerprint feature extraction is carried out respectively based on the signal decomposition results of the traditional VMD algorithm and the improved VMD algorithm provided by the invention, and the mixed identification rate obtained by different signal decomposition algorithms is counted. Table 1 is a comparison table of the mixed recognition rates obtained by different signal decomposition methods in this embodiment.
Collecting distance The invention Traditional VMD algorithm
1 meter (10dB) 93.4％ 86.4％
5 m (5dB) 90.2％ 70.1％
TABLE 1
And then, adopting simulation software, adding Gaussian white noise into the acquired pure signal simulation, setting the range of the signal-to-noise ratio to be-10 to 10dB, and forming a signal sample set with different signal-to-noise ratios. Fig. 6 is a statistical chart of the hybrid identification rate of the present invention under different signal-to-noise ratios in this embodiment. Then, the VMD-HOS algorithm (see the literature, "Alghannai Aghnaiya, Yaser Dalveren, and Ali Kara. on the Performance of spatial Mode Decomposition-Based Radio Frequency recognition of Bluetooth Devices [ J ]. Sensors (base, Switzerland)20.6(2020):1704. Web.") was used as a comparison method to count the hybrid recognition rates obtained by the present invention and the comparison method under different SNR conditions. Table 2 is a comparison table of the mixed recognition rates obtained by the present invention and the comparison method under different snr conditions in this example.
Figure GDA0003501120990000111
TABLE 2
By comparing table 1 with fig. 6, it is found from the experimental result analysis of the collected data with different signal-to-noise ratios and the simulated signal-to-noise ratio that the recognition rate of the actual signal-to-noise ratio is not much different from the recognition rate under simulation. The validity of the scheme is verified. As can be seen from the table 2, under the environment with a low signal to noise ratio of-5 dB to 0dB, the identification rate of the invention is higher than that of the existing VMD-HOS algorithm, and the effectiveness of the invention in the radio frequency fingerprint extraction under the environment with a low signal to noise ratio is verified.
Although illustrative embodiments of the present invention have been described above to facilitate the understanding of the present invention by those skilled in the art, it should be understood that the present invention is not limited to the scope of the embodiments, and various changes may be made apparent to those skilled in the art as long as they are within the spirit and scope of the present invention as defined and defined by the appended claims, and all matters of the invention which utilize the inventive concepts are protected.

Claims (1)

1. A radiation source identification method based on a VMD-Teager energy operator is characterized by comprising the following steps:
s1: for M radiation sources to be identified, respectively acquiring a plurality of signal samples for each radiation source, and recording the d signal sample as xd(t)＝Id(t)+jQd(t)，Id(t)、Qd(t) respectively represent signal samples xd(t) in-phase and quadrature components, D ═ 1,2, …, D, t, time;
s2: for the in-phase component I of each signal sampled(t), orthogonal component Qd(t) respectively adopting improved VMD algorithm to carry out signal decomposition to obtain N modal components
Figure FDA0003501120980000011
The specific steps of signal decomposition by adopting the improved VMD algorithm comprise:
s2.1: the number n of initialization modes is 2;
s2.2: recording the signal to be decomposed as x (t), and performing VMD decomposition on the signal to be decomposed x (t) by adopting a VMD algorithm to obtain n modal components imfi′(t)，i′＝1,2,…,n；
S2.3: calculating each signal x (t) to be decomposed and each modal component imfi′(t) coefficient of correlation Ri′；
S2.4: judging n correlation coefficients Ri′Whether the mode number is decreased along with the increase of the mode number i', if so, the step S2.5 is carried out, otherwise, the step S2.6 is carried out;
s2.5: judging whether N is less than N, wherein N represents the preset number of modes needing to be decomposed, if so, entering a step S2.6, otherwise, entering a step S2.7;
s2.6: making n equal to n +1, and returning to the step S2.2;
s2.7: judging whether N is larger than N, if so, entering step S2.8, otherwise, entering step S2.10;
s2.8: calculating correlation coefficients R corresponding to the n modal componentsi′Selecting 2 correlation coefficients with the minimum difference value from the difference values of the adjacent 2 correlation coefficients, combining the corresponding 2 modal components, and summing the 2 correlation coefficients to obtain the correlation coefficient corresponding to the modal componentCounting;
s2.9: making n equal to n-1, and returning to step S2.7;
s2.10: taking the current N modal components as a signal decomposition result;
s3: for the in-phase component I of each signal sampled(t), orthogonal component Qd(t) N modal components obtained by decomposition
Figure FDA0003501120980000012
Selecting the first N' modal components to be combined to obtain a reconstructed in-phase component
Figure FDA0003501120980000013
Reconstructing the orthogonal components
Figure FDA0003501120980000014
The value of N' is set according to actual needs;
s4: reconstructing the in-phase components using Teager energy operator pairs, respectively
Figure FDA0003501120980000015
Reconstructing the orthogonal components
Figure FDA0003501120980000016
Demodulating to obtain corresponding instantaneous amplitude
Figure FDA0003501120980000017
Then respectively calculating and reconstructing in-phase components
Figure FDA0003501120980000018
Reconstructing the orthogonal components
Figure FDA0003501120980000021
Energy of
Figure FDA0003501120980000022
Figure FDA0003501120980000023
Figure FDA0003501120980000024
S5: for each signal sample, its in-phase component I is calculated separatelyd(t), orthogonal component QdN modal components of (t)
Figure FDA0003501120980000025
Respectively, of
Figure FDA0003501120980000026
And corresponding reconstructed in-phase component
Figure FDA0003501120980000027
Reconstructing the orthogonal components
Figure FDA0003501120980000028
Respectively, of
Figure FDA0003501120980000029
S6: constructing and obtaining a feature vector by adopting the energy and box dimension of each signal sample
Figure FDA00035011209800000210
Taking the characteristic vector as input, and taking the radiation source serial number corresponding to the signal sample as a label to construct and obtain D training samples;
s7: training the M classification model by using the training sample constructed in the step S6 to obtain a radiation source identification model;
s8: when radiation source identification is needed, radiation source signals x '(t) ═ I' (t) + jQ '(t) are acquired, and I' (t) and Q '(t) respectively represent an in-phase component and a quadrature component of the signals x' (t)Performing signal decomposition by the same method as in step S2 to obtain N modal components of the in-phase component I '(t) and the quadrature component Q' (t)
Figure FDA00035011209800000211
Figure FDA00035011209800000212
Reconstructing by the same method in the step S3 to obtain a reconstructed in-phase component
Figure FDA00035011209800000213
Reconstructing the orthogonal components
Figure FDA00035011209800000214
The reconstructed in-phase component is calculated by the same method in step S4
Figure FDA00035011209800000215
Reconstructing the orthogonal components
Figure FDA00035011209800000216
Energy E ofI′、EQ′N modal components of the in-phase component I ' (t) and the quadrature component Q ' (t) of the signal x ' (t) are calculated by the same method as in step S5
Figure FDA00035011209800000217
Of box dimension Yi I′、Yi Q′And reconstructing the in-phase component
Figure FDA00035011209800000218
Reconstructing the orthogonal components
Figure FDA00035011209800000219
Of (2) box dimension
Figure FDA00035011209800000220
Construct the obtainerEigenvector
Figure FDA00035011209800000221
And inputting the data into the radiation source recognition model trained in the step S7 to obtain a radiation source recognition result.
CN202011144744.2A 2020-10-23 2020-10-23 Radiation source identification method based on VMD-Teager energy operator Active CN112287796B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011144744.2A CN112287796B (en) 2020-10-23 2020-10-23 Radiation source identification method based on VMD-Teager energy operator

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011144744.2A CN112287796B (en) 2020-10-23 2020-10-23 Radiation source identification method based on VMD-Teager energy operator

Publications (2)

Publication Number Publication Date
CN112287796A CN112287796A (en) 2021-01-29
CN112287796B true CN112287796B (en) 2022-03-25

Family

ID=74425029

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011144744.2A Active CN112287796B (en) 2020-10-23 2020-10-23 Radiation source identification method based on VMD-Teager energy operator

Country Status (1)

Country Link
CN (1) CN112287796B (en)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114757236B (en) * 2022-06-13 2022-10-11 广东工业大学 Electroencephalogram signal denoising optimization method and system based on TQWT and SVMD
CN115169406B (en) * 2022-07-15 2023-04-07 中国人民解放军国防科技大学 Instantaneous phase fingerprint feature enhancement method based on empirical mode decomposition
CN115510924B (en) * 2022-11-17 2023-04-07 中铁第一勘察设计院集团有限公司 Radio frequency fingerprint identification method based on improved variational modal decomposition

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108090462A (en) * 2017-12-29 2018-05-29 哈尔滨工业大学 A kind of Emitter Fingerprint feature extracting method based on box counting dimension
CN109164489A (en) * 2018-10-15 2019-01-08 西南石油大学 A kind of earthquake fluid prediction technique based on VMD Yu TK energy operator
CN109307862A (en) * 2018-07-05 2019-02-05 西安电子科技大学 A kind of target radiation source individual discrimination method
CN110147848A (en) * 2019-05-24 2019-08-20 哈尔滨工业大学 A kind of radiation source personal feature Enhancement Method based on time-variable filtering theory
CN110163040A (en) * 2018-03-28 2019-08-23 信阳农林学院 Radar emitter signal identification technology in non-gaussian clutter
CN110197209A (en) * 2019-05-15 2019-09-03 电子科技大学 A kind of Emitter Recognition based on multi-feature fusion
CN110542828A (en) * 2019-07-17 2019-12-06 佛山电力设计院有限公司 Distribution network fault positioning method and system

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP1920359A2 (en) * 2005-09-01 2008-05-14 Astragroup AS Post-recording data analysis and retrieval
CN106772417B (en) * 2016-12-31 2017-11-14 华中科技大学 A kind of multiple dimensioned infrared spectrum pattern measurement method and system of moving-target various dimensions

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108090462A (en) * 2017-12-29 2018-05-29 哈尔滨工业大学 A kind of Emitter Fingerprint feature extracting method based on box counting dimension
CN110163040A (en) * 2018-03-28 2019-08-23 信阳农林学院 Radar emitter signal identification technology in non-gaussian clutter
CN109307862A (en) * 2018-07-05 2019-02-05 西安电子科技大学 A kind of target radiation source individual discrimination method
CN109164489A (en) * 2018-10-15 2019-01-08 西南石油大学 A kind of earthquake fluid prediction technique based on VMD Yu TK energy operator
CN110197209A (en) * 2019-05-15 2019-09-03 电子科技大学 A kind of Emitter Recognition based on multi-feature fusion
CN110147848A (en) * 2019-05-24 2019-08-20 哈尔滨工业大学 A kind of radiation source personal feature Enhancement Method based on time-variable filtering theory
CN110542828A (en) * 2019-07-17 2019-12-06 佛山电力设计院有限公司 Distribution network fault positioning method and system

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
Novel Method for Identifying Fault Location of Mixed Lines;Lei Wang et al.;《www.mdpi.com/journal/energies》;20180612;1-19 *
基于VMD和SVM的舰船辐射噪声特征提取及分类识别;李余兴 等;《国防科技大学学报》;20190228;第41卷(第1期);89-94 *
通信辐射源个体识别技术研究;李秋雪;《万方数据》;20211201;1-67 *

Also Published As

Publication number Publication date
CN112287796A (en) 2021-01-29

Similar Documents

Publication Publication Date Title
CN112287796B (en) Radiation source identification method based on VMD-Teager energy operator
CN108845250B (en) On-load tap-changer fault identification method based on vibration signal feature extraction
Wang et al. Fault diagnosis of diesel engine based on adaptive wavelet packets and EEMD-fractal dimension
Yan et al. Improved Hilbert–Huang transform based weak signal detection methodology and its application on incipient fault diagnosis and ECG signal analysis
CN109307862A (en) A kind of target radiation source individual discrimination method
CN103220241B (en) Signal box counting dimension feature extracting method under a kind of Low SNR
CN107392123B (en) Radio frequency fingerprint feature extraction and identification method based on coherent accumulation noise elimination
CN106533577B (en) A kind of non-Gaussian noise suppressing method based on energy measuring
CN111189638B (en) HMM and QPSO optimization algorithm-based bearing fault degree identification method
CN112101245A (en) Short-time Fourier transform mechanical impact feature extraction method based on frequency domain window function
CN108009584A (en) Deficient based on the detection of single source point determines blind source separation method
CN109506938A (en) Based on the failure extracting method for improving wavelet packet threshold denoising and local mean value decomposition
CN111160317B (en) Weak signal blind extraction method
CN107196661A (en) The non-homogeneous lack sampling blind reconstructing method of harmonic combinations signal based on sparse theory
CN105578480A (en) Undersampling frequency spectrum perception pre-decision method orienting broadband modulation converter
CN113325277A (en) Partial discharge processing method
CN106534014A (en) Accurate detection and separation method for multi-component LFM signal
CN106209703A (en) A kind of Frequency Hopping Signal Blind Parameter Estimation and device
CN110909480A (en) Denoising method and denoising device for water turbine vibration signal
CN112666084A (en) SF based on Raman spectrum6Noise reduction method for decomposition characteristic component detection
CN103245518A (en) Rotary mechanical fault diagnosis method based on differential local mean decomposition
CN107340133A (en) A kind of bearing condition monitoring method based on fitting Lifting Wavelet and higher order cumulants analysis
CN106301755B (en) A kind of noise-reduction method and system of the energy leakage signal based on wavelet analysis
CN109165545A (en) A kind of new signal characteristic extracting methods based on fractal dimension
CN112689288A (en) Radio frequency fingerprint extraction and identification method based on WANN

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant