Parameter estimation of frequency shift keying radar signal intercepted by Nyquist folding receiver using periodic linear frequency modulation local oscillator

Frequency shift keying (FSK) signal has been widely used in modern radar systems. To intercept the FSK radar signal with a high interception probability, radar reconnaissance receiver should have an ultra ‐ wideband receiving ability. The Nyquist folding receiver (NYFR) is a novel non ‐ cooperative ultra ‐ wideband receiver and it requires a small amount of equipment. When the FSK radar signal is intercepted by the NYFR, the corresponding output is a hybrid modulation signal. Therefore, it is necessary to investigate the parameter estimation of FSK radar signal intercepted by the NYFR. Herein, based on the NYFR using sinusoidal frequency modulation local oscillator (LO), the drawbacks of several existing parameter estimation methods are first analysed. To overcome these drawbacks, the NYFR architecture is improved by changing the LO modulation type to periodic linear frequency modulation. Then, based on the improved NYFR, a novel parameter estimation approach is proposed according to the LO modulation characteristics and the template matching rule. The proposed method is


| INTRODUCTION
With the development of electronic warfare, the frequency space of modern radar systems is becoming wider [1]. To intercept radar signals spread across the frequency domain, radar reconnaissance receiver should have an ultra-wideband receiving ability.
During the past decades, several non-cooperative wideband receiver architectures have been proposed. The channelization structure is a typical wideband receiver and it can achieve real time processing [2]. However, its instantaneous receiving bandwidth is restricted by the sample rate of analogue to digital converter (ADC). To expand the instantaneous receiving bandwidth, some new wideband receiving structures [3,4] have been proposed in recent years. The Nyquist folding receiver (NYFR) is one of the novel wideband receiving architectures and it can break through the restriction of the Nyquist rate [5]. Since the NYFR only uses a local oscillator (LO) and an ADC, it can achieve a high interception probability with a small amount of receiving resources.
In the meantime, because frequency shift keying (FSK) radar signal is easy to generate and its range resolution is only dependent on the code rate, it has been used in many modern radar systems. The FSK radar signal is a typical low probability interception (LPI) radar signal [6] and the common FSK radar signal code schemes include the M-FSK (e.g. 2FSK) code and the Costas code [6][7][8]. Considering the noncooperative wideband receiver, it is necessary to investigate the parameter estimation of FSK radar signal intercepted by the NYFR.
Recently, relevant studies are mainly focused on the NYFR architectures analysis and some parameter estimation methods of radar signals are intercepted by the NYFR. For example, the synchronous NYFR (SNYFR) is proposed to simplify the NYFR output analysis [9]. In addition, the photonic NYFR structure [10] and the NYFR output noise properties [11] are investigated. In terms of processing the radar signal intercepted by the NYFR, the mono-pulse radar signal, the binary phase coded signal, the linear frequency modulation signal and the hybrid modulation radar signal are investigated [12][13][14][15][16]. However, for the FSK radar signal intercepted by the NYFR, there are only a few studies in the literature.
Thus, herein, we will focus on the parameter estimation of FSK radar signal intercepted by the NYFR. Firstly, considering the existing parameter estimation approaches based on the NYFR, the drawbacks of existing methods for the FSK radar signal estimation are analysed mathematically. Then, the LO modulation type of the NYFR is improved to make the receiving structure become more suitable for the FSK signal processing. Furthermore, a parameter estimation method of FSK radar signal intercepted by the NYFR using the improved LO is given, which can estimate the parameters of FSK radar signal with different code schemes. Finally, several simulations are given to show the merits of the proposed method.
Here, Section 2 gives the NYFR architecture and the model of FSK radar signal intercepted by the NYFR. Section 3 investigates the shortcomings of existing methods when they estimate the FSK radar signal intercepted by the NYFR. In Section 4, the LO is improved to overcome the drawbacks of existing methods and the corresponding NYFR output signal model is given. Then, the parameter estimation approach of FSK signal intercepted by the NYFR using the improved LO is proposed in Section 5. Section 6 shows the simulation results and the conclusion is given in Section 7.

