Research on Fast Time–Frequency Reconstruction Algorithm for Wideband Compressive Spectrum Sensing

TL;DR

This paper introduces a new algorithm for efficiently reconstructing time-frequency signals in cognitive radio systems, improving performance even in noisy environments.

Contribution

The novel FTFR algorithm reduces computational complexity for wideband spectrum sensing in cognitive radio.

Findings

The FTFR algorithm uses multi-coset sampling to acquire sub-Nyquist samples for signal reconstruction.

The algorithm effectively restores time-frequency characteristics even in low SNR conditions.

It outperforms existing methods in reconstructing temporal and frequency distributions.

Abstract

Cognitive Radio (CR) is widely acknowledged as a pivotal technology for mitigating the scarcity of spectrum resources, with Transform Domain Communication Systems (TDCSs) regarded as one of the primary candidate technologies for CR. However, conventional Wideband Spectrum Sensing (WBSS) techniques utilized in TDCS exhibit limitations and are insufficient for adapting to the current complex electromagnetic environment. This paper tackles the time–frequency reconstruction challenge in WBSS by proposing a fast time–frequency reconstruction (FTFR) algorithm. The proposed algorithm acquires sub-Nyquist samples through the introduction of a Multi-Coset Sampling structure and reconstructs the autocorrelation of signals across various windows through a series of low-complexity operations. It captures the dynamic variations of signals by integrating spectra from adjacent time windows. In comparison to existing time–frequency reconstruction algorithms in WBSS, the proposed algorithm demonstrates reduced computational complexity. Simulation experiments indicate that the FTFR algorithm can effectively reconstruct the time–frequency characteristics of signals and significantly restore the primary temporal and frequency distributions, even in low Signal-to-Noise Ratio (SNR) environments.

Figures

Full text

1. Introduction

With the advancement of wireless communication technologies, the scarcity of spectrum resources has emerged as a critical constraint limiting the performance improvement of wireless communication systems. Consequently, Cognitive Radio (CR)—a technology capable of dynamically sensing and utilizing spectrum resources—has garnered considerable attention [1]. CR enhances spectral efficiency by adaptively learning from the communication system to sense the spectrum usage in the surrounding environment, fully exploiting spectral holes, and dynamically adjusting its operating frequency to minimize conflicts as much as possible [2]. Transform Domain Communication Systems (TDCS) are regarded as highly promising candidate technologies within the realm of cognitive radio, owing to their flexible spectrum utilization, robust interference resistance, and low probability of being intercepted [3]. To date, the majority of research on spectrum sensing has primarily concentrated on Narrowband Spectrum Sensing (NBSS). However, since TDCS operates over a wide band, NBSS proves inadequate for addressing the inherent complexity and dynamism of wideband signals; thus, it fails to achieve higher throughput or satisfy the demands of contemporary broadband communication systems [4]. Furthermore, TDCS continues to utilize traditional power-spectrum estimation methods for spectrum sensing, including AR spectrum estimation, which fail to account for the increasingly complex variations in the electromagnetic environment. Consequently, research into Wideband Spectrum Sensing (WBSS) that is suitable for TDCS holds substantial practical significance.

Transitioning from Narrowband Spectrum Sensing (NBSS) to Wideband Spectrum Sensing (WBSS), an intuitive approach involves dividing the entire wideband into multiple narrowbands and applying NBSS methods to sequentially sense each sub-band, thereby achieving comprehensive sensing across the broadband. However, this method presents significant drawbacks in practical applications: it results in high latency, fails to meet real-time requirements, and operates inefficiently [5]. An alternative solution, as proposed in [6,7], employs a multi-resolution sensing strategy that optimizes detection accuracy through the adjustment of wavelet-transform basis function precision. This process occurs in two stages: initially, a coarse-grained detection is performed across the entire band to filter out several potentially active narrowbands, followed by fine-grained detection focused exclusively on these promising narrowbands to minimize sensing costs. However, as the spectrum range that requires detection continues to expand, the initial detection may characterize each narrowband as potentially active, necessitating sequential fine-grained detection of every narrowband, which inevitably results in substantial increases in both time and energy consumption.

According to the Nyquist sampling theorem, conventional Wideband Spectrum Sensing (WBSS) methods are constrained by the sampling rate of the front-end Analog-to-Digital Converter (ADC) when the spectrum width that requires detection is substantial. This situation creates a dilemma, if the sampling rate requirements are satisfied, the costs associated with hardware, energy consumption, and scalable deployment may become excessively high, rendering large-scale deployment unfeasible. Conversely, if these costs are controlled, meeting the real-time sensing requirements becomes impractical. In this context, the advent of Compressed Sensing (CS) technology [8] presents a novel solution for WBSS. According to the theory of CS, if a signal exhibits sparse or compressible characteristics, it can be efficiently reconstructed from a limited number of random measurements, thereby achieving low-dimensional representation and a high probability of recovery of the original signal [9,10]. A substantial number of scholars have researched and developed various sub-Nyquist-based spectrum sensing algorithms grounded in CS principles. Based on the nature of the information extracted from sub-Nyquist samples, these methods can primarily be categorized into two categories: CS-based approaches and compressive covariance sensing (CCS)-based approaches.

