CFAR based NOMP for Line Spectral Estimation and Detection

The stopping criterion proposed in [19] assumes that the noise variance is known and constant. This stopping criterion sets the threshold accurately to guarantee a specified probability of false alarm. In practice, the noise variance is unknown and can often be variable. To provide predictable detection and false alarm behaviour in realistic scenarios, CFAR detection or “adaptive threshold detection ” is developed.

The idea of CFAR detector is to estimate the noise variance from the data in real time, so that the detection threshold can be adjusted to maintain the desired $\bar{\rm P}_{\rm FA}$. Below we present the details of the CFAR design.

We consider a binary hypothesis testing problem

and we choose between the null hypothesis ${\mathcal{H}}_{0}$ and the alternative hypothesis ${\mathcal{H}}_{1}$, where $\mathcal{\varepsilon}$ is the AWGN with its elements $\mathcal{\varepsilon}_{\mathbf{n}}$ being $\mathcal{\varepsilon}_{\mathbf{n}}\sim{\mathcal{CN}}({0},\sigma^{2})$. With the knowledge of noise variance $\sigma^{2}$, a typical detector deciding ${\mathcal{H}}_{1}$ is

where $\tilde{\mathcal{Y}}$ denotes the normalized $D$ dimensional DFT of ${\mathcal{Y}}$ defined as

${\tilde{\mathbf{n}}}_{\rm peak}$ is the peak localization of $|\tilde{\mathcal{Y}}|^{2}$ given by

where $\tilde{\mathcal{N}}=\{\tilde{\mathbf{n}}\in{\mathbb{N}}^{D},0\leq k_{d}\leq N_{d}-1,d=1,2,\cdots,D\}$. In other words, the detector (11) decides that the signal is present if the peak value of the spectrum exceeds a threshold. It is also worth noting that (11) can be evaluated efficiently through $D$ dimensional DFT. Given the false alarm probability $\bar{\rm P}_{\rm FA}$, the threshold multiplier $\alpha^{\prime}$ can be obtained by [19]

Motivated by the conventional CFAR approach, the noise variance $\sigma^{2}$ is estimated with the data samples as

where ${\tilde{\mathbf{n}}}_{\rm peak}$ (13) is the index of the cell under test (CUT) $\tilde{\mathcal{Y}}_{{\tilde{\mathbf{n}}}_{\rm peak}}$, ${\mathcal{T}}_{{\tilde{\mathbf{n}}}_{\rm peak}}$ is the index set of the reference cells, $N_{r}$ denotes the number of reference cells and $N_{r}=|{\mathcal{T}}_{{\tilde{\mathbf{n}}}_{\rm peak}}|$. The reference cells are averaged to estimate the noise variance. In practice, we set guard cells between the reference cells and the CUT. The reason is that a frequency not on the DFT grid exactly might straddle frequency cells. In this case, the energy in the cell adjacent to $\tilde{\mathcal{Y}}_{{\tilde{\mathbf{n}}}_{\rm peak}}$ would contain both noise and signal energy. The extra energy from the signal would tend to raise the estimate of the noise variance, resulting in a higher threshold and a lower $\bar{\rm P}_{\rm FA}$ and $\bar{\rm P}_{\rm D}$ than intended. The details of selecting reference cells and guard cells can be referred to [20]. While in multiple target scenario, the construction of reference cells from [20] can not be directly borrowed. The reason is that there may some other signals lie on the reference cells. Once NOMP-CFAR estimates these signals, eliminates the effects of these signals, and detects a new signal, directly constructing the reference cells via the CFAR approach tends to make the estimate of the noise variance lower, resulting in a lower threshold and a higher $\bar{\rm P}_{\rm FA}$ than intended, see Subection III-B for further details. Consequently, the cells nearby the other targets should be excluded from the reference cells.

With $\hat{\sigma}^{2}$ (15), the detector in the case of unknown noise variance is

