在短时傅里叶变换中求特定频率的值

Question

我试图使用C++重新创建Matlab使用的谱图函数。该函数使用短时傅里叶变换 (STFT)。我找到了一些执行STFT的C++代码这里。这个代码似乎对所有频率都很有效，但我只想要几个。我在这上找到了一个类似问题的帖子，答案如下：

只需在感兴趣的频率上以复杂指数表示数据的内积即可。如果g是你的数据，那就用f代替你想要的频率的值(例如，1，3，10，.)

没有数学背景，我想不出该怎么做。从维基百科页面看，内部积部分看起来很简单，但我完全不知道他的意思(关于DFT的公式)。

在感兴趣的频率上的复指数

有人能解释一下我是怎么做到的吗？STFT之后的数据结构是一个充满复数的矩阵。我只是不知道怎么提取我想要的频率。

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

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算法或过滤器代码，该代码使用内部积的计算等效值来度量信号中特定的平稳正弦频率的存在或幅度。性能或分辨率将取决于滤波器的长度和信号。