| NYFR ARCHITECTURE AND SIGNAL MODEL
The NYFR uses an ADC and an LO to realise wideband receiving [17]. Generally, the LO of the NYFR contains a certain type of modulation and one of the LO modulation types is the sinusoidal frequency modulation (SFM) [5]. The NYFR architecture using SFM LO [5,17] is shown in Figure 1.
In Figure 1, f s is the LO carrier frequency, k is an integer, p SFM ðtÞ is the pulse template and p SFM ðtÞ is the sampling pulse [17]. The NYFR uses zero crossing rising (ZCR) voltage time of the LO to generate the radio frequency (RF) sample clock. Because the LO modulation type is SFM, the RF sample clock can provide a non-uniform sample rate [5]. The SFM modulation of the LO in Figure 1 can be expressed as [5,9] where f sin is the modulation frequency and m f is the modulation coefficient.
In Figure 1, the intercepted analogue signal with high frequency is filtered by the pre-selected filter Hð f Þ firstly, which allows multiple Nyquist zones (NZs) to be sampled. Then, the intercepted signal is sampled by the non-uniform RF pulse-based sampling p SFM ðtÞ and filtered by an interpolation filter H 0 ðf Þ. After the interpolation filtering, the intercepted signal was folded to the intermediate frequency (IF) that is located in the first NZ, which means the folded signal satisfies the Nyquist rate [16,17]. Finally, the interpolation filter output is discretised/digitalised by the ADC and the discrete NYFR output x os ðnÞ is obtained [5]. Because of the non-uniform RF pulse-based sampling and the folding process, the intercepted signal with high frequency can be sampled by an ADC with low sample rate, and the NYFR output will contain the LO modulation part that indicates the high frequency of the intercepted signal [5].
To simplify the derivation here, the signal is assumed as a complex signal and the FSK radar signal can be written as where A is the signal amplitude, t ∈ ½0; T �, T is the signal duration, f i ¼ f 1 ; f 2 ; ⋯; f C are the carrier frequencies in each code, f i > 0, f i > f s , C means the number of different carrier frequencies, N c is the total number of codes in the code scheme, T c ¼ T =N c is the duration of each code, mod means modulo operator, and the code rate is r c ¼ 1=T c . Usually, for the intercepted FSK radar signal, we have N c ≥ C. When the intercepted signal is complex, the NZ bandwidth is equal to f s and the FSK signal spectrum here is assumed in one NZ [16]. Based on (1), when the FSK radar signal in (2) is intercepted by the NYFR, the corresponding NYFR output can be expressed as [5,9].
where wðtÞ is the noise, and k NZ ¼ roundðf i =f s Þ is the NZ index, which indicates the original carrier frequency of the intercepted signal. Compared with the intercepted FSK signal in (2), the NYFR output contains an added SFM part and it is an FSK/SFM hybrid modulation signal. Besides, the FSK signal with high frequency has been folded to the IF (3). According to (3), the F I G U R E 1 NYFR architecture using SFM LO QIU ET AL.
-457 ultra-wideband receiving of the NYFR is at the expense of its output modulation complexity.
Herein, the parameters of the intercepted FSK signal to be estimated are the carrier frequencies f i ; i ¼ 1; 2; ⋯; C and the code rate. In (3), due to the non-cooperative receiver, the NZ index k NZ determined by the signal carrier frequency is also an unknown parameter. If the NZ index k NZ can be estimated first, the NYFR output can be transformed into an FSK signal by LO demodulation, and the parameter estimation can be simplified. Thus, the NZ index estimation is a critical step for the parameter estimation of the FSK signal intercepted by the NYFR. In the next section, several existing NZ index estimation methods will be analysed first.

| ANALYSIS OF EXISTING METHODS
Because the unknown parameter in the added SFM component in (3) is the NZ index, estimating the NZ index firstly can simplify the hybrid modulation signal processing. Presently, the main existing methods of NZ index estimation are the multiple channel method [12] and the matching component function (MCF) [13,18]. In addition, the local high-order phase function (LHPF) [19] has the optimal transformation kernel for SFM component, and it can be used to estimate the signal that contains SFM part as well. Thus, focusing on the FSK signal intercepted by the NYFR using SFM LO, the three methods are investigated.

| Multiple channel
The NZ index estimation using the multiple channel method is illustrated in Figure 2.
In Figure 2, the multiple channel method constructs M channels to demodulate the SFM part in the NYFR output, where M means the number of the monitored NZs [16]. The group of SFM demodulation signals can be expressed as s SFM ðt; kÞ ¼ e jπkm f sinð2πf sin tÞ ð4Þ where the NZ searching index k ¼ 0; 1; ⋯; M − 1. Based on (3) and (4), the estimated NZ index can be written as [12] k NZ ¼ arg where ℑð⋅Þ represents the Fourier transform operator. In (5), the essence of the multiple channel method is using the maximum peak value of the demodulated signal spectra with different NZ searching indexes to estimate the NZ index in (3). However, unlike the mono-pulse radar signal, the FSK radar signal (e.g. Costas FSK signal) in frequency domain has multiple peaks and it is difficult to identify the maximum peak with different SFM demodulation signals under noisy condition.
The Costas FSK signal is used as an example. The signal amplitude is 1, the signal duration is 1.2 μs, the Costas carrier frequency array is [4.16, 4.14, 4.22, 4.18, 4.20, 4.12] GHz, the number of different carrier frequencies C is 6 and the total number of codes N c is 6. The NYFR LO modulation type is SFM, the carrier frequency of the LO is 1 GHz, the LO modulation frequency is 10 MHz, the LO modulation coefficient is 1, and the ADC sampling frequency is 2 GHz. Thus, the NZ index is roundð4:16 GHz=1 GHzÞ ¼ 4. The output signal-to-noise ratio (SNR) is 5 dB. Figure 3 gives the demodulated signal with different NZ searching indexes k in frequency domain. , there are several peaks in the frequency domain and the corresponding maximum peak values are almost the same, which means the difference between the peak values of the spectra with different NZ searching indexes is small. Therefore, the multiple channel method is hard to obtain the NZ index of Costas FSK signal intercepted by the NYFR under noisy conditions.

