自动信道检测算法

Posted 2023-03-14

【中文标题】自动信道检测算法

【英文标题】：Algorithm for automatic channel detection

【发布时间】：2017-08-19 09:42:02

【问题描述】：

我目前正在从事一个业余项目，使用 SDR 对无线电信号（更准确地说，我对 L 波段卫星频道感兴趣）执行自动调制分类 (AMC)。我希望它能够发现功率谱中的通道以及它们的频率和带宽，这样我就可以将 AMC 应用到此步骤的输出上。

我的第一个（简单）方法相当简单：在每 N 个样本（比如 1024 个）应用一个窗函数后，对最后一个 N 执行 FFT，计算功率谱的估计值并应用一些指数平滑来减少噪声.然后，在平滑的功率谱中，找到最大和最小信号电平，根据两个电平的加权平均值计算某个阈值，并使用该阈值来确定哪些频率区间属于一个通道。

这在我的单元测试中效果很好（一个 QPSK 通道 + 高斯噪声）。然而，在现实生活中，我要么得到一些渠道，要么得到很多误报。当然，我可以通过微调阈值计算中的权重来解决这个问题，但它就不再是自动的了。

我一直在对 Google 进行一些研究，但也许我没有使用正确的搜索关键字，或者对这个主题没有真正的兴趣（这很奇怪，因为频率扫描仪必须以某种方式执行此任务）。

如何找到平均权重的适当值？也许有比计算整个频谱的阈值更好的方法？

【问题讨论】：

您多久采样一次？ （或者，等效地：您的 1024 个样本跨越了多少时间？）我假设您已经熟悉 Nyquist 采样定理？ （我问是因为 L 波段信号高于 1 GHz，所以如果你采集足够频繁的样本以避免违反该定理，那么 1024 个样本的跨度将小于一微秒，所以看起来你可能会得到一些通过在更长的时间内采集更多样本来改进。）

250000 对于这种特殊情况，我的单元测试没有考虑采样率，它们只是以标准化频率运行。无论如何，不​​要太担心，我使用的 SDR（Nuand 的 BladeRF）能够调谐到 VHF 和 6 GHz 之间的任何频率，所以我不必对整个频谱进行采样（那太疯狂了！）。最后，它就像在 -125000Hz 和 125000Hz 之间对基带信号进行采样。

平滑听起来适得其反。没有理由噪声看起来像频域中的峰值。

【参考方案1】：

您的平滑方法似乎有点适得其反：为什么足够长的 DFT 中的噪声会形成类似于“尖锐”形状的东西？您很可能会以这种方式抑制窄带载波。

然后：确实有很多信号检测器，其中许多只是基于频谱中的能量检测（就像您的方法一样）。

但是，估算器必须越好，您需要的信号信息就越多。您是在寻找 20 kHz 宽窄带信道，还是寻找数十兆赫的高速 QPSK 卫星下行链路？你对通道/脉冲形状有什么了解吗？自相关？

这是一个有点广泛的领域，所以我建议你看看一些已经有效的东西：

Sebastian Müller 的gr-inspector 是一款经过验证且非常棒的自动通道检测器，对于某些类型的传输，还可以推断调制参数。

查看demo video（dubstep 警告！）；你画的好像是这样的：

但是：这只是它可以做的事情之一。

更重要的是：基于 DFT 的能带检测实际上只是检测信号的众多方法之一。例如，它不会检测扩频传输，例如 GPS。为此，您需要了解传播技术或至少需要一个大型的基于自相关的检测器。

【讨论】：

我不是在频率上平滑，而是在时间上。对于功率谱的每次计算（比如spec），我会像这样更新我的谱图：for (i = 0; i < LEN; ++i) spectrogram[i] += alpha * (spec[i] - spectrogram[i]); 我正在寻找任何看起来像 PSK 通道的东西，任何带宽都可以适合我的谱图，这第一步不不必太准确。在获得通道位置的提示后，我可以将频谱居中，进行一些过滤/抽取并应用其他技术来猜测实际的通道参数（波特率等）。

嗯，在这种情况下，只有指数衰减而不是线性 IIR 低通的东西在检测分组链接时会不会表现得更好？

恐怕我不关注你了。为什么指数衰减在这里有意义？能详细点吗？

在现实世界中，你得到的东西是打包的——所以它们只是短暂的“昙花一现”。因此，您的滤波器应该让频率仓中的功率短暂“上升”。然而，如果它是一个合适的线性滤波器，输出会以同样快的速度下降——我们不希望这样，相反，我们希望能够快速上升，然后缓慢下降（首先，然后加速如果真的只是一次爆发）

查看这个 gr-fosphor 可视化的上半部分，以图形方式理解我的意思：youtube.com/watch?v=mjD-l3GAghU

【参考方案2】：

在考虑了 Marcus Müller 的建议后，我开发了一种算法，该算法仍然依赖于全局能量阈值，但也依赖于对本底噪声的估计。它可以在许多方面进行改进，但由于其最简单的实现已经提供了可接受的真实世界捕获结果（同样，250 ksps 的 L 波段捕获）我相信它可能是其他人开始的良好基础。

