teresa_zp 2015-04-14 05:24 采纳率: 100%
浏览 3320
已采纳

fftw使用的问题 频域抽取失败 高手请进~

int FFTw_IFFTw_Fun2(IN dDataArray* pRRIData, IN float fmin, IN float fmax, OUT dDataArray& VLFData)
{
try
{
int nSamples = pRRIData->m_nSamples;
int N = pRRIData->m_nLen;
int Nout = floor(N/2)+1;//实数据的DFT具有 Hermitian对称性,所以只需要计算n/2+1(向下取整)个输出就可以了。比如对于r2c,输入in有n个数据,输出out有floor(n /2)+1个数据。

    fftw_complex* out1;
    fftw_complex* out2;
    double* out3;
    out1 =  (fftw_complex *)fftw_malloc(sizeof(double)*Nout * 2);
    out2 =  (fftw_complex *)fftw_malloc(sizeof(double)*Nout * 2);
    out3 =  (double *)fftw_malloc(sizeof(double)*N);

    /////////////////////////fft变换/////////////////////////////////////////////////
    fftw_plan p = fftw_plan_dft_r2c_1d(N, pRRIData->m_pData, out1,  FFTW_ESTIMATE);     
    fftw_execute(p);// 执行变换 ... 
    fftw_destroy_plan(p);

    /////////////////////////频域抽取/////////////////////////////////////////////////
    float f1 = fmin;
    float f2 = fmax;

    float K = 0.0f;
    for (int m=0; m<Nout-1; m++)
    {
        K=m/(Nout/nSamples);
        if ((K < f1 || K > f2) || (K > (nSamples - f2) && K < (nSamples - f1)))//1/dt为一个频率周期
        {
            out2[m][0] = 0;
            out2[m][1] = 0;
        }
        else 
        {
            out2[m][0] = out1[m][0];
            out2[m][1] = out1[m][1];
        }
    }

    //////////////////////////////////IFFT变换////////////////////////////////////////
    fftw_plan  c2cP;
    c2cP = fftw_plan_dft_c2r_1d( N, out2, out3,  FFTW_ESTIMATE);//需注意FFTW的逆变换没有除以N,即数据正变换再反变换后是原始数据的N倍。
    fftw_execute(c2cP);
    fftw_destroy_plan(c2cP);

    VLFData.m_nLen = N;
    VLFData.m_nSamples = 2;
    VLFData.m_pData = new double[VLFData.m_nLen];
    memset(VLFData.m_pData,0,sizeof(double)*VLFData.m_nLen);

    for (int i=0; i< N; i++)
    {
        VLFData.m_pData[i] = out3[i]/N;//取实数
    }

    fftw_free(out1);
    fftw_free(out2);
    fftw_free(out3);
    TRACE(_T("int  FFTw_IFFTw_Fun2 函数返回1~"));
    return 1;
}
catch (...)
{
    TRACE(_T("int  FFTw_IFFTw_Fun2 函数catch到错误~"));
    return -1;
}

}

        float fmin = 0.0f;
        float fmax = 0.0033f;
        dataOutput.VLFData.Init();
        FFTw_IFFTw_Fun2(&dataOutput.RRIData, fmin, fmax, dataOutput.VLFData);

出来的波形跟原始数据波形极其相似,我用matlab频域抽取的话,是正确的,用C++实现就出了问题。

  • 写回答

2条回答 默认 最新

  • teresa_zp 2015-04-16 06:13
    关注

    找到了!这一行出了问题,K=m/(Nout/nSamples); C++默认=右侧为整形数据 整形赋值给float类型,float没有起到相应的作用,仍然为整形数据。所以抽取频域数据的时候出现了问题~改为K=1.0*m/(Nout/nSamples); 即可

    本回答被题主选为最佳回答 , 对您是否有帮助呢?
    评论
查看更多回答(1条)

报告相同问题?

悬赏问题

  • ¥15 lammps拉伸应力应变曲线分析
  • ¥15 C++ 头文件/宏冲突问题解决
  • ¥15 用comsol模拟大气湍流通过底部加热(温度不同)的腔体
  • ¥50 安卓adb backup备份子用户应用数据失败
  • ¥20 有人能用聚类分析帮我分析一下文本内容嘛
  • ¥15 请问Lammps做复合材料拉伸模拟,应力应变曲线问题
  • ¥30 python代码,帮调试,帮帮忙吧
  • ¥15 #MATLAB仿真#车辆换道路径规划
  • ¥15 java 操作 elasticsearch 8.1 实现 索引的重建
  • ¥15 数据可视化Python