| Matching component function
In this subsection, the MCF is analysed. The MCF can estimate the NZ indexes directly for several types of radar signals intercepted by the NYFR [13,18]. We assume the ADC sampling frequency is f ADC , the sampling interval is T ADC ¼ 1=f ADC , the number of points in one code is N 1 ¼ T c =T ADC , and the total number of signal points is N ¼ T =T ADC . Hence, the discrete expression of (3) can be written as The MCF definition is [16].
Therefore, according to (7), (9) and the Jacobi-anger expansion, the MCF of FSK signal intercepted by the NYFR can be calculated as where J m i ð⋅Þ; i ¼ 1; 2 are the Bessel functions with m i order. The third summation term in (10) is assumed as X and it is written as To further analyse (12), we use the FSK signal segment that contains two carrier frequencies. After instantaneous autocorrelation, the corresponding result in (12) whose length is N 1 can be written as Focusing on the second term in (13), it can be calculated as In (14), the main contribution of the MCF peak value comes from the non-exponential part. Considering f cp ≠ f cq , the value of the non-exponential part in (14) is possibly small or negative. Therefore, according to (14), the value of X 1 may be small because the carrier frequencies of the intercepted FSK (a) signal are unknown. Moreover, the value of X will be small, which leads to a decline in the value of Pðk NZ Þ. Consequently, the MCF will lose its NZ index estimation ability.
The simulations of the above MCF analysis are given and two 2FSK signals are chosen for an instance. The durations and the amplitudes of both signals are 1 μs and 1, respectively. The code schemes of the two 2FSK signal are the same and they are [1, 0, 0, 1, 0, 1, 1, 0, 1, 0]. For the first 2FSK signal, its carrier frequencies are 4.16 and 4.18 GHz. Hence, the number of different carrier frequencies C is 2 and the total number of codes N c is 10. The carrier frequencies of the second 2FSK signal are 4.15 and 4.18 GHz. The NYFR parameters are the same as those in Figure 3, the output SNR is 5 dB, and τ M ¼ 100. Figure 4 shows the MCFs of the two 2FSK signals intercepted by the NYFR.
In Figure 4 (a), the MCF shows a peak when the NZ searching index k ¼ 4, which means the MCF can estimate the NZ index of the first intercepted 2FSK signal. However, the MCF in Figure 4 (b) shows no peak when the NZ searching index k ¼ 4. Noticing the carrier frequencies of the second 2FSK signal and (14), the non-exponential part in (14) is sinð3π=2Þ=sinð0:03π=2Þ and it is negative. Thus, the value of X in (12) is small and the MCF loses its NZ index estimation ability.
Based on the above simulations, the MCF can only estimate the NZ index when the value of X in (12) is large. Because the carrier frequencies of the intercepted FSK signal are unknown, the MCF is not the most suitable parameter estimation method for the FSK signal intercepted by the NYFR.