In CS-based approaches, the essence lies in fully exploiting the sparsity of signals within the frequency domain, acquiring corresponding data samples through sub-Nyquist sampling, and subsequently reconstructing the wideband analog signal using various reconstruction algorithms [11,12]. Although the application of CS technology can indeed reduce the sampling rate and provide certain advantages for WBSS, complex reconstruction algorithms have emerged as a bottleneck in enhancing real-time processing capabilities. To address this impasse, studies [13,14] have developed algorithms that simplify signal reconstruction without necessitating additional samples, thereby providing new directions for future research. Reference [15] utilizes model selection methods to estimate the sparsity of the signal, thereby eliminating reliance on prior knowledge of signal sparsity as configured in earlier research, thus enhancing the algorithm’s adaptability and accuracy. Reference [16] proposes a rapid compressive Wideband Spectrum Sensing method that further reduces computational complexity while preserving high-precision spectrum reconstruction performance through the joint utilization of low-rank and sparse properties. This method operates without the necessity for prior knowledge of sparsity, rendering it compatible with low-complexity platforms for real-time sub-Nyquist rate sensing.

When applying compressive sensing to reconstruct wideband analog signals, several challenges can arise, including the difficulty of achieving a sparse representation of actual signals, high computational complexity, noise interference that degrades the quality of the reconstructed signals and more. To address these challenges, numerous scholars have undertaken a series of in-depth studies. Dyonisius Dony Ariananda was one of the first researchers to propose a wideband power-spectrum estimation algorithm grounded in compressive sampling [17]. This algorithm does not enforce sparsity constraints on the power spectrum and can effectively recover the power spectrum of wideband sensed signals from sub-Nyquist rate samples, eliminating the necessity to reconstruct the original signal. Subsequently, reference [18] employed the asymptotic theory of circulant matrices to introduce a novel dimensionality reduction technique that simplifies existing structured covariance estimation algorithms, achieving comparable performance at a reduced computational cost. Reference [19] leverages the statistical structure of stochastic processes for signal compression, thereby providing an alternative perspective for inference that is independent of sparsity. Another approach involves utilizing existing sampling structures to directly reconstruct the power spectrum of wideband analog signals from sub-Nyquist samples through optimized algorithms [20,21,22,23,24]. This can significantly reduce the computational burden of power-spectrum reconstruction while providing new insights for the design of fast time–frequency reconstruction algorithms in this study.

In practical applications, such as wireless communication, electronic warfare, and radar, traditional spectrum sensing and time–frequency reconstruction schemes often have limitations. On the one hand, traditional methods usually focus on static power-spectrum analysis and have difficulty effectively capturing the dynamic changes of signals over time. This results in an inability to timely and accurately sense critical information like spectrum holes in rapidly changing signal environments. On the other hand, despite research by many scholars on the time–frequency reconstruction of wideband signals [25,26,27,28], existing algorithms commonly struggle with real-time sensing. Further analysis reveals two main reasons for this issue: first, some algorithms lose temporal information during processing, preventing a complete and accurate representation of the signal’s time–frequency characteristics; second, some algorithms retain temporal information but analyze multiple short time segments, which not only reduces accuracy but also significantly increases computational load, severely limiting real-time sensing capabilities. To address these shortcomings of traditional schemes and inspired by fast compressive covariance sensing [20], this paper proposes a fast time–frequency reconstruction algorithm. It aims to effectively solve the real-time sensing problems of existing algorithms and provide a more efficient and accurate solution for the time–frequency reconstruction of wideband signals. This will better meet the needs of signal dynamic-change monitoring and spectrum resource management in practical applications.

The remainder of this paper is organized as follows. Section 2 summarizes the WBSS model and the fundamental principles of wideband compressive covariance sensing. Section 3 focuses on the proposed FTFR algorithm, detailing its procedure and analyzing its computational complexity. Section 4 presents experimental validations of the algorithm’s reconstruction performance. Finally, Section 5 concludes the paper with a summary.

2. Principles of Compressed Spectrum Sensing

In this section, we explore the principles of wideband compressive covariance sensing. We begin by presenting the model of Wideband Spectrum Sensing to provide a comprehensive review of its theoretical foundations. This is succeeded by a concise overview of existing algorithms for wideband compressive power-spectrum estimation.

2.1. Wideband Compressive Spectrum Sensing Model

In wireless communication environments, wideband signals generally exhibit spectral sparsity, indicating that within the entire frequency band only a limited number of sub-bands contain non-zero power due to signal activity. Building upon this concept, wideband compressive spectrum sensing harnesses the theory of CS, effectively exploiting signal sparsity to significantly reduce the sampling rate requirements for wideband signals while preserving the accuracy of spectrum sensing. Subsequently, we systematically review the comprehensive process of wideband compressive spectrum sensing.

As illustrated in Figure 1, the process of wideband compressive spectrum sensing primarily consists of five steps: sparse representation, sparse measurement, sparse reconstruction, spectrum analysis, and spectrum sensing.

1.Sparse Representation: Let us consider a signal x containing N observations that demonstrates sparsity in a specific domain, with a sparsity level of $[eqn]$ . If the signal does not exhibit sparsity in the designated domain, an appropriate transformation basis can be selected for projection, thereby rendering the signal sparse in the transformed domain. This process can be mathematically represented by Equation (1). In this context, $[eqn]$ represents the transformation basis; if x is already sparse, then $[eqn]$ becomes the identity matrix $[eqn]$ .

