快速傅立叶变换段错误 [C++]

Question

我的快速傅立叶变换出现段错误，我完全不知道为什么。 下面的代码生成一个正弦波，通过 FFT 运行它，试图找到正弦波的频率，但它会抛出一个错误。

#include <iostream>
#include <iomanip>
#include <string>
#include <cmath>
#include <vector>   
#include <complex>

# define M_PI 3.14159265358979323846

using namespace std;
using namespace std::complex_literals;

/*
int main()
{
    ofstream outfile;
    outfile.open("data.dat",ios::trunc | ios::out);
    for(int i=0;i<201;i++)
    {
        float rads = M_PI/180;
        outfile << (float)(3276*sin(1.8*i*rads)+32767) << endl;
    }
    outfile.close();
    return 0;
*/
/*
y[t] = AMPLITUDE * sin (2 * M_PI * 0.15 * t + 0) + ZERO_OFFSET;
f = 15 cycles / NUM_POINTS = 0.15 Hz
To have one full-cycle, loop from y[0:t) where t 
is the time or number of points it takes to have a full cycle (i.e. wavelength)*/

vector<float> FFT (vector<float> p) {

    int n = p.size();
    if (n == 1) { // n is a power of 2
        return p;
    }
    //cout << p.size() << ",";
    //for (auto &z:p) {
        //cout << z << ",";
    //}
    // -1 represents i
    std::complex<double> complex(M_PI*2, 1);
    float w = exp(complex.real()/double(n)); // eulers number raised to argument
    
    vector<float> PE; // evens & odds
    vector<float> PO;

    for (int i = 0; i < p.size(); i++) {
        if (int(p[i]) % 2 == 0) {
            PE.push_back(p[i]);
        } else {
            PO.push_back(p[i]);
        }
    }
    if (PE.size() == 0) {
        vector<float> PE = {0};
    }
    if (PO.size() == 0) {
        vector<float> PO = {0};
    }

    vector<float> YE,YO = FFT(PE), FFT(PO);
    
    vector<float> y(n,0);

    for (int j = 0; j < n/2; j++) {
        y[j] = YE[j] + pow(w,j) * YO[j];
        y[j + n/2] = YE[j] - pow(w, j) * YO[j];
    }

    return y;
}

vector<float> get_frequencies(vector<float> vec, int SR){
    vector<float> frequencies = FFT(vec);

    return frequencies;
}

int main() {

    /*vector<int> vec(500, 0);

    int i = 0;
    for (i; i < 100; i+=10) {
        vec[i] = 10;
    }

    for (i; i < 200; i+=20) {
        vec[i] = 10;
    }

    for (i; i < 300; i+=5) {
        vec[i] = 10;
    }

    for (i; i < 400; i+=2) {
        vec[i] = 10;
    }

    for (i; i < 500; i+=7) {
        vec[i] = 10;
    }*/
    vector<float> vec;
    
    int base_freq = 100;

    #define NB_OF_SAMPLES 10
    #define OFFSET        50
    #define AMPLITUDE     5000


    double angle = 0.0;

    for(int i=0; i < NB_OF_SAMPLES; i++)
    {
        angle += (2 * M_PI * base_freq / NB_OF_SAMPLES);
        vec.push_back(AMPLITUDE * sin(angle) + OFFSET);
    }

    //cout << "RAW:\n";
    //for (auto &z:vec) {
        //cout << z << ",";
        //if (i != vec.size()) {
            //cout << ",";
        //}
    //}
    vector<float> frequencies = get_frequencies(vec, NB_OF_SAMPLES);
    //vector<float> frequencies = {1,0,1,0,1,0,1,0,1,0};
    cout << "\nFREQUENCIES:\n";
    for (auto &s:frequencies) {
        cout << s << ",";
        
    }
    return 0;
}

应该注意的是，我没有数学背景，所以一步一步的解释是最有帮助的。

非常感谢您的帮助， 波士顿

Answer 1

我意识到代码可能有问题的一些事情是： 首先，在您的 FFT 函数中，在声明您编写的复数时，这有点不正确，因为在处理浮点数时，如果您使用常数对其执行任何操作（在本例中为 2 和 1），您需要在常量中使用浮点数，如 2.0 和 1.0。 另外，在声明向量 YE 和 YO 时，您已编写

complex(M_PI*2, 1);

，这是无效的，您应该编写

vector<float> YE,YO = FFT(PE), FFT(PO);

这些只是发现的语法错误。 您可以进行更改并检查错误是否已解决。