| Local high-order phase function
Since the NYFR output in (3) contains the SFM component and the LHPF can analyse the SFM [19], the LHPF is another existing method to process the signal intercepted by the NYFR using SFM LO.
For the discrete NYFR output in (6), we consider two FSK codes containing frequency hopping and assume n ¼ nT ADC to simplify the derivation. The FSK signal containing two codes can be written as The LHPF of (15) can be computed as where τ represents the delay and we have The Taylor expansion is Substituting (18) into (17), it yields Based on (19,16) can be calculated as Thus, the time-frequency distribution (TFD) of (20) is Focusing on (21), it contains e 2jπðf cp −f cq Þτ term that comes from the frequency hopping. Therefore, jH L ðn; ψÞj cannot completely accumulate the FSK/SFM signal energy along τ 2 . Besides, the term e 2jπðf cp −f cq Þτ in (21) also means the information of FSK carrier frequency difference is lost and (21) is unable to show the different carrier frequencies of FSK signal in the TFD.
The second 2FSK signal in Figure 4 is used as an example, and the NYFR parameters are the same as those in Figure 4. The output SNR is 10 dB. Figure 5 gives the TFD of the NYFR output calculated by the LHPF.
In Figure 5, the TFD shows the SFM component of the 2FSK/SFM signal. Each red circle in Figure 5 represents the time at which frequency hopping is supposed to occur. However, it can be observed that the frequency hopping time in each red circle is difficult to be identified, and the reason is that the LHPF cannot completely accumulate the energy of FSK signal. When the SNR is lower, the frequency hopping time will be more difficult to find. In addition, as shown in Figure 5, the difference between the carrier frequencies is disappeared in the TFD, which proves the analysis of (21). As a result, the LHPF method is unable to give an accurate estimation result and it is not a suitable method for the FSK radar signal intercepted by the NYFR using SFM LO. Furthermore, the computational complexity of the LHPF is generally large due to the two dimensional searching.
In this section, based on the NYFR using SFM LO, three existing NZ index estimation methods are analysed mathematically. However, each approach has its drawback for different types of FSK radar signal code schemes. Thus, it is hard to obtain the NZ index of unknown intercepted FSK radar signal by using the present approaches. The simulations in this section also prove the corresponding analysis. Additionally, Table 1 summarizes the advantages [12,18,19] and the disadvantages of the three existing methods when they estimate the FSK radar signal intercepted by the NYFR using SFM LO.
In the next section, we will find another way to estimate the FSK radar signal intercepted by the NYFR.

| NYFR USING PERIODIC LINEAR FREQUENCY MODULATION LOCAL OSCILLATOR AND SIGNAL MODEL
When the FSK radar signal is intercepted by the NYFR using SFM LO, the NYFR output is an FSK/SFM hybrid modulation signal. Based on the analysis in Section 3, the existing methods are difficult to process the SFM and the FSK modulation jointly. Because the LO modulation type is known, a feasible way to solve the problem of FSK/SFM signal parameter estimation is to change the NYFR output modulation type, which means the LO modulation type can be changed to simplify the NYFR output processing.
Considering the LO ZCR signal generation [20,21], the periodic linear frequency modulation (PLFM) can be used as the LO modulation type. For the NYFR architecture in Figure 1, the LO modulation is replaced by PLFM and it can be expressed as where μ 0 is the LO chirp rate, T LO is the period duration of the PLFM signal, and other parameters of the NYFR are the same as those in Figure 1. The LO still uses the ZCR voltage time ZCRðsin½2πf s t þ θ PLFM ðtÞ�Þ to generate the non-uniform RF sample clock and the sampling pulse can be written as According to the NYFR architecture and the FSK radar signal in (2), when the LO modulation is the PLFM in (22), the corresponding NYFR output can be written as [5,9].
where wðtÞ is the output noise and other parameters are the same as those in (3) and (22). The IF carrier frequency is still assumed as (23) can be rewritten as In (24), it can be seen that the NYFR output is an FSK/ PLFM hybrid modulation signal. The parameters to be estimated are the FSK code rate, the carrier frequencies and the NZ index.

| PARAMETER ESTIMATION APPROACH
In this section, focusing on the FSK radar signal intercepted by the NYFR using PLFM LO and based on the PLFM characteristics, a novel parameter estimation method is given.