$$[eqn]$$

In this equation, s represents the sparse representation of the signal x within the transformed domain.2.Sparse Measurement: Upon completion of the sparse representation, the sparse signal is multiplied by a measurement matrix $[eqn]$ of size $[eqn]$ , where the number of measurements M is associated with the sparsity level, the properties of the measurement matrix, and the chosen reconstruction method, with  $[eqn]$ . The resulting measurement vector is given by Equation (2)

$$[eqn]$$

In this context, y represents the compressed measurements obtained from the sparse signal s, and  $[eqn]$ denotes the measurement matrix used for acquiring these measurements. To realize this process, various sampling structures have been proposed, including the Analog-to-Information Converter (AIC) [17], Modulated Wideband Converter (MWC) [29], and Multi-Coset Sampling (MCS) [30]. These structures can effectively reduce the number of samples required while preserving the key features of the signal.3.Sparse Reconstruction: The low-dimensional signal obtained from sparse measurements cannot be directly utilized for spectrum analysis; consequently, it is essential to employ reconstruction algorithms to approximate the original sparse spectrum signal. Commonly employed reconstruction algorithms include Basis Pursuit (BP) [8], Orthogonal Matching Pursuit (OMP) [14], Compressive Sampling Matching Pursuit (CoSaMP) [31], and Simultaneous Orthogonal Matching Pursuit (SOMP) [32]. These algorithms are fundamentally grounded in optimization theory and rely on the sparsity of the signal to accurately reconstruct the original spectral information.4.Spectrum Analysis: After obtaining the reconstructed signal, techniques such as the Fast Fourier Transform (FFT) can be employed to transform the signal into the frequency domain, facilitating the analysis of its spectral characteristics, which provides the foundation for spectrum sensing.5.Spectrum Sensing: Finally, in the frequency domain, unoccupied spectrum bands, referred to as spectrum holes, can be identified by establishing an energy threshold. This step typically utilizes energy detection [5], which entails scanning the reconstructed spectrum to identify regions where the energy falls below the established threshold, thereby indicating potential available spectrum resources.

2.2. Wideband Compressive Power-Spectrum Estimation

In 2012, Dr. Geert Leus and his research team were among the pioneers in proposing wideband compressive power-spectrum estimation [17]. This method establishes the autocorrelation relationship between sub-Nyquist samples and Nyquist-rate samples through periodic sampling processes in both the time and frequency domains. By not imposing any sparsity constraints, it leverages the cross-correlations of the compressed signal and utilizes simple least squares to reconstruct the signal’s autocorrelation, from which the power spectrum is subsequently calculated. The algorithm for reconstructing the power spectrum in the time domain has significantly informed our research efforts. A brief overview of this algorithm is presented in the following section.

Consider a wide-sense stationary analog signal $[eqn]$ with a bandwidth of $[eqn]$ .The signal $[eqn]$ is sampled below the Nyquist rate using AIC, and the structure of AIC is shown in Figure 2. AIC is a sampling structure with M channels. In the i-th channel, the signal $[eqn]$ is modulated by a complex-valued periodic waveform $[eqn]$ with period $[eqn]$ . In practice, $[eqn]$ can be chosen as a periodic square wave:

$$[eqn]$$

where $[eqn]$ denotes the modulation frequency of the i-th channel. Alternatively, a random modulation approach can be employed, generating $[eqn]$ according to a specific probability distribution, such as complex Gaussian sampling or randomly selecting values from the set $[eqn]$ for binary sampling. After modulation by the complex periodic waveform $[eqn]$ , the signal is integrated by an integrator with period $[eqn]$ . Obviously, the sampling rate of AIC is $[eqn]$ times that of the Nyquist sampling rate. The output at the k-th sampling index of the i-th channel can be expressed as:

$$[eqn]$$

where $[eqn]$ denotes one period of $[eqn]$ , that is

$$[eqn]$$

Assume that $[eqn]$ is a piecewise constant function, with constant values in each interval of length T, i.e., $[eqn]$ ; for $[eqn]$ , we have:

$$[eqn]$$

Figure 3 illustrates the digital interpretation of the AIC sampling device. The process begins with an integrator that produces $[eqn]$ , followed by a set of channels. In the i-th channel, filtering is performed by a filter $[eqn]$ of length N, and then N-times decimated to obtain the final sampled result, i.e.,  $[eqn]$ , where

$$[eqn]$$

• denotes the convolution operator. This completes the compressed sampling operation of the wide-sense stationary analog signal $[eqn]$ . Next, the autocorrelation $[eqn]$ of the original signal is reconstructed using the cross-correlation $[eqn]$ between $[eqn]$ and $[eqn]$ .

According to the definition of cross-correlation and the relationship $[eqn]$ , the cross-correlation between $[eqn]$  and  $[eqn]$ can be expressed as a N-times sampled version of the cross-correlation function between $[eqn]$  and  $[eqn]$ , given by:

$$[eqn]$$

where $[eqn]$ represents statistical expectation, which is used to describe the mean characteristics of random signals. Substituting expression (7) into expression (8), $[eqn]$ can be written as:

$$[eqn]$$

where

$$[eqn]$$

represents the deterministic cross-correlation function between $[eqn]$  and $[eqn]$ . Combining expressions (8)–(10), we obtain:

$$[eqn]$$

where

