解决使用傅里叶变换开源库fftw分析音频频谱结果与matlab或audacity不一致的问题

发布时间:2023年12月22日

找的一些demo输出结果与实际结果相差巨大,修复后效果如下:

采用一个采样率48000,精度16bit,单通道的46Hz,振幅为32767的正弦波测试(理论上应该得输出一个一模一样的正弦波)。输出如下图,可以看到和matlab或audacity差不多。

fftw测试结果,

audacity输出结果:

源码如下:

#include <iostream>
#include <fstream>
#include <fftw3.h>
#include <math.h>
#include <stdio.h>
#include <sys/types.h>
#include <sys/stat.h>
#include <unistd.h>
#include <stdlib.h>
#include <string.h>
#include <errno.h>
#include <stdio.h>
#include <stddef.h>
#include <sys/statvfs.h>

#include <time.h>
#include <linux/input.h>
#include <sys/types.h>
#include <sys/stat.h>
#include <fcntl.h>
uint64_t bag_get_boot_time() {
? struct timespec ts;

? uint64_t ret = 0;
? if(clock_gettime(CLOCK_MONOTONIC, &ts)!=0){
? // ?LOGE(TAG,"CLOCK_GETTIME fatal ERR %s",strerror(errno));
? ? return 0;
? }
? ret ?= (uint64_t)ts.tv_sec * 1000000U + ts.tv_nsec / 1000U;

? return ret;
}
uint64_t st;
uint64_t et1;
uint64_t et2;
// 读取WAV文件,并返回采样数据
bool readWavFile(const std::string& filename, double*& samples, int& sampleRate, int& numSamples)
{
? ? // 打开WAV文件
? ? std::ifstream file(filename, std::ios::binary);
? ? if (!file)
? ? {
? ? ? ? std::cerr << "无法打开WAV文件" << std::endl;
? ? ? ? return false;
? ? }
? ? // WAV文件头部格式(44字节)
? ? char header[44];
? ? file.read(header, 44);
? ? // 检查是否为PCM格式的WAV文件
? ? if (header[20] != 1)
? ? {
? ? ? ? std::cerr << "不支持非PCM格式的WAV文件" << std::endl;
? ? ? ? return false;
? ? }
? ? // 获取采样率和采样点数
? ? sampleRate = *(int*)&header[24];
? ? numSamples = *(int*)&header[40];
? ? numSamples = 4096;
? ? printf("--- %d %d\n",sampleRate,numSamples);
? ? // 分配内存并读取采样数据
? ? samples = new double[numSamples];
? ??
?? ?// 读取16位有符号整型采样数据
?? ?int16_t* buffer = new int16_t[numSamples];
?? ?file.read((char*)buffer, sizeof(int16_t) * numSamples);
? ? st = bag_get_boot_time();
?? ?// 转换为浮点型数据,并归一化到 [-1.0, 1.0] 范围内
?? ?for (int i = 0; i < numSamples; ++i){
?? ??? ?samples[i] = static_cast<double>(buffer[i]) / 1.f;
? ? ? // ?printf("%f ",samples[i]);
? ? }
? ? printf("\n");
?? ?delete[] buffer;
? ? return true;
}
// 使用FFTW进行频域分析
void performFFT(const double* samples, int numSamples, int sampleRate)
{
? ? // 创建FFTW输入数组
? ? fftw_complex* in = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * numSamples);
? ? // 创建FFTW输出数组
? ? fftw_complex* out = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * numSamples);
? ? // 创建FFTW计划
? ? fftw_plan plan = fftw_plan_dft_1d(numSamples, in, out, FFTW_FORWARD, FFTW_ESTIMATE);
? ? // 将采样数据复制到输入数组中
?? ?for (int i = 0; i < numSamples; ++i)
?? ?{
? ? ? ? in[i][0] = samples[i]; ?// 实部
?? ??? ?in[i][1] = 0.0; ? ? ? ? // 虚部(初始化为0)
?? ?}
? ? // 执行频域变换
?? ?fftw_execute(plan);
?? ?// 输出频域结果(幅度谱)
?? ?for (int i = 0; i < numSamples/2 ; ++i)
?? ?{
? ? ? ? double frequency = static_cast<double>(i) * sampleRate / numSamples;
? ? ? ? double amplitude = sqrt(out[i][0] * out[i][0] + out[i][1] * out[i][1])*2/numSamples;
? ? ? ? std::cout << "频率:" << frequency << "Hz,振幅:" << amplitude ;
? ? ? ? printf(" db=%f %f %f\n",10*log10(amplitude/32768),out[i][0],out[i][1]);
?? ?}
? ? // 清理内存和释放资源
?? ?fftw_destroy_plan(plan);
?? ?fftw_free(in);
?? ?fftw_free(out);
}
int main()
{
? ? std::string filename = "/home/cf/myprj/shock_2.wav";
? ? double* samples;
? ? int sampleRate, numSamples;
? ?
? ? if (!readWavFile(filename, samples, sampleRate, numSamples))
? ? ? ? return 1;
? ? et1 = bag_get_boot_time();
? ? performFFT(samples, numSamples, sampleRate);
? ? et2 = bag_get_boot_time();
? ? printf("c1=%d c2=%d\n",et2-st,et1-st);
?? ?delete[] samples;
? ??
?? ?return 0;
}

文章来源:https://blog.csdn.net/cf125313/article/details/135146463
本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。