| NZ index estimation
Comparing the FSK signal in (2), the unknown parameter of the added PLFM part in (24) is the NZ index. If the NZ index is obtained firstly, the FSK/PLFM hybrid modulation signal in (24) can be transformed into an FSK signal and the following parameter estimation will be simplified. Thus, the NZ index is estimated first in this subsection.
According to the known LO parameters, a group of PLFM signals can be constructed and they can be written as Then, the instantaneous autocorrelation of (26) in time domain can be computed as r l ðt; kÞ ¼ y l ðt; kÞ y l * ðt þ τ l ; kÞ where the delay τ l is small, τ l < T LO and τ l < T c . The Fourier transform of the instantaneous autocorrelation in (27) can be calculated as Noticing the result in (28), two situations should be considered, and they are.
Situation 1: T c ≥ T LO Situation 2: T c < T LO When T c ≥ T LO (i.e. situation 1), we assume roundðT c =T LO Þ¼ P, and (28) can be computed as.
According to (29), we have where SaðxÞ ¼ sin x=x. Moreover, when k ¼ k NZ and f ¼ 0, (30) can be computed as In (31), based on the characteristics of Sað⋅Þfunction, several terms of jR l ð0; k NZ Þj show peaks at f ¼ 0. Furthermore, according to (30), if k ≠ k NZ , jR l ð0; kÞj has no peak at f ¼ 0 and when T c < T LO (i.e. situation 2), one LO period duration is considered to simplify the analysis. Assuming roundðT LO =T c Þ ¼ Q, the corresponding result in (28) can be calculated as The last term in (33) is written as two different terms and they can be expressed as In addition, the explanations of (33) and (34) are demonstrated in Figure 6. Figure 6 shows the instantaneous autocorrelation result in one LO period. When the FSK signal in one LO period contains different carrier frequencies, the corresponding instantaneous autocorrelation result in (33) will contain the specific signal pieces whose frequencies are f ci − f ciþ1 , and they can be regarded as the cross term pieces. In other words, the cross term pieces come from the instantaneous autocorrelation of different carrier frequencies. Besides, the duration of the cross term piece is τ l , which can be seen in Figure 6.
Because the FSK code scheme and code duration are unknown, the end of instantaneous autocorrelation result in one LO period may contain the cross term piece. In detail, Figure 6 (a) shows the end of instantaneous autocorrelation result has no cross term piece. Hence, the corresponding Fourier transform of the instantaneous autocorrelation in one LO period duration can be expressed by (33) and R l11 ðf ; kÞ in (34). Meanwhile, Figure 6 (b) demonstrates that the end of instantaneous autocorrelation result contains the cross term piece whose frequency is f ci − f ciþ1 and the length of the last cross term piece in one LO period is T LO − ðQT c − τ l Þ. Thus, the Fourier transform of the instantaneous autocorrelation in one LO period duration can be expressed by (33) and R l12 ðf ; kÞ in (34). According to (33) and (34), we have when k ¼ k NZ and f ¼ 0, (35) can be calculated as In (36), several terms of jR l1 ð0; k NZ Þj show peaks at f ¼ 0. In addition, if k ≠ k NZ , jR l1 ð0; k NZ Þj in (35) has no peak at -463 f ¼ 0. Furthermore, the whole signal duration can be proved in the same way as the analysis of one LO period from (33) to (36). Thus, when T c < T LO , jR l ð0; k NZ Þj achieves its maximum as well.
In summary, according to the analysis from (29) to (36), for the FSK radar signal intercepted by the NYFR using PLFM LO, jR l ðf ; kÞj will reach its maximum at f ¼ 0 when k ¼ k NZ . Therefore, the NZ index can be estimated by searching the maximum value of jR l ð0; kÞj and the estimated NZ indexk NZ can be written aŝ Besides, focusing on (30) and (35), when k NZ ≠ k and if we have According to (39), the values of the terms in (30) and (35) that contain Safπτ l ½ðk NZ − kÞμ 0 τ l þ ðf ci − f ciþ1 Þ − f �g will increase, which leads to the growth of jR l ð0; kÞj; k ≠ k NZ . In other words, the difference between the peak value and the non-peak value of jR l ð0; kÞj will be declined. As a result, the NZ index estimation performance will be deteriorated in a noisy environment.