$$[eqn]$$ $$[eqn]$$ $$[eqn]$$

Concatenating $[eqn]$ distinct cross-correlation functions $[eqn]$ forms the vector $[eqn]$ , expressed as:

$$[eqn]$$

Here, $[eqn]$  and $[eqn]$ are $[eqn]$ matrices, given by:

$$[eqn]$$ $$[eqn]$$

Due to the band-limited nature of $[eqn]$ , the range of $[eqn]$ is restricted to $[eqn]$ , where L can be set according to practical needs. The expressions for $[eqn]$  and $[eqn]$ are shown in Equations (18) and (19), respectively.

$$[eqn]$$ $$[eqn]$$

By examining the ranges of $[eqn]$  and $[eqn]$ , Equation (11) can be expressed as a cyclic convolution, leading to the vector form:

$$[eqn]$$

Here, $[eqn]$ is a $[eqn]$ fast cyclic matrix given by:

$$[eqn]$$

where

$$[eqn]$$ $$[eqn]$$

If $[eqn]$ is a full-rank matrix and $[eqn]$ , Equation (20) can be theoretically solved using the least-squares method. However, due to the potential ill-conditioning of $[eqn]$ , direct least squares may amplify measurement noise and lead to error accumulation. To mitigate this, we introduce the Tikhonov regularization approach:

$$[eqn]$$

where $[eqn]$ is a regularization parameter selected via the L-curve criterion to balance data fidelity and solution smoothness.

Furthermore, the solution to the matrix Equation (20) is complicated by the error-accumulating nature of the special matrix $[eqn]$ in $[eqn]$ . This may introduce extra error and complicates the solution. To address this, in addition to adopting the Tikhonov regularization approach, another strategy involves considering the structural characteristics of $[eqn]$ into a diagonal block matrix $[eqn]$ , which is expressed as:

$$[eqn]$$

We define the vectors $[eqn]$  and $[eqn]$ as the Discrete Fourier Transform (DFT) of $[eqn]$  and $[eqn]$ , respectively, as shown below:

$$[eqn]$$ $$[eqn]$$

Based on these definitions, Equation (20) can be reformulated as:

$$[eqn]$$

By expressing $[eqn]$  and $[eqn]$ as frequency components, Equation (28) can be rewritten as a set of matrix equations, as shown in Equation (29):

$$[eqn]$$

If $[eqn]$ is full-rank for all l, then $[eqn]$ can be solved using the least-squares method. Subsequently, $[eqn]$ can be reconstructed via the Inverse Discrete Fourier Transform (IDFT), and finally, the power spectrum of the wide-sense stationary analog signal $[eqn]$ can be reconstructed by performing a DFT on $[eqn]$ .

3. Fast Time–Frequency Reconstruction Algorithm

In Wideband Spectrum Sensing, capturing the dynamic changes of signals over time is a key challenge. Traditional methods often focus on static power-spectrum analysis, which fails to effectively capture these dynamic changes. This can lead to misjudgments in complex electromagnetic environments. To address this, inspired by the efficient power-spectrum reconstruction algorithms utilizing Multi-Coset Sampling (MCS) referenced in [17,20], this paper introduces an effective time–frequency analysis algorithm. The core of this algorithm resides in the integration of sub-Nyquist sampling technology for signal compression, along with time-window segmentation and fast autocorrelation estimation, which facilitates the rapid construction of the time–frequency distribution in both the time and frequency domains. By reconstructing the time–frequency map, the algorithm accurately captures the dynamic changes in the signal, enhancing the real-time capabilities of spectrum sensing and effectively mitigating misjudgment issues in power-spectrum reconstruction caused by a lack of temporal information in complex electromagnetic environments.

3.1. Multi-Coset Sampling

According to the statistical data presented in [4], the Multi-Coset Sampling (MCS) scheme requires the least computational effort and is also the simplest to implement, establishing it as the predominant compressed sampling structure currently available. MCS utilizes a periodic non-uniform sampling approach. In this framework, a continuous set of L sampling points is regarded as a single sampling block, within which p sampling points are selected, with  $[eqn]$ , maintaining a period of L, where T denotes the Nyquist sampling interval. During implementation, MCS is typically designed as a sampling structure that contains M channels, each responsible for processing a portion of the sampled data, as illustrated in Figure 4. According to compressive sensing theory, the parameter M must satisfy the theoretical lower bound $[eqn]$ , where K is the sparsity level of the signal, to ensure the Restricted Isometry Property (RIP) of the measurement matrix and guarantee unique reconstruction. However, in practical sparse signal recovery, moderate redundancy (e.g., $[eqn]$ ) beyond the theoretical minimum is required to address non-ideal sparsity and noise interference.

In the $[eqn]$ channel, the signal is initially subjected to a fixed delay $[eqn]$ . Subsequently, it is sampled by a low-speed analog-to-digital converter (ADC) at a sampling rate of $[eqn]$ , where N denotes the downsampling factor that represents the compression ratio of MCS relative to Nyquist sampling. A smaller value of N corresponds to higher resource utilization. The coefficient $[eqn]$ is an integer less than N, and this set is given by

$$[eqn]$$

where

$$[eqn]$$

When designing $[eqn]$ , for any integer n that satisfies $[eqn]$ , there exist two indices $[eqn]$ and $[eqn]$ such that $[eqn]$ , where $[eqn]$ denotes the floor function. A detailed mathematical proof can be found in [17].

