| 网站首页 | 业界新闻 | 小组 | 威客 | 人才 | 下载频道 | 博客 | 代码贴 | 在线编程 | 编程论坛
欢迎加入我们,一同切磋技术
用户名:   
 
密 码:  
共有 1732 人关注过本帖
标题:快速傅立叶变换FFT C程序的疑问
取消只看楼主 加入收藏
岁月如刀
Rank: 6Rank: 6
来 自:冰冻星球
等 级:侠之大者
威 望:7
帖 子:165
专家分:477
注 册:2013-7-21
结帖率:50%
收藏
已结贴  问题点数:20 回复次数:14 
快速傅立叶变换FFT C程序的疑问
最近在看《数字信号处理C语言程序集》,其中1.2 快速傅立叶变换的程序fft.c,按照书上的程序和例题,运行后的结果与书上的不一致。请高手帮忙看看是咋回事。
fft.c如下:
#include <math.h>

void fft(x, y, n, sign)
int n, sign;
double x[], y[];
{
    int i, j, k, l, m, n1, n2;
    double c, c1, e, s, s1, t, tr, ti;

    for (j=1,i=1; i<16; i++)
    {
        m = i;
        j = 2 * j;
        if (j == n)
        {
            break;
        }
    }

    n1 = n - 1;
    for (j=0,i=0; i<n1; i++)
    {
        if (i < j)
        {
            tr = x[j];
            ti = y[i];
            x[j] = x[i];
            y[j] = y[i];
            x[i] = tr;
            y[i] = ti;
        }

        k = n / 2;
        while (k < (j+1))
        {
            j = j - k;
            k = k / 2;
        }

        j = j + k;
    }

    n1 = 1;
    for (l=1; l<=m; l++)
    {
        n1 = 2 * n1;
        n2 = n1 / 2;
        e = 3.14159265359 / n2;
        c = 1.0;
        s = 0.0;
        c1 = cos(e);
        s1 = -sign * sin(e);

        for (j=0; j<n2; j++)
        {
            for (i=j; i<n; i+=n1)
            {
                k = i + n2;
                tr = c * x[k] - s * y[k];
                ti = c * y[k] + s * x[k];
                x[k] = x[i] - tr;
                y[k] = y[i] - ti;
                x[i] = x[i] + tr;
                y[i] = y[i] + ti;
            }
        

            t = c;
            c = c * c1 - s * s1;
            s = t * s1 + s * c1;
        }
    }

    if (sign == -1)
    {
        for (i=0; i<n; i++)
        {
            x[i] /= n;
            y[i] /= n;
        }
    }
}

主程序如下:
#include <math.h>
#include <stdio.h>
#include "fft.c"

main()
{
    int i, j, n;
    double a1, a2, c, c1, c2, d1, d2, q1, q2, w, w1, w2;
    double x[32], y[32], a[32], b[32];

    n = 32;
    a1 = 0.9;
    a2 = 0.3;
    x[0] = 1.0;
    y[0] = 0.0;

    for (i=1; i<n; i++)
    {
        x[i] = a1 * x[i-1] - a2 * y[i-1];
        y[i] = a2 * x[i-1] + a1 * y[i-1];
    }

    printf("\nCOMPLEX INPUT SEQUENCE\n");
    for (i=0; i<n/2; i++)
    {
        for (j=0; j<2; j++)
        {
            printf("  %10.7f + J %10.7f",x[2*i+j],y[2*i+j]);
        }
        printf("\n");
    }

    q1 = x[n-1];
    q2 = y[n-1];
    for (i=0; i<n; i++)
    {
        w = 6.28318530718 / n * i;
        w1 = cos(w);
        w2 = -sin(w);
        c1 = 1.0 - a1 * w1 + a2 * w2;
        c2 = a1 * w2 + a2 * w1;
        c = c1 * c1 + c2 * c2;
        d1 = 1.0 - a1 * q1 + a2 * q2;
        d2 = a1 * q2 + a2 * q1;
        a[i] = (c1 * d1 + c2 * d2) / c;
        b[i] = (c2 * d1 - c1 * d2) / c;
    }   

    printf("\nTHEORETICAL DFT\n");
    for (i=0; i<n/2; i++)
    {
        for (j=0; j<2; j++)
        {
            printf("  %10.7f + J%10.7f",a[2*i+j],b[2*i+j]);
        }
        printf("\n");
    }

    fft(x, y, n, 1);

    printf("\nDISCRETE FOUTIER TRANSFORM\n");
    for (i=0; i<n/2; i++)
    {
        for (j=0; j<2; j++)
        {
            printf("  %10.7f + J %10.7f",x[2*i+j],y[2*i+j]);
        }
        printf("\n");
    }

    fft(x, y, n, -1);

    printf("\nINVERSE DISCRETE FOUTIER TRANSFORM\n");
    for (i=0; i<n/2; i++)
    {
        for (j=0; j<2; j++)
        {
            printf("  %10.7f + J %10.7f",x[2*i+j],y[2*i+j]);
        }
        printf("\n");
    }
}
搜索更多相关主题的帖子: include double C语言 
2013-07-25 15:34
岁月如刀
Rank: 6Rank: 6
来 自:冰冻星球
等 级:侠之大者
威 望:7
帖 子:165
专家分:477
注 册:2013-7-21
收藏
得分:0 
From: TonyDeng
你已经发现了问题,麻烦你尊重一下别人,把发现的现象描述一遍,不要让人浪费时间重复你的发现,基本的礼貌。此帖锁定!