| FSK parameter estimation
Once the NZ index is estimated, a PLFM demodulation signal r PLFM ðtÞ can be constructed based on the known LO parameters and it is expressed as If the NZ index estimation resultk NZ is correct, the NYFR output in (24) can be demodulated by (40) and the LO demodulation result is where w d ðtÞ is the noise after demodulation. The signal in (41) is an FSK signal, which means the NYFR output hybrid modulation signal has been simplified. To estimate the parameter of FSK signal, the approaches such as the TFD and the cyclic spectral analysis [6,22] are commonly used. In this subsection, a fast parameter estimation method is proposed.
Because the signal x FSK ðtÞ contains different frequencies, its Fourier transform result shows several peaks and the IF carrier frequencies f ci can be estimated by finding the peak positions of ℑ½x FSK ðtÞ� in frequency domain. Consequently, the estimated carrier frequencies of the intercepted FSK radar signalf i can be computed aŝ wheref ci are the estimated IF carrier frequencies.
Considering the code rate estimation of FSK signal, the template matching rule (TMR) [23] is used here. The TMR is normally used to obtain the code rate of phase shift keying (PSK) signal, whereas its effectiveness for the code rate estimation of the FSK signal will be proved in rest of this subsection.
In terms of the TMR, we can cut a signal segment from the FSK signal in (41) as a matched template. The matched template duration is T sg and it can be expressed as where t ∈ ½0; T sg �, w sg ðtÞ is the noise, T sg should be small and satisfy T sg < T c . To simplify our derivation, the noise parts in (41) and (43) are omitted. For the FSK signal in (41), the signal section to be matched is where t hop means the hopping frequency time and t hop ∈ ½0; T sg �.
Then, the signal similarity can be written as x sg ðtÞy sg ðtÞdt ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi q ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffiffi According to (43), (44) and (45), the corresponding signal similarity can be computed as when t hop ¼ 0, the signal section in (44) contains no frequency hopping and (46) can be calculated as In (47), if the carrier frequencies of the matched template and the signal segment are different, rð0Þ is a value smaller than 1. Meanwhile, if the carrier frequencies of the matched template and the signal segment are same, rðt hop Þ meets it to the maximum and max½rðt hop Þ� ¼ 1. In addition, if t hop ≠ 0, rðt hop Þ < 1.
When y sg ðtÞ moves along the time axis, the similarity curve for the FSK signal can be defined as Based on the analysis of (46) and (47), if the signal section contains no frequency hopping and the carrier frequencies of the matched template and the signal segment are the same, the value of the similarity curve C s ðtÞ equals 1. If the signal section has no frequency hopping and the carrier frequencies are different, the value of C s ðtÞ is a constant and it is smaller than 1. According to (47), when the value of C s ðtÞ remains unchanged, it represents a carrier frequency of the FSK signal. Furthermore, if the signal section contains frequency hopping, the value of C s ðtÞ will be constantly changed over time. Because T sg is small, the transient duration of C s ðtÞ is small and the middle of the transient duration is the frequency hopping time. Thus, for the FSK signal, the similarity curve can be regarded as a function consisting piecewise constant and the minimum interval of the piecewise constant function represents the FSK code duration.
Moreover, the Haar wavelet transform of the similarity curve can be used to indicate the frequency hopping time [24] and the wavelet transform result is where W T ð⋅Þ is the Haar wavelet transform operator and m is the translation. Each peak position of W C ðmÞ in time domain is the frequency hopping time. Therefore, the estimated code rater c is equal to the maximum peak position of non-directcurrent (non-DC) component of ℑ½W C ðmÞ� in frequency domain [23,24]. Focusing on the FSK radar signal intercepted by the NYFR using PLFM LO, the proposed parameter estimation method steps are given.
Step. 1 Construct a group of PLFM signals in (25) and estimate the NZ indexk NZ based on (37).
Step. 2 Demodulate the PLFM component of the NYFR output and estimate the FSK signal carrier frequencies according to (42).
Step. 3 Based on the TMR, compute the similarity curve in (48). Estimate the code rate by using the Haar wavelet transform in (49) and its non-DC component peak position in frequency domain.
Finally, the computational complexity of the proposed method is considered. For the NZ index estimation, the instantaneous autocorrelation in (27) and the corresponding Fourier transform are calculated. Then, the searching process to find the maximum peak position in (37) is conducted. Thus, the computational complexity is OðNlog 2 NÞ, where N is the total number of signal points. In terms of the FSK signal parameter estimation, the Fourier transform, the similarity curve, the Haar wavelet transform and one dimensional peak searching in frequency domain are conducted to estimate the carrier frequencies and the code rate. Thus, the corresponding computational complexity is OðNN sg Þ, where N sg ¼ T sg f ADC is the number of the points in the matched template duration. As a result, the computational complexity of the proposed method is OðNlog 2 NÞ þ OðNN sg Þ. Compared with the TFD method using two dimensional searching, the computational load of the proposed method is found to be lower.

| SIMULATION ANALYSIS
In this section, several simulations are given to show the merits of the proposed method.

| NZ index estimation
Firstly, the simulations of the proposed NZ index estimation method are given. The Costas FSK radar signal and the 2FSK radar signal are considered. The parameters of the Costas FSK signal are the same as those in Figure 3. For the 2FSK signal, the carrier frequencies are 4.15 and 4.18 GHz, and other parameters are the same as those in Figure 4. In terms of the NYFR parameters, the LO modulation type is PLFM, the LO carrier frequency f s is 1 GHz, the LO chirp rate μ 0 is 400 MHz/μs, the LO period duration T LO is 0.0625 μs, the ADC sampling frequency is 2 GHz and the monitored NZ number M is 10. The noise is white Gaussian noise and the output SNR is 8 dB. Figure 7 and Figure 8 illustrate the NZ index estimation results of the Costas FSK signal and the 2FSK signal, respectively.
Because the LO carrier frequency is 1 GHz, the NZ indexes of the two NYFR outputs are 4. In Figures 7 and 8, each NZ index estimation result achieves its maximum when the NZ searching index k ¼ 4. Obviously, the proposed NZ index estimation method can obtain the NZ indexes of the two common FSK radar signals intercepted by the NYFR using PLFM LO.
In addition, focusing on the parameters of the 2FSK signal and the NYFR in Figure 8, the LO period duration is 0.0625 μs and the 2FSK code duration is 0.1 μs. These parameters coincide with the situation 1 in Section 5.1. Considering the situation 2 in Section 5.1, the LO period duration is changed to prove our analysis. The LO period duration is assumed as 0.125 μs and other parameters remain the same. Figure 9 gives the NZ index estimation result under condition 2.
In Figure 9, when the NZ searching index is 4, the estimation result achieves its maximum, which means the proposed method can obtain the correct NZ index in this situation. Thus, the simulations in Figures 8 and 9 show that the proposed NZ index estimation approach is available for the two situations in Section 5.1. QIU ET AL.