At the sampling moment l, the output of the $[eqn]$ branch is

$$[eqn]$$

For all channels’ sampling points, Equation (32) can be written in matrix form as

$$[eqn]$$

where $[eqn]$ is a selection matrix with only one non-zero element per row. From the above analysis, it is evident that the sampling points of MCS constitute a subset of the Nyquist sampling points.

3.2. Fast Time–Frequency Reconstruction Algorithm

The overall flowchart of the FTFR algorithm is illustrated in Figure 5. First, the received signal $[eqn]$ undergoes MCS to acquire sub-Nyquist samples. Next, utilizing the delay matrix $[eqn]$ and the compressed sampling results $[eqn]$ , an indicator sequence $[eqn]$ and data sequence $[eqn]$ are constructed. Short-Time Fourier Transform (STFT) is applied to $[eqn]$ and $[eqn]$ , respectively, to derive a matrix containing time–frequency information. Subsequently, the Inverse Fast Fourier Transform (IFFT) is computed row-wise for this matrix, and the elements at corresponding positions within the matrices are divided to compute the signal’s autocorrelation across different time windows. The resulting matrices, which encapsulate time–frequency information, are subsequently flattened after concatenation. Finally, STFT is applied to the flattened autocorrelation matrix to conduct time–frequency analysis. It is noteworthy that all operations depicted in the diagram are matrix operations, where $[eqn]$ denotes element-wise division. A detailed derivation of the algorithm will be provided in the subsequent sections.

To provide a clearer elucidation of the algorithmic process, Figure 5 is decomposed into a three-dimensional diagram as shown in Figure 6. The introduction of the STFT in place of the FFT fundamentally incorporates the use of a time window. Its primary objective is to integrate temporal information into the signal-processing framework. Theoretically, this implies that the signal undergoes initial processing through segmentation using time windows. When receiving the signal $[eqn]$ , the segment of the signal within each time window is represented by

$$[eqn]$$

where $[eqn]$ denotes the time window, and  $[eqn]$ represents the center position of the $[eqn]$ time window. Each window has a width W, and it moves with a step size $[eqn]$ , which is generally set to $[eqn]$ . After applying MCS to the windowed signal, the output is

$$[eqn]$$

To establish the relationship between the sub-Nyquist sampled sequence and the Nyquist sampled sequence, we define the data sequence $[eqn]$ and indicator sequence $[eqn]$ within the $[eqn]$ time window as follows

$$[eqn]$$ $$[eqn]$$

It is straightforward to determine that

$$[eqn]$$

According to the Wiener–Khinchin theorem, the power-spectral density of a wide-sense stationary random process is given by the Fourier transform of its autocorrelation function. Therefore, the power spectrum within the $[eqn]$ window is

$$[eqn]$$

To achieve continuous time–frequency analysis, it is essential to fuse spectra from adjacent time windows, which can be accomplished through weighted averaging. Assuming that two adjacent time windows have the power spectrum $[eqn]$ and $[eqn]$ , these can be fused using a weighting function $[eqn]$

$$[eqn]$$

where $[eqn]$ represents the transition relationship between $[eqn]$ and $[eqn]$ . In this algorithm, the weighting function is characterized by a linear function, specifically $[eqn]$ . By fusing the complete time–frequency representation of the signal according to the aforementioned method, the final time–frequency analysis formula is derived as

$$[eqn]$$

where $[eqn]$ is the weighting function for the $[eqn]$ time window.

At this point, our objective is to estimate the autocorrelation function $[eqn]$ within the $[eqn]$ time window, for which we utilize the algorithm presented in [20].

$$[eqn]$$

where $[eqn]$ , $[eqn]$ and $[eqn]$ . Here, $[eqn]$ represents the set of valid indices for computing the autocorrelation at lag k, $[eqn]$ is the subset of those indices where both $[eqn]$ and $[eqn]$ are non-zero. The normalization factor $[eqn]$ corresponds to the size of the set $[eqn]$ .

First, we calculate $[eqn]$ , according to the definition of autocorrelation

$$[eqn]$$

To solve for $[eqn]$ , we define two new sequences $[eqn]$ and $[eqn]$ ,

$$[eqn]$$ $$[eqn]$$

Based on the convolution property, Equation (43) can be simplified to

$$[eqn]$$

where ⊛ represents circular convolution. According to the circular convolution theorem, we obtain

$$[eqn]$$

where $[eqn]$ and $[eqn]$ . $[eqn]$ and $[eqn]$ denote the FFT and IFFT of $[eqn]$ points, respectively.

Subsequently, we calculate $[eqn]$ ,which is defined as

$$[eqn]$$

Similarly, we define two new sequences $[eqn]$ and $[eqn]$ ; that is,

$$[eqn]$$ $$[eqn]$$

By analogy, based on the convolution property, Equation (48) can also be simplified to

$$[eqn]$$

According to the circular convolution theorem, we obtain

$$[eqn]$$

where $[eqn]$ and $[eqn]$ .

To perform time–frequency analysis, we first compute the autocorrelation function of the signal within each time window, as described in Equation (42). Subsequently, based on the Wiener–Khinchin theorem, we derive the power spectrum for each time window from the computed autocorrelation functions. Finally, by applying Equation (41), we obtain a continuous time–frequency analysis result through the weighted fusion of power-spectral information from each time window. The entire algorithm process is summarized and presented in Algorithm 1.