Given this detector $T({{\mathcal{Y}}})$ (16), its false alarm probability $\bar{\rm P}_{\rm FA}$ and detection probability $\bar{\rm P}_{\rm D}$ are analyzed. Firstly, we need to calculate the probability density function (PDF) of $\hat{\sigma}^{2}$ (15). Note that the normalized DFT is a unitary matrix, and the DFT of ${\mathcal{Y}}$ under the null hypothesis ${\mathcal{H}}_{0}$ is still the AWGN, i.e.,

In general, the exact PDF of $\hat{\sigma}^{2}$ (15) is hard to obtain as the reference cells depend on the peak localization ${\tilde{\mathbf{n}}}_{\rm peak}$ of the spectrum. Provided that the reference cells are not too near the peak, the PDF of $\tilde{\mathcal{Y}}_{\tilde{\mathbf{n}}}$, ${\tilde{\mathbf{n}}}\in{\mathcal{T}}_{{\tilde{\mathbf{n}}}_{\rm peak}}$ is approximated as a Gaussian distribution, i.e.,

Based on (18), we proceed to calculate the average false alarm $\bar{\rm P}_{\rm FA}$ and the result is summarized as Proposition 1.

Note that we have provided two formulas (19) and (20) for $\bar{\rm P}_{\rm FA}$ with the required threshold multiplier $\alpha$. Eq. (19) is simple for calculating via numerical integration. In contrast, eq. (20) is hard to calculate as the range of the terms in the summation can be very large. However, eq. (20) indeed provides some insight into the relationship between the establishing result and the previous results, as shown in the following Remarks.

Define

where $\hat{T}$ is (21). Note that $\Delta=0$ implies that the integrated amplitude of the detected signal just cross the corresponding threshold. For $\Delta>0$, $\Delta$ is the excess in which the amplitude of the detected signal above the corresponding threshold. While for $\Delta<0$, $-\Delta$ is the excess in which the amplitude of the detected signal below the corresponding threshold. For the traditional CFAR detector, it makes a binary decision based on the sign of $\Delta$. In contrast, we provide a “soft” CFAR detector which outputs $\Delta$. The “soft” CFAR detector is summarized in Algorithm 1.

It is also meaningful to obtain the detection probability $\bar{\rm P}_{\rm D}$ for the required threshold multiplier $\alpha$. Under the alternative hypothesis ${\mathcal{H}}_{1}$, we declare that a frequency $\bm{\omega}$ exists if the spectrum $\tilde{\mathcal{Y}}$ at the DFT frequency ${\bm{\omega}}_{\hat{\mathbf{n}}}$ exceeds its corresponding threshold, where the $d$th element of $\hat{\mathbf{n}}$ is

where $\bar{\omega}_{n_{d}}\triangleq(n_{d}-1)2\pi/N_{d}$ denotes the $d$th dimensional the $n_{d}$th DFT frequency grid. Straightforward calculation shows that under the alternative hypothesis ${\mathcal{H}}_{1}$,

$|\omega_{d}-\bar{\omega}_{\hat{n}_{d}}|\leq\pi/N_{d}$ and $|\beta_{\hat{\mathbf{n}}}|\leq 1$. It can be known that $\tilde{\mathcal{Y}}_{\hat{\mathbf{n}}}\sim\mathcal{CN}(\sqrt{N}x\beta_{\hat{\mathbf{n}}},\sigma^{2})$. Conditioned on $\hat{T}$, the detection probability ${\rm P}_{\rm D}(\hat{T})$ is

where $Q_{1}(\cdot,\cdot)$ is the Marcum Q-function. The average detection probability can be calculated as

Substituting (22) and (III-A) into (41), one obtains

where the $\stackrel{{\scriptstyle a}}{{=}}$ is due to the substitution $T=\frac{\hat{T}}{\sigma^{2}}$, $\stackrel{{\scriptstyle b}}{{\approx}}$ is due to ${\rm E}[\beta_{\hat{\mathbf{n}}}]=0.88^{D}$ by assuming a uniform distribution for a frequency within a DFT grid interval, see [19] for further details.