| Code rate estimation based on TMR
Then, several simulations of the code rate estimation based on the TMR method are conducted. The parameters of the FSK signals and the NYFR are the same as those in Figures 7 and 8. The SNR is still 8 dB. According to (41), the NYFR outputs are demodulated by the LO demodulation signals with the estimated NZ indexes in Figures 7 and 8. Figures 10 and 11 give the similarity curves of the demodulated Costas FSK signal and the demodulated 2FSK signal based on the TMR and the corresponding Haar wavelet transform results, respectively.
In Figure 10  is about 0.2 μs and it equals the FSK code duration. Figure 10 (b) shows the Haar wavelet transform of Figure 10 (a). There are five peaks in Figure 10 (b) and the position of each peak represents the frequency hopping time. For the data in Figure 10 (b), its first non-DC component peak position in frequency domain can indicate the code rate estimation result [24].
As shown in Figure 11 (a), the similarity curve of the 2FSK signal shows two values in vertical axis, which coincides with the number of the 2FSK signal carrier frequencies. In addition, the minimum interval of the piecewise constant function in Figure 11 (a) is about 0.1 μs, and it equals the 2FSK code duration. Figure 11 (b) also demonstrates the Haar wavelet transform of Figure 11 (a) and it can be seen that there are seven peaks. Each peak position represents the frequency hopping time of the 2FSK signal. Similarly, the first non-DC component peak position of the Fourier transform of the data in Figure 11 (b) can give the code rate estimation result. Figures 10 and 11 prove the analysis in Section 5.2.