3.3. Computational Complexity

The FTFR algorithm proposed in this paper swiftly reconstructs the power spectrum within each time window by applying a width-W time window to the signal $[eqn]$ , followed by the fusion of the final time–frequency information. Assuming that the total signal length is K, and that each window moves with a step size $[eqn]$ , a total of $[eqn]$ windows are required, where $[eqn]$ denotes the ceiling function for rounding up. It can be observed that within each time window, the computations mainly involve simple multiplication and division operations as well as FFT/IFFT operations, while the final time–frequency data-fusion stage entails only weighted averaging. Algorithm 1 Fast time–frequency reconstruction algorithm

• Require:   signal $[eqn]$ , time window $[eqn]$ , window width W, window center position $[eqn]$ , sub-sampling factor N, sampling block length L
• Ensure:   Time–frequency representation $[eqn]$ of the signal $[eqn]$
• 1:Initialize an empty array of signal segments $[eqn]$
• 2:for each time window $[eqn]$ do
• 3:    Extract signal segment $[eqn]$ from $[eqn]$ , centered at window center $[eqn]$ with width W
• 4:    Segment $[eqn]$ and append each segment to $[eqn]$
• 5:end for
• 6:Initialize arrays for sub-Nyquist sampling sequence $[eqn]$ and Nyquist sampling sequence $[eqn]$
• 7:for each signal segment $[eqn]$ in $[eqn]$ do
• 8:    Apply Multi-Coset Sampling (MCS) to obtain sub-Nyquist sampling output $[eqn]$
• 9:    Define $[eqn]$ and $[eqn]$
• 10:end for
• 11:Initialize an empty array for power spectra $[eqn]$
• 12:for each time window $[eqn]$ do
• 13:    Calculate the autocorrelation $[eqn]$ for the i-th time window
• 14:    Use autocorrelation estimation and Fast Fourier Transform (FFT) to derive $[eqn]$
• 15:    Append $[eqn]$ to $[eqn]$
• 16:end for
• 17:Initialize the time–frequency representation $[eqn]$
• 18:for adjacent time windows i and $[eqn]$ do
• 19:    Fuse $[eqn]$ and $[eqn]$ using a weighting function
• 20:    Update $[eqn]$ to include the fused power spectrum
• 21:end for
• 22:return the final time–frequency representation $[eqn]$ of the signal

Traditional time–frequency reconstruction algorithms based on compressive sensing typically employ the OMP algorithm, with a single iteration complexity of $[eqn]$ , where D is the signal dimension and O is the number of dictionary atoms. For I time windows, the total complexity becomes $[eqn]$ . Classic covariance-sensing-based methods avoid sparse recovery but require matrix inversion for autocorrelation reconstruction, with a per-window complexity of $[eqn]$ , where L is the correlation matrix dimension), leading to a total complexity of $[eqn]$ for I windows. In contrast, FTFR reduces computational complexity from polynomial to logarithmic operations via FFT/IFFT, also lowering ADC sampling rate requirements. Specifically, to compute the power spectrum of the signal $[eqn]$ within each time window at a spectral resolution of $[eqn]$ , it is necessary to calculate the correlation function $[eqn]$ and the normalization factor $[eqn]$ , followed by performing a $[eqn]$ -points FFT. Consequently, the complexity within each time window is

$$[eqn]$$

The overall complexity for calculating the power spectrum across all time windows is

$$[eqn]$$

During the time–frequency information fusion phase, the weighted fusion of adjacent time windows requires weighted averaging calculations for each frequency point, resulting in a complexity of $[eqn]$ . By combining both aspects, the overall complexity of the FTFR algorithm is

$$[eqn]$$

Since the FFT operations constitute the primary component of the computation, they can be executed in parallel efficiently, rendering this method particularly suitable for real-time applications. Through optimized implementation on high-performance hardware, such as Field-Programmable Gate Arrays (FPGAs), the reconstruction of the time–frequency map can be achieved in a short timeframe. Currently, Intel’s Agilex 7 FPGA can perform 38 trillion floating-point operations per second, enabling it to generate the time–frequency map for a bandwidth of 1 GHz sampled at a spectral resolution of 50 kHz in approximately 644 microseconds, which effectively meets the demands of most real-time sensing applications.

4. Simulation Results

In this section, we assess the performance of the proposed FTFR algorithm through a series of simulation experiments. To test and validate the fundamental performance characteristics of the algorithm, we generated broadband stationary signals for the experiments. Specifically, in our experiments, we processed zero-mean, unit-variance Gaussian white noise using a Finite Impulse Response (FIR) filter to generate five wide-sense stationary signals within the frequency range of [0, 500] MHz. The carrier frequencies of these signals are [55, 135, 220, 315, 420] MHz, respectively, each with a bandwidth of 15 MHz. According to the Nyquist sampling theorem, the sampling frequency must be at least $[eqn]$ . Our objective is to determine the presence of these signals at various moments and frequency positions; therefore, we collected data samples over a sampling duration of 1 s, with the duration of signals on different carriers constrained as follows:

$$[eqn]$$ $$[eqn]$$ $$[eqn]$$ $$[eqn]$$ $$[eqn]$$

