| 网站首页 | 业界新闻 | 小组 | 威客 | 人才 | 下载频道 | 博客 | 代码贴 | 在线编程 | 编程论坛
欢迎加入我们,一同切磋技术
用户名:   
 
密 码:  
共有 514 人关注过本帖
标题:请教FFTW3的应用
只看楼主 加入收藏
Lidong005
Rank: 1
等 级:新手上路
帖 子:34
专家分:1
注 册:2010-9-1
结帖率:90%
收藏
 问题点数:0 回复次数:0 
请教FFTW3的应用
各位,我写了一段FFTW3做二维DFT的代码,但是结果不正确,请帮我找找原因何在?

#include "stdafx.h"
#include <iostream>
#include <string>
#include <iomanip>
#include <fstream>
#include <conio.h>
#include <math.h>
#include <complex>
#include <fftw3.h>
#include <stdio.h>
#define M 64
#define N 64
#define Pi 3.14159265359
#define complex_d complex<double>
#define ELEM(r,c) (r*M+c)
double data[M][2*N];


using namespace std;



void dataread( )                      //读取文件数据   经过检验读取的数据是正确的
{
    using namespace std;
    int a[2];
    double b[N+1];
    double c[1];
    streampos weizhi;
    streampos WZhi;


    ifstream infile("ISAR2D_1.dat",ios::in);

    if(!infile)
    {
        cerr<<"open error!"<<endl;
        exit(1);
    }
    for(int i=0;i<2;i++)
    {
        infile>>a[i];
        

    }
    cout<<endl;
    streampos weizhi1 =infile.tellg();                          //weizhi1 是 64 64
   

    infile.seekg(weizhi1);

    for(int i=0;i<N+1;i++)
    {
        infile>>b[i];
    }
    streampos weizhi2 =infile.tellg();                          //weizhi2 第二行结束
      

    infile.seekg(weizhi2);
    infile>>c[0];
    streampos weizhi3=infile.tellg();                         //weizhi3 第三行第二个数开始位置

    weizhi=weizhi3-weizhi2;                                   //weizhi 求出第一个数所占的空间
  ;
    WZhi= weizhi3;
    for(int i=0 ,k=1;i<M;i++,k++)
    {   
        infile.seekg(WZhi);
        for(int j=0;j<2*N;j++)
        {
            infile>>data[i][j];
        }
        streampos weizhi4=infile.tellg();  
        infile.seekg(weizhi4);
        infile>>c[0];
        //cout<<c[k]<<endl;
        streampos weizhi5=infile.tellg();  
        WZhi=weizhi5;
    }
       //c[65]=0.0;
   
}



int showresult(fftw_complex* in, fftw_complex* out)
{


    FILE* shibu = NULL;
    fopen_s(&shibu, "shibu.txt", "w");
    if(shibu == NULL)
    {
        throw "创建数据文件失败。";
    }

    for(int i=0;i<M;i++)
    {
        for(int j=0;j<N;j++)
        {
            fprintf_s(shibu, "%f\t",out[ELEM(i, j)][0]);
        }
        fprintf_s(shibu, "\n");
    }

    FILE* xubu = NULL;                                    
    fopen_s(&xubu, "xubu.txt", "w");
    if(xubu == NULL)
    {
        throw "创建数据文件失败。";
    }


    for(int i=0;i<M;i++)                 
    {
        for(int j=0;j<N;j++)
        {
            fprintf_s(xubu, "%f\t", out[ELEM(i, j)][1]);
        }
        fprintf_s(xubu, "\n");
    }

    fclose(shibu);
    fclose(xubu);
    return 1;
}




int main()
{
    fftw_complex *in, *out;
    fftw_plan p;
    int i, j;
    //complex_d G[M][N];


    //读取数据
    dataread( );

    // 分配存储空间
    in  = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * M* N);
    out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * M *N);

    // 设置变换计划
    p = fftw_plan_dft_2d(M, N, in, out, FFTW_FORWARD, FFTW_ESTIMATE);

    //加窗函数
    for (i=0; i<M; i++)
    {
        for (j=0; j<2*N; j=j+2)
        {
          data[i][j]=data[i][j];//*(0.5-0.5*cos(2*Pi*28/64));
        }
    }

    // 设置测试数据  
     
    for (i=0; i<M; i++)
    {
       int k=0;
      for (j=0; j<2*N; j=j+2)                                      
        {   
            
            in[ELEM(i,k)][0] =data[i][j];
            in[ELEM(i,k)][1] =data[i][j+1];
            //cout<<in[ELEM(i, k)][0]<<endl;
            k=k+1;
        }
        
    }
    // 执行变换
    fftw_execute(p);

    // 显示测试结果
    showresult(in, out);

    // 释放内存
    fftw_destroy_plan(p);
    fftw_free(in);
    fftw_free(out);
    return 1;
}
搜索更多相关主题的帖子: 应用 
2010-09-16 15:15
快速回复:请教FFTW3的应用
数据加载中...
 
   



关于我们 | 广告合作 | 编程中国 | 清除Cookies | TOP | 手机版

编程中国 版权所有,并保留所有权利。
Powered by Discuz, Processed in 0.016951 second(s), 9 queries.
Copyright©2004-2024, BCCN.NET, All Rights Reserved