Now. Game over.
2013-07-25 15:35
岁月如刀
Rank: 6Rank: 6
来 自:冰冻星球
等 级:侠之大者
威 望:7
帖 子:165
专家分:477
注 册:2013-7-21
收藏
得分:0 
From: rjsp
运行后的结果与书上的不一致
------ 运行后的结果是什么?书上的结果是什么?尤其是后者,你不说别人知道个球呀!

Now. Game over.
2013-07-25 15:35
岁月如刀
Rank: 6Rank: 6
来 自:冰冻星球
等 级:侠之大者
威 望:7
帖 子:165
专家分:477
注 册:2013-7-21
收藏
得分:0 
From: beyondyf
这贴没什么问题,锁贴的理由很牵强。如果这贴都要锁,那这论坛里没几个值得打开的贴子了。

这段代码没有语法错误,所以也没什么现象可以描述。楼主说的也很清楚,程序的输出与书中的结果不一致。如果非要个现象,这就是现象。

输出结果与标准结果不一致这一定是代码中存在逻辑错误,这种错误往往很难发现。这题有标准结果可对比还好,而实际应用中多数没有标准结果可对比,或者是只有部分样本而程序的逻辑错误又只出现在部分问题域内,也就是说大部分情况下它是对的,只在少数点上会出错。说所有的软件都存在这样的错误有点绝对,但这样的错误绝不在少数,这也是我们时常要给操作系统打补丁的原因。

书归正传,仔细阅读了一下楼主的代码,其中FFT函数中存在一处逻辑错误


Now. Game over.
2013-07-25 15:35
岁月如刀
Rank: 6Rank: 6
来 自:冰冻星球
等 级:侠之大者
威 望:7
帖 子:165
专家分:477
注 册:2013-7-21
收藏
得分:0 
From: beyondyf
你介不介意我删掉你3楼的发言?

Now. Game over.
2013-07-25 15:36
岁月如刀
Rank: 6Rank: 6
来 自:冰冻星球
等 级:侠之大者
威 望:7
帖 子:165
专家分:477
注 册:2013-7-21
收藏
得分:0 
From: love云彩

帖子还好吧,就是问题描述不太好,至少让新手可以借鉴一下
一段程序代码的简练,取决于算法够不够严谨!

Now. Game over.
2013-07-25 15:38
岁月如刀
Rank: 6Rank: 6
来 自:冰冻星球
等 级:侠之大者
威 望:7
帖 子:165
专家分:477
注 册:2013-7-21
收藏
得分:0 
From: beyondyf
你想借鉴什么?这份代码写的很一般,当然这不是楼主的问题,他也是从书上抄的代码想试一下。原书作者是信号处理方面的专家,但不是编程的高手。

所以这贴没什么可借鉴的编程技巧。想看懂这份代码光有C语言的知识是不够的,更重要的是懂傅里叶变换。

这个论坛里的众人大部分是大学生吧,关于傅里叶级数的知识在高数下册中有讲解。你们需要的是重新复习一下这部分知识,然后找本数值分析或信号处理的书再看看快速傅里叶变换的算法原理。

Now. Game over.
2013-07-25 15:38
岁月如刀
Rank: 6Rank: 6
来 自:冰冻星球
等 级:侠之大者
威 望:7
帖 子:165
专家分:477
注 册:2013-7-21
收藏
得分:0 
From: rjsp
我模仿问一题:
int main() { return 0; }
运行后的结果与我猜想的不一致。请高手帮忙看看是咋回事。

Now. Game over.
2013-07-25 15:39
岁月如刀
Rank: 6Rank: 6
来 自:冰冻星球
等 级:侠之大者
威 望:7
帖 子:165
专家分:477
注 册:2013-7-21
收藏
得分:0 
From: beyondyf
我想你误会了,你对他的质疑那是你的观点,我并不想干涉。只是用词让我觉得有点有悖你留给我的印象,所以寻问你的意见。

至于6楼,我则纯是表达对他挑事的不满,呵呵说的明确点免得再误会。

楼主的描述其实很完整,只不过需要多做几个步骤。

复制他提到的书的名字到百度里搜索找到那本书,下载打开他所说的章节就能知道“书上的结果”。

复制楼主的代码贴到编辑器里,编译执行就能知道他的结果。对比一下就能知道哪里不一致了。

呵呵,是费事点,如果了解傅里叶变换可以省去这些步骤。懒得做这几个步骤也可以不做,但不能说人家提供的信息不足。

至于你的例子,呵呵与楼主的类比性不强。

Now. Game over.
2013-07-25 15:39
岁月如刀
Rank: 6Rank: 6
来 自:冰冻星球
等 级:侠之大者
威 望:7
帖 子:165
专家分:477
注 册:2013-7-21
收藏
得分:0 
以下是引用tlliqi在2013-7-25 15:37:49的发言:

此贴应移走

是啊,我移到这里来了

Now. Game over.
2013-07-25 15:41
快速回复:快速傅立叶变换FFT C程序的疑问
数据加载中...
 
   



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

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