The original time–frequency representation of the simulated signals is illustrated in Figure 7a. We employed a Hamming window as the analysis window, selecting a width of 1024 to balance time resolution and frequency resolution. To implement the FTFR algorithm, Multi-Coset Sampling (MCS) was initially employed to obtain sub-Nyquist samples. In this arrangement, we utilized 8 channels, specifically $[eqn]$ , with each channel’s sampling rate set at 60 MHz and different delays $[eqn]$ , which satisfy the requirements outlined in [17]. The frequency resolution was established at 50 kHz, corresponding to a power-spectrum length of 20,000. To achieve this resolution, each sampling channel was required to collect at least $[eqn]$ samples. To more closely simulate real-world scenarios, we added noise to the signals. Figure 7b illustrates the time–frequency representation with a Signal-to-Noise Ratio (SNR) set at 0 dB. Unless otherwise specified, subsequent simulations will also use an $[eqn]$ dB. After applying the FTFR algorithm to reconstruct the time–frequency analysis, we set a threshold to clearly identify spectral holes based on the signal’s power-spectrum characteristics, i.e.,

$$[eqn]$$

where k is a constant less than 1 used to adjust the sensitivity of the threshold, and $[eqn]$ denotes the average value of the power spectrum.

We first compare and analyze the results of reconstructing the time–frequency representation using different window functions, with the fixed window width set to 1024, including the application of a threshold. As illustrated in Figure 8, experiments were conducted using the Hann window, Hamming window, Blackman window, and Gaussian window. From the figures, it is evident that in the time–frequency representations of Figure 8a Hann window, Figure 8b Hamming window, and Figure 8c Blackman window, the energy distribution of signals across different carriers appears relatively blurred, resulting in lower time and frequency resolutions. Conversely, in the time–frequency representation of Figure 8d Gaussian window, the signal’s energy distribution is clearer, yielding higher time and frequency resolutions. This difference may be attributed to the wider main lobes of the Hann, Hamming, and Blackman windows, which result in lower time resolution and a reduced ability to accurately capture rapid changes in the signal. In contrast, the Gaussian window offers improved localization in both the time and frequency domains, thereby simultaneously providing superior time and frequency resolutions. It is important to note that the computational complexity of the Gaussian window is higher than that of other windows, and it requires the adjustment of the shape parameter and standard deviation to achieve optimal performance. Relevant simulation experiments will be discussed in detail in future research.

The width of the window function plays a crucial role in time–frequency analysis, directly influencing the resolution of the signal in both the time and frequency domains. According to the Uncertainty Principle, it is impossible to simultaneously achieve high time resolution and high frequency resolution in time–frequency analysis. Consequently, we further investigate the effect of various time-window widths on the reconstruction of the signal’s time–frequency representation. As illustrated in Figure 9, panels (a–d) depict the reconstructed time–frequency representations after applying a threshold for window widths of 256, 512, 1024, and 2048, respectively. A window width of 256 offers high time resolution, allowing for the capture of short-term dynamic changes in the signal; however, it results in lower frequency resolution, which makes it challenging to accurately distinguish between signals with close frequencies, potentially leading to misjudgments of spectral holes. Conversely, a time-window width of 2048 provides higher frequency resolution; however, this comes at the expense of time resolution, risking the oversight of many briefly appearing spectral holes. In comparison, window widths of 512 and 1024 strike a better balance between frequency and time resolution. Considering the application context of spectrum sensing in TDCS, these window widths enable the precise and efficient identification of spectral holes.

To comprehensively validate the performance of the FTFR algorithm, we compare the time–frequency reconstruction results under different SNR conditions in Figure 10. An analysis of this figure reveals that at $[eqn]$ dB, the FTFR algorithm can effectively restore the time and frequency distribution of the signals. As the SNR decreases, the time–frequency representation of the original signal gradually becomes affected by noise, resulting in less distinct signal characteristics. Nonetheless, the FTFR algorithm remains capable of relatively clear reconstruction of the primary time and frequency distributions. When the SNR further drops to −15 dB, the original signal is nearly overwhelmed by noise, resulting in increased interference in the reconstructed time–frequency representation. Despite this, it is still possible to approximately reconstruct the main time and frequency distributions. Therefore, the FTFR algorithm proposed in this paper can achieve satisfactory reconstruction results even under low SNR conditions. Nonetheless, further improvements are necessary to filter out additional noise for improved performance.

In order to further evaluate the performance of the proposed fast time–frequency reconstruction (FTFR) algorithm under various Signal-to-Noise Ratio (SNR) conditions, we have plotted Receiver Operating Characteristic (ROC) curves. True-Positive Rate (TPR) and False-Positive Rate (FPR) are key metrics for evaluating detection performance, defined in this context as follows:

• •True-Positive Rate (TPR): The proportion of actually occupied frequency bands that are correctly identified as occupied by the FTFR algorithm. It is calculated as $[eqn]$ , where TP represents true positives, i.e., the number of frequency bands that are actually occupied and correctly detected; FN represents false negatives, i.e., the number of frequency bands that are actually occupied but not detected.
• •False-Positive Rate (FPR): The proportion of actually unoccupied frequency bands that are incorrectly identified as occupied by the FTFR algorithm. It is calculated as $[eqn]$ , where FP represents false positives, i.e., the number of frequency bands that are actually unoccupied but incorrectly marked as occupied by the algorithm; TN represents true negatives, i.e., the number of frequency bands that are actually unoccupied and correctly identified as such.