The average detection probability $\bar{\rm P}_{\rm D}$ in the single target scene is (III-A). In the following, we discuss the detection probability $\bar{\rm P}_{\rm all}$ for multiple targets scene. Let ${\rm E}_{k}$ denote the event of the $k$th target being detected, and let ${\rm P}\left({\rm E}_{k}\right)$ denote the probability of the $k$th target being detected. As a result, the probability of all targets being detected can be calculated as

where ${\rm P}\left({\rm E}_{i}|{\rm E}_{1}\cdots{\rm E}_{i-1}\right)$ denote the probability of the $i$th target being detected conditioned on all the previous $i-1$ detected targets. Obviously, ${\rm P}\left({\rm E}_{i}|{\rm E}_{1}\cdots{\rm E}_{i-1}\right)\leq\tilde{\rm P}_{{\rm D}}({\rm SNR}_{i})$. Consequently, ${\rm P}\left({\rm E}_{1}{\rm E}_{2}\cdots{\rm E}_{K}\right)$ can be upper bounded as

In Subsection IV-B, $\bar{\rm P}_{\rm all}^{\rm upper}$ is evaluated to provide an upper bound of the detection probability of all the targets.

Subsection III-A has proposed the CFAR detector for a single target scenario and provided (19) for calculating the threshold multiplier $\alpha$ for a specified $\bar{\rm P}_{\rm FA}$. Now we focus on the multiple target detection problem and design NOMP-CFAR algorithm.

For the target detection and estimation problem, one should try to find a minimal representation of the line spectral such that the residue can be well approximated as the AWGN. Note that the conventional CFAR based approach consists of two steps: Firstly, performing the fast Fourier transform (FFT) on $\mathcal{Y}$ to obtain the spectrum $\tilde{\mathcal{Y}}$ ¹¹1For range estimation, range Doppler estimation, the spectrum corresponds to range spectrum and range-Doppler spectrum, respectively.. Secondly, implementing the CFAR detection such as CA or OS methods for all the cells by specifying the guard length, training length and the average false alarm rate $\bar{\rm P}_{\rm FA}$, the target detection results are obtained. Obviously, the FFT based approach is suboptimal as it suffers from grid mismatch and inter-target intergerence. Besides, the resolution of the FFT-based approach is limited by the Rayleigh criterion. Therefore, we propose NOMP-CFAR to overcome the drawbacks of the traditional CFAR based approach, while inherit the advantages of the CFAR behaviour.

The main steps of the NOMP-CFAR can be divided into two steps: Initialization and MOE or Detection step. The initialization step provides a good initial point of the NOMP-CFAR. While for the MOE step, it is iterative and it includes an activation step and deactivation step, i.e., by activating the target or deactivating the target. For the initialization step: It is assumed that the number of targets $K$ is upper bounded by a known constant $K_{\rm max}$, i.e., $K\leq K_{\rm max}$. This is reasonable as one could choose $K_{\rm max}=N$ as the number of targets must be less than the number of measurements $N$. Then we run NOMP consisting of coarse detection, single refinement, cyclic refinement, least squares (LS) step to obtain the candidate amplitudes and frequencies set ${\mathcal{K}}_{K_{\rm max}}=\{(\hat{x}_{k},\hat{{\bm{\omega}}}_{k}),k=1,2,\cdots,K_{\rm max}\}$. Note that this step is necessary as it alleviates the intertarget interference and provides fine results for the MOE step. For further details please see [19]. For the MOE step, the CFAR approach is implemented to perform target detection, i.e., activate and deactivate the components of the frequencies. In addition, all the steps such as coarse detection, single refinement, cyclic refinement, LS in NOMP are also incorporated to enhance the target estimation and detection accuracy. Below we present the details of the MOE step.