| Parameter estimation performance
At last, the parameter estimation performance is given. In this subsection, three FSK signals are considered. The first FSK signal is the Costas FSK signal in Figure 10, and the second FSK signal is the 2FSK signal in Figure 11. The 4FSK signal [25] is used as the third signal, the signal duration is 1 μs, the signal amplitude is 1, the carrier frequencies of the 4FSK signal are 4.2540 GHz, 4.2902 GHz, 4.2265 GHz, 4.2665 GHz and the code scheme is [0, 1, 1, 2, 3, 3, 0, 1, 2, 1]. The parameters of the NYFR using the PLFM LO are the same as those in Figure 7. Besides, the NYFR using SFM LO is also used as the comparison. For the NYFR using SFM LO, the carrier frequency of the LO is 1 GHz, the LO modulation frequency is 20 MHz, the LO modulation coefficient is 1, and the ADC sampling frequency is 2 GHz.
Considering the NZ index estimation, the multiple channel method in Section 3.1 is used as the comparison. As with the NYFR using PLFM LO, the demodulation signal group of the multiple channel method is changed to the PLFM signal in (25). Meanwhile, the modulation type of the demodulation signal group is still SFM when the LO modulation type is SFM.
The output SNR is from −10 to 10 dB, the SNR interval is 1 dB and 500 Monte-Carlo experiments are conducted for each SNR. The delay τ l is 0.02 μs and the accuracy rate r a is defined as where N total is the total number of Monte Carlo experiments and N correct is the correct number of Monte Carlo experiments. Thus, the NZ index accuracy rates of the three FSK signals using the proposed method and the comparison are given in Figure 12.
In Figure 12, based on the NYFR using PLFM LO, the 2FSK signal and the 4FSK signal accuracy rates of the proposed method and the comparison are almost the same when the SNR is greater than 1 dB, and the accuracy rates of both methods are higher than 90%. When SNR < −5 dB, the 2FSK signal and the 4FSK signal accuracy rates of the comparison are around 60%, while the accuracy rates of the proposed method are about 22%. However, the accuracy rate of the Costas FSK signal using the proposed method is noticeably higher. When SNR > −1 dB, the accuracy rate of the proposed method can achieve 95%, whereas the accuracy rate of the comparison is nearly greater than 90% even when the value of SNR is 10 dB. Because the multiple channel method uses the peak value of the signal spectrum in each channel to estimate the NZ index, its NZ index accuracy rate is worse for the signal (i.e. Costas FSK signal) whose spectrum is spread out. By contrast, the proposed method can overcome this drawback  and it is suitable for the variety of common FSK radar signal code schemes.
Meanwhile, when the NYFR LO modulation type is SFM, the 2FSK accuracy rate is better than the others because the bandwidth of 2FSK is narrowest. The 4FSK accuracy rates of the proposed method and the comparison are almost the same. However, the multiple channel method fails to estimate the NZ index of the Costas FSK signal when the LO modulation type is SFM, which proves the analysis in Section 3.1.
In addition, two 2FSK signals are considered. The parameter settings of the first 2FSK signal and the NYFR are the same as those in Figure 12. The carrier frequencies of the second 2FSK signal are 4.150 and 4.166 GHz, and other parameters of the second 2FSK signal are the same as those of the first 2FSK signal. Figure 13 illustrates the NZ index estimation accuracy rates of the two 2FSK signals.
In Figure 13, the NZ index estimation accuracy rate of the first 2FSK signal is better than that of the second 2FSK signal. When SNR > 0 dB, the accuracy rate of the second 2FSK signal reaches 90%, while the first 2FSK signal can meet the same performance when SNR > −4 dB. Focusing on the carrier frequencies of the second 2FSK signal, the estimation process meets the situation in (39) when the NZ searching index k ¼ 5. Therefore, the NZ index estimation performance is deteriorated. Nevertheless, the proposed method still works.
Furthermore, considering the code rate and the carrier frequency estimation, the root mean squaring errors (RMSEs) are given. The three FSK radar signals in Figure 12 are still used as the intercepted signals and the parameters of the NYFR using PLFM LO are the same as those in Figure 12. The output SNR is from −10 to 10 dB, the SNR interval is 1 dB and 500 Monte-Carlo experiments are conducted for each SNR. Because the TFD method can estimate the code rate [6], the short time Fourier transform (STFT) is chosen as the comparison. The NYFR using SFM LO in Figure 12 is as well used as the comparison. Besides, the proposed FSK signal parameter estimation method is also used to obtain the corresponding code rates and carrier frequencies when the LO modulation is SFM. Figure 14 and Figure 15 demonstrate the RMSEs of the code rate estimation and the carrier frequency estimation, respectively.
In Figure 14, the TMR method has a better code rate estimation performance than the STFT method when the LO modulation type is PLFM. In detail, when SNR > −10 dB, the RMSEs of three FSK signals using the TMR methods are smaller than 1 MHz, whereas the STFT can achieve the same performance when SNR > 0 dB. Because the resolution of the STFT is low, its estimation accuracy is worse when the SNR is high. Besides, since the TMR method only requires one dimensional searching, its computational load is lower. Meanwhile, the TFD method needs two dimensional searching, which brings a heavy computational load. Furthermore, the code rate estimation performances of the 2FSK and 4FSK signals are the same when the LO modulation types are different. However, when the LO modulation type is SFM, the NZ index estimation of the Costas FSK signal in Figure 12 is failed. Thus, the corresponding code rate estimation is unsatisfactory, and the RMSE is greater than 1 MHz even the SNR is higher than 9 dB.
As to the carrier frequency estimation in Figure 15, when SNR > −1 dB and the LO modulation type is PLFM, the RMSEs of three FSK signals are smaller than 1 MHz. Focusing on the NYFR using SFM LO, because the NZ index estimation of the 2FSK signal in Figure 12 is the best, the corresponding carrier frequency estimation performance in Figure 15 is better, and the RMSE is smaller than 1 MHz when SNR > −10 dB. However, since the NZ index estimation of the Costas FSK signal is failed, the corresponding RMSE is still greater than 10 MHz when the SNR is high. Therefore, Figure 15 proves that the NYFR using PLFM LO is more suitable for the parameter estimation of the intercepted FSK radar signal.
In summary, for the FSK radar signal intercepted by the NYFR using PLFM LO, when the SNR is greater than −1 dB, the proposed parameter estimation approach can achieve a good performance. Herein, focusing on the FSK radar signal intercepted by the NYFR, a novel parameter estimation method is proposed. Firstly, based on the NYFR using SFM LO, three existing parameter estimation methods are investigated and the corresponding drawbacks of these approaches for estimating the NYFR output are analysed. Then, to overcome the shortcomings, the NYFR LO modulation type is improved and it is changed to the PLFM. According to the NYFR using PLFM LO, the novel NZ index estimation method is proposed, which is suitable for different common FSK radar signal code schemes. In addition, different signal parameter situations of the NZ index estimation are discussed. Furthermore, the FSK signal code rate estimation based on the TMR method is analysed mathematically. The proposed parameter estimation method only requires one dimensional searching and its computational load is low. Finally, several simulations are conducted to show the merits of the proposed approach. Herein, it is focused on a single FSK radar signal intercepted by the NYFR, and the LO modulation optimization for the multiple radar signals' interception and processing is a future study.