As shown in Figure 11. These curves illustrate the detection performance at SNR levels of 0 dB, −5 dB, −10 dB, and −15 dB. In Figure 11, as the SNR decreases from 0 dB to purple −15 dB, we observe a shift towards the lower left corner, indicating that the balance between TPR and FPR deteriorates with increasing noise levels, showing a decrease in signal detection accuracy. Notably, even under challenging low SNR conditions, such as −15 dB, the FTFR algorithm maintains a certain level of detection capability, which can be attributed to its effective noise-reduction capabilities and accurate capture of signal characteristics. Additionally, by analyzing the Area Under the Curve (AUC) of the ROC curves, one can quantify the overall performance of the algorithm at different SNR levels. The closer the AUC value is to 1, the better the discriminative ability of the algorithm. In this study, even at the lowest SNR level, the FTFR algorithm demonstrates relatively high AUC values, proving its good applicability in complex electromagnetic environments.

Subsequently, we investigated the impact of different down-sampling factors N on the time–frequency reconstruction performance under an $[eqn]$ . In the experiment, we set the ADC sampling rates to 100 MHz, 60 MHz, and 50 MHz, corresponding to down-sampling factors N of 10, 17, and 20, with compression ratios of 0.8, 0.48, and 0.4, respectively. We present the reconstructed time–frequency representations at different sampling rates in Figure 12. From this figure, it is evident that as the sampling rate increases, the performance of time–frequency reconstruction also improves. Higher sampling rates facilitate a more accurate capture of the finer details of the signal, thereby enhancing the quality of the reconstruction.

We went on to analyze the ROC curves corresponding to different sampling rates, as shown in Figure 13. At this noise level, with the increase in sampling rate, the balance between TPR and FPR gradually shifts towards the upper right, indicating that a higher sampling rate contributes to improved signal-detection accuracy. Specifically, at a 100 MHz sampling rate, the ROC curve is closest to the top left corner, demonstrating the highest discrimination and optimal detection performance; whereas at a 50 MHz sampling rate, the position of the ROC curve is lower, reflecting relatively weaker detection capability. By comprehensively analyzing Figure 12 and Figure 13, although a higher sampling rate can achieve better performance, increasing the sampling rate to improve reconstruction quality comes with certain costs. On one hand, data volume increases with the rise in sampling rate, which not only raises storage and transmission costs but may also lead to increased processing delays. On the other hand, higher sampling rates require more advanced hardware capabilities, including greater computational power and faster data-processing speeds. Therefore, in practical applications, it is necessary to weigh the selection of sampling rates based on specific needs to achieve the best balance between performance and cost.

Finally, we compare the effect of different MCS compression channel number M on the reconstruction results, as shown in the Figure 14. The values of M corresponding to (a–d) are 5, 6, 8, and 10, respectively. The results show that as M increases, the details in the time–frequency representation become richer, and the dynamic changes of the signal are captured more clearly. Specifically, when $[eqn]$ , that is, it is in a critical state, the time–frequency representation is relatively blurry, making it difficult to discern the details and dynamic changes of the signal. When $[eqn]$ , the details in the time–frequency representation improve, but there is still some blurriness, and the dynamic changes of the signal are not clear enough. When $[eqn]$ , the time–frequency representation is very detailed, with clear dynamic changes of the signal, allowing for accurate capture of the time–frequency characteristics. When $[eqn]$ , the details in the time–frequency representation further increase, but the improvement over $[eqn]$ is not significant, and the computational complexity and data volume also increase. From the previous complexity analysis, although $[eqn]$ can achieve slightly better results than $[eqn]$ , it also brings higher computational complexity. Overall, $[eqn]$ is a good choice as it provides excellent performance in capturing the details and dynamic changes in the time–frequency representation while maintaining lower computational complexity and data volume.

5. Conclusions

In this paper, we address the challenge of time–frequency reconstruction in WBSS by proposing a FTFR algorithm. The objective of this algorithm is to efficiently and accurately recover the time–frequency characteristics of signals from sub-Nyquist samples, thereby facilitating the dynamic selection of spectra in TDCS. By incorporating the MCS structure, the FTFR algorithm acquires sub-Nyquist samples, significantly reducing the required sampling rate. By leveraging low-complexity operations such as STFT, IFFT, and matrix computations, the algorithm achieves precise calculation of the autocorrelation for signals across different windows. Through the fusion of spectra from adjacent time windows, the algorithm effectively captures the dynamic changes in the signals. The proposed algorithm demonstrates lower computational complexity compared to previous methods for reconstructing time–frequency analyses, rendering it suitable for practical applications. A series of simulation experiments validate the computational efficiency and effectiveness of the algorithm, indicating its capacity to approximate the primary time and frequency distributions even under low SNR conditions.

However, during the time–frequency fusion phase, this study primarily considered fusion along the temporal dimension, excluding the incorporation of phase information. If a discontinuity occurs in the phase during the fusion process, it may introduce noise into the frequency spectrum. This limitation may result in misjudgments in the reconstructed time–frequency representation, potentially affecting the accurate capture of transient signal features. In future research, we will concentrate on improving the time–frequency fusion strategy and exploring methods to incorporate phase information, further enhancing the accuracy and reliability of the time–frequency representation.