Given the candidate set ${\mathcal{K}}=\{(\hat{x}_{k},\hat{{\bm{\omega}}}_{k}),k=1,2,\cdots,\hat{K}\}$ and $\hat{K}=|{\mathcal{K}}|\leq K_{\rm max}$, one could use the CFAR criterion to detect the $k$th frequency. To implement the CFAR detector, we proceed as follows: The residue of the observation $\mathcal{Y}_{\rm r}({\mathcal{K}})$ after eliminating all the targets ${\mathcal{K}}$ is

Provided that all the amplitudes and the frequencies are detected and estimated in high accuracy, $\mathcal{Y}_{\rm r}({\mathcal{K}})$ is approximately AWGN, i.e.,

and ${\rm vec}(\varepsilon_{\rm r})\sim{\mathcal{CN}}({\mathbf{0}},\sigma_{\rm r}^{2}{\mathbf{I}}_{N})$. To detect the $k$ target, we obtain the pseudo measurement $\mathcal{Y}_{{\rm r}}({\mathcal{K}}\backslash(\hat{x}_{k},\hat{{\bm{\omega}}}_{k}))$ by adding the contribution of the $k$th frequency to the residue $\mathcal{Y}_{\rm r}({\mathcal{K}})$, i.e.,

where $\stackrel{{\scriptstyle a}}{{\approx}}$ is due to (46). We now input $\mathcal{Y}_{{\rm r}}({\mathcal{K}}\backslash(\hat{x}_{k},\hat{{\bm{\omega}}}_{k}))$ and the required threshold multiplier $\alpha$ to the CFAR detector to provide a soft detection of the $k$th frequency, $k=1,2,\cdots,\hat{K}$, i.e., calculate $\Delta_{k}$ (37). Note that $\Delta_{k}\geq 0$ means that the power of the CUT exceeds the estimated threshold, and this target should be detected if the CFAR is adopted. If we radically implement the CFAR detector to detect all the targets and then stop, miss event may happen in certain scenarios. For example, in the Initialization step, a strong target may be detected as two closely-spaced frequencies with distinguishable amplitudes. This splitting phenomenon degrades the SNR of the frequency and the CFAR detector may miss the two targets at the same time, which misses the true frequency. To overcome this problem, we only activate or deactivate one frequency in each iteration. In this way, even after splitting, the frequency corresponding to the weaker amplitude may be eliminated and the frequency one with stronger amplitude can be refined, and the true frequency may be detected. To activate or deactivate one frequency in each iteration, the greedy approach is adopted by calculating

and ${\Delta}_{\hat{k}}$. Depending on the sign of ${\Delta}_{\hat{k}}$, we choose two different actions:

• A1

  ${\Delta}_{\hat{k}}<0$: In this setting, deactive the $\hat{k}$th frequency by removing the $\hat{k}$th frequency, i.e., updating $\hat{\mathcal{K}}^{\prime}={\mathcal{K}}\backslash(\hat{x}_{\hat{k}},\hat{{\bm{\omega}}}_{\hat{k}})$. Perform the single refinement, cyclic refinement and LS for the remaining targets and obtain $\hat{\mathcal{K}}^{\prime\prime}$. This refinement step is necessary to ensure the high estimation accuracy. Then for the updated parameter set $\hat{\mathcal{K}}^{\prime\prime}$, perform CFAR detection as done before.

• A2

  ${\Delta}_{\hat{k}}\geq 0$: In this setting, the frequency set is kept. Then, perform CFAR detection on the residual $\mathcal{Y}_{\rm r}({\mathcal{K}})$ to determine whether a new amplitude frequency pair should be added on the list ${\mathcal{K}}$ or not. If the CFAR does not detect any target or the number of targets is greater than $K_{\rm max}$, then the algorithm stops. Otherwise, add the new amplitude frequency pair into the list and iterate.

Note that compared to the traditional CFAR, our procedure can be viewed as a soft CFAR as we proceed very gently. The NOMP-CFAR is summarized as Algorithm 2.

The above NOMP-CFAR can be easily extended to deal with the MMV setting and the compressive setting.