该算法的中心思想是基于对噪声功率的连续估计来计算能量阈值，并随着频谱图的每次更新而更新它。此估计跟踪在所有 FFT 运行期间每个 FFT bin 达到的最大和最小级别，使用它们来估计该 bin 中的 PSD 并丢弃异常值（即通道、杂散...）。该算法将每隔固定数量的样本执行一次。更正式地说：

算法参数：

int N        /* spectrogram length (FFT bins) */
complex x[N] /* last N samples */
float alpha  /* spectrogram smoothing factor between 0 (infinite smoothing,
                spectrogram will never be updated) and 1 (no smoothing at all) */
float beta   /* level decay factor between 0 (levels will never decay) and 1
                (levels will be equal to the spectrogram). */
float gamma  /* smooth factor applied to updates of the current noise estimation
                between 0 (no updates allowed) and 1 /* no smoothing */
float SNR    /* detection threshold above noise, in dB */

alpha、beta 和 gamma 的推荐值分别为 1e-2、1e-3 和 .5。

变量：

complex dft[N]  /* Fourier transform of the last N samples */
float spect[N]  /* smoothed spectrogram */
float spmin[N]  /* lower levels for spectrogram bins */
float spmax[N]  /* upper levels for spectrogram bins */
int runs        /* FFT run counter (initially zero) */
float N0        /* Current noise level estimation */
float new_N0    /* New noise level estimation */
float min_pwr   /* Minimum power density found */
int min_pwr_bin /* FFT bin where min_pwr is */
int valid       /* Number of valid bins for noise estimation */
int min_runs    /* Minimum number of runs required to detect channels */

算法：

min_runs = max(2. / alpha, 1. / beta)
dft = FFT(x);
++runs;
if (runs == 1) then /* First FFT run */
    spect = dft * conj(dft) /* |FFT(x)|^2 */
    spmin = spect /* Copy spect to spmin */
    spmax = spect /* Copy spect to spmax */
    N0    = min(spect); /* First noise estimation */
else
    /* Smooth spectrogram w.r.t the previous run */
    spect += alpha * (dft * conj(dft) - spect)

    /* Update levels. This has to be performed element-wise  */
    new_N0 = 0
    valid = 0
    min_pwr = INFINITY
    min_pwr_bin = -1
    for (int i = 0; i < N; ++i)
        /* Update current lower levels or raise them */
        if (spect[i] < spmin[i]) then
            spmin[i] = spect[i]
        else
            spmin[i] += beta * (spect[i] - spmin[i]);
        end

        /* Update current upper levels or decrease them */
        if (spect[i] > spmax[i]) then
            spmax[i] = spect[i]
        else
            spmax[i] += beta * (spect[i] - spmax[i]);
        end

        if (runs > min_runs) then
            /* Use previous N0 estimation to detect outliers */
            if (spmin[i] < N0 or N0 < spmax[i]) then
                new_N0 += spect[i]
                ++valid
            end

            /* Update current minimum power */
            if (spect[i] < min_pwr) then
                min_pwr = spect[i]
                min_pwr_bin = i
            end
        end
    end

    /*
     * Check whether levels have stabilized and update noise
     * estimation accordingly
     */
    if (runs > min_runs) then
        /*
         * This is a key step: if the number of valid bins is
         * 0 this means that our previous estimation was
         * absolutely wrong. We reset it with a cruder estimation
         * based on where the minimum value of the current
         * spectrogram was found
         */

        if (valid == 0) then
            N0 = .5 * (spmin[min_pwr_bin] + spmax[min_pwr_bin])
        else
            N0 += gamma * (new_N0 / valid - N0)
        end

        /*
         * Detect channels based on this threshold (trivial,
         * not detailed)
         */

        detect_channels(spect, 10^(SNR / 10) * N0)
    end
end

尽管该算法强烈假设本底噪声是平坦的（在大多数情况下这是错误的，因为在现实世界的收音机中，调谐器输出通过响应不平坦的低通滤波器），它甚至可以工作如果这个条件不成立。这些是不同 alpha 值的一些算法结果，N = 4096 和 SNR = 3 dB。噪声估计用黄色标记，通道阈值用绿色标记，高电平用红色标记，频谱图用白色标记，低电平用青色标记。我还提供了每次 FFT 运行后 N0 估计的演变：

alpha = 1e-1 的结果：

alpha = 1e-2 的结果。请注意，随着频谱图变得更加清晰，有效 bin 的数量是如何减少的：

alpha = 1e-3 的结果。在这种情况下，电平非常紧凑，本底噪声明显不平坦，以至于从一个 FFT 运行到另一个 FFT 存在无效的 bin。在这种情况下，我们回退到寻找具有最低功率密度的 bin 的粗略估计：

min_runs 计算至关重要。为了防止噪声电平上升（即跟随通道而不是本底噪声），我们必须至少等待2. / alpha FFT 运行，然后才能信任信号电平。这个值是通过实验发现的：在我之前的实现中，我直观地使用了 1. / alpha，但对于 alpha = 1e-3 却惨遭失败：

我还没有在其他场景（如突发传输）上对此进行测试.不过，考虑到我正在使用的渠道类型，目前这不是优先事项。