Skip to content

这个fft是怎么使用的? (DSP-156) #100

Open
@WhiteDoveBuct

Description

@WhiteDoveBuct

Answers checklist.

  • I have read the ESP-DSP documentation and the issue is not addressed there.
  • I have read the documentation ESP-IDF Programming Guide and the issue is not addressed there.
  • I have updated my ESP-DSP branch (master or release) to the latest version and checked that the issue is present there.
  • I have searched the issue tracker for a similar issue and not found a similar issue.

General issue report

我使用kissfft库的代码:

Matrix compute_spectrum(const Matrix &strided_input) {
  AUTO_TIME(FFT,"fft")
  int num_frames = strided_input.size();
  int nfft = strided_input[0].size();  // FFT 长度
  int fft_out_size = nfft / 2 + 1;     // rFFT 输出大小

  // 结果存储
  Matrix magnitude_spectrum(num_frames, vector<float>(fft_out_size, 0.0f));

  // 初始化 KissFFT 配置
  kiss_fftr_cfg cfg = kiss_fftr_alloc(nfft, 0, nullptr, nullptr);
  if (!cfg) {
    return {magnitude_spectrum};
  }

  // 计算每一行的 FFT
  for (int i = 0; i < num_frames; ++i) {
    vector<kiss_fft_cpx> fft_output(fft_out_size);
    kiss_fftr(cfg, strided_input[i].data(), fft_output.data());

    // 计算幅度谱或功率谱
    for (int j = 0; j < fft_out_size; ++j) {
      float real = fft_output[j].r;
      float imag = fft_output[j].i;
      float magnitude = sqrt(real * real + imag * imag);
      magnitude_spectrum[i][j] = use_power ? magnitude * magnitude : magnitude;
    }
  }

  // 释放 FFT 资源
  kiss_fft_free(cfg);

  return magnitude_spectrum;
}

我使用kissfft库的结果:

Image

输入是 298 个 512 维的矩阵,512维时域转到(512 /2+1)=257维频率

使用esp-dsp测试的代码:

#define N_SAMPLES 512
int N = N_SAMPLES;
// Input test array
__attribute__((aligned(16)))
float x1[N_SAMPLES];
__attribute__((aligned(16)))
float x2[N_SAMPLES];
// Window coefficients
__attribute__((aligned(16)))
float wind[N_SAMPLES];
// working complex array
__attribute__((aligned(16)))
float y_cf[N_SAMPLES * 2];
// Pointers to result arrays
float *y1_cf = &y_cf[0];
float *y2_cf = &y_cf[N_SAMPLES];

// Sum of y1 and y2
__attribute__((aligned(16)))
float sum_y[N_SAMPLES / 2];

extern const        uint8_t input_start []      asm("_binary_input_bin_start");
extern const        uint8_t input_end   []      asm("_binary_input_bin_end");
extern "C" void app_main()
{
    AUTO_TIME(ALL, "all");

    float *input = (float *)input_start;
    int input_size = (int)(input_end - input_start) / sizeof(float);
    ESP_LOGI(TAG, "input_size = %d", input_size);

    int rows = 298;
    int cols = N;

    esp_err_t ret;
    ESP_LOGI(TAG, "Start Example.");
    ret = dsps_fft2r_init_fc32(NULL, CONFIG_DSP_MAX_FFT_SIZE);
    if (ret  != ESP_OK) {
        ESP_LOGE(TAG, "Not possible to initialize FFT. Error = %i", ret);
        return;
    }
    unsigned int start_b;
    unsigned int end_b;
    {
        AUTO_TIME(fft, "FFT");

        for (int i = 0 ; i < N ; i++) {
            x1[i] = input[i];
            if(i<5){
                std::cout << x1[i] << " ";
            }
        }
        std::cout << std::endl;

        for (int i = 0 ; i < N ; i++) {
            x2[i] = input[i + N];
            if(i<5){
                std::cout << x2[i] << " ";
            }
        }
        std::cout << std::endl;

        // Convert two input vectors to one complex vector
        for (int i = 0 ; i < N ; i++) {
            y_cf[i * 2 + 0] = x1[i];
            y_cf[i * 2 + 1] = 0;
        }
        // FFT
        start_b = dsp_get_cpu_cycle_count();
        dsps_fft2r_fc32(y_cf, N);
        end_b = dsp_get_cpu_cycle_count();
        // Bit reverse
        dsps_bit_rev_fc32(y_cf, N);
        // Convert one complex vector to two complex vectors
        dsps_cplx2reC_fc32(y_cf, N);
    }

    std::vector<float> y1_cf_vec( N / 2);
    std::vector<float> y2_cf_vec( N / 2);
    for (int i = 0 ; i < N / 2 ; i++) {
        y1_cf_vec[i] = sqrt((y1_cf[i * 2 + 0] * y1_cf[i * 2 + 0] + y1_cf[i * 2 + 1] * y1_cf[i * 2 + 1]));
        y2_cf_vec[i] = sqrt((y2_cf[i * 2 + 0] * y2_cf[i * 2 + 0] + y2_cf[i * 2 + 1] * y2_cf[i * 2 + 1]));
        
    }
    for (int i = 0 ; i < 5 ; i++) {
    
        std::cout << y1_cf_vec[i] * y1_cf_vec[i] << " ";
    }
    std::cout << std::endl;
    for (int i = 0 ; i < 5 ; i++) {
    
        std::cout << y2_cf_vec[i] *  y2_cf_vec[i]<< " ";
    }
    std::cout << std::endl;

    ESP_LOGI(TAG, "FFT for %i complex points take %i cycles", N, end_b - start_b);

    ESP_LOGI(TAG, "End Example.");
}

使用esp-dsp测试的结果:

Image

红色 1 是输入数据的一部分
红色 2 是输出数据的一部分
我们可以看到只有 黄色标记的第一个直流通量才是和kissfft的结果一致,其它的都不对,这是为什么呢?还有就是我们默认是257维的结果,为啥本项目中的fft这个demo的结果是256维,我要怎么做才能生成257维的结果

Metadata

Metadata

Assignees

No one assigned

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions