c++ - 在短时傅里叶变换中获取特定频率的值

Question

我正在尝试使用 C++ 重新创建 spectrogram Matlab 使用的函数。该函数使用 Short Time Fourier Transform (STFT)。我找到了一些 C++ 代码 here执行 STFT。该代码似乎适用于所有频率，但我只想要一些。我找到了 this发布具有以下答案的类似问题:

Just take the inner product of your data with a complex exponential at the frequency of interest. If g is your data, then just substitute for f the value of the frequency you want (e.g., 1, 3, 10, ...)

没有数学背景，我不知道该怎么做。从 Wikipedia page 看，内积部分似乎很简单但我完全不知道他的意思(关于 DFT 的公式)

a complex exponential at frequency of interest

有人可以解释一下我是如何做到这一点的吗？我在 STFT 之后的数据结构是一个充满复数的矩阵。我只是不知道如何提取我想要的频率。

相关函数，其中 window 是 Hamming，所需频率的 vector 还不是输入，因为我不知道如何处理它们:

Matrix<complex<double>> ShortTimeFourierTransform::Calculate(const vector<double> &signal,
    const vector<double> &window, int windowSize, int hopSize)
{
    int signalLength = signal.size();
    int nOverlap = hopSize;
    int cols = (signal.size() - nOverlap) / (windowSize - nOverlap);
    Matrix<complex<double>> results(window.size(), cols);

    int chunkPosition = 0;
    int readIndex;
    // Should we stop reading in chunks? 
    bool shouldStop = false;
    int numChunksCompleted = 0;
    int i;
    // Process each chunk of the signal
    while (chunkPosition < signalLength && !shouldStop)
    {
        // Copy the chunk into our buffer
        for (i = 0; i < windowSize; i++)
        {
            readIndex = chunkPosition + i;
            if (readIndex < signalLength)
            {
                // Note the windowing! 
                data[i][0] = signal[readIndex] * window[i];
                data[i][1] = 0.0;
            }
            else
            {
                // we have read beyond the signal, so zero-pad it!
                data[i][0] = 0.0;
                data[i][1] = 0.0;
                shouldStop = true;
            }
        }
        // Perform the FFT on our chunk
        fftw_execute(plan_forward);

        // Copy the first (windowSize/2 + 1) data points into your spectrogram.
        // We do this because the FFT output is mirrored about the nyquist 
        // frequency, so the second half of the data is redundant. This is how
        // Matlab's spectrogram routine works.
        for (i = 0; i < windowSize / 2 + 1; i++)
        {               
            double real = fft_result[i][0];
            double imaginary = fft_result[i][1];
            results(i, numChunksCompleted) = complex<double>(real, imaginary);
        }
        chunkPosition += hopSize;
        numChunksCompleted++;
    } // Excuse the formatting, the while ends here.
    return results;
}

Answer 1

查找 Goertzel algorithm或 filter例如，使用内积与复指数的计算等效的代码来测量信号中特定固定正弦频率的存在或幅度。性能或分辨率将取决于滤波器的长度和您的信号。