傅立叶变换

傅立叶变换最开始是在大学时学过的,当初的印象是它跟音频图像有关系,各种时领频域的名词让人很容易混淆,而且教材解说得很生涩,学过之后大都忘得一干二净。偶然看见知乎上有一篇文章详细生动地具体描述了傅立叶变换的含义,看过后令人印象深刻。傅立叶变换具体内容都在这个链接里面,在这篇文章里面就不再复述,这里只是说它的应用和代码实现,为了保证运算速度代码用的是C语言实现。
现实生活中很多信息或者说是信号大都不是连续的,而是离散的。因此这篇文章说的傅立叶变换只针对离散傅立叶变换。

傅立叶变换(DFT)公式如下所示:(简书上不能使用 MathJax 编辑公式 ,有点遗憾,只能截图)


DFT公式.png

x(n)为输入函数,x(k) 由实数和虚数组成( x(k) ( 0 < k < N -1) 是x(n) 中任意某一个点) X(k)则为傅立叶变换后的函数
该公式的算法复杂度是O(N^2)
为了方便用代码写出来可以将这个公式拆分:


公式1.png

然后DFT公式可转换为:


公式2.png
公式3.png

由公式三可以写出DFT公式

/**
 离散傅立叶变换
 
 @param x 原始实数
 @param y 原始虚数
 @param a 变换后实数
 @param b 变换后虚数
 @param N 变换的个数
 @param sign 为 -1时 表示进行逆变换
 */
void dft(double x[],double y[],double *a, double *b, int N,int sign){
    float q ,W,d,c,s;
    q = 2 * M_PI/N;
    for (int k = 0; k < N ; k++) {
        a[k] = b[k] = 0;
        W = q * k;
        for (int  n = 0; n < N; n++) {
            d = W * n;
            c = cos(d);
            s =  sin(d) * sign;
            a[k] += x[n]*c + y[n]*s;
            b[k] += y[n]*c - x[n]*s;
        }
    }
    if(sign == -1){
        c = 1.0/N;
        for (int k = 0; k < N ; k++) {
            a[k] = c * a[k];
            b[k] = c * b[k];
        }
    }
}

通常来说,一般应用傅立叶变换处理的数据量是比较大的,复杂度O(N^2)这个方法无法保证程序高效率运行。因而出现了快速傅立叶变换(FFT)。
FFT是在DFT的基础上进行的快速变换:
x(n)可以根据奇偶数拆分:

公式4.png

结合公式4和DFT公式得出

公式5.png

由于权系数W有对称性


公式6.png

因而有以下公式:

公式7.png

将x(n)不断按照上面方式拆分直至无法拆分为止,这样就有如下图形算法的出现

FFT.jpeg

上图算法称为蝶形算法,图中 W(k,N)为旋转因子
每个蝶形算拆分如下:

FFT_SPLIT.png

FFT算法复杂度是O(N*log2N)
FFT算法的前提是N必须保证是N的2次幂,如果不够则补零直到条件成立为止。
公式看起来麻烦其实理解后也很简单。算法介绍完毕就开始写代码。

快速傅立叶变换(FFT)

/**
 离散快速傅立叶变换
 @param x 原始实数 函数结束后保存变换后的实数数值
 @param y 原始虚数 函数结束后保存变换后的虚数数值
 @param N 变换的个数
 @param sign 为 -1时 表示进行逆变换
 */
void FFT(double x[],double y[], int N,int sign){
    
    int i,j,k,L,m = 0,n1,n2;
    double c,c1,e,s,s1,t,tr,ti;
    double M = log2(N);
    if(floor(M) != ceil(M)){
        double *a =  (double *)malloc(sizeof(double) * N);
        double *b =  (double *)malloc(sizeof(double) * N);
        dft(x, y, a, b, N, 1);
        for (int i = 0; i < N; i++) {
            x[i] = a[i];
            y[i] = b[i];
        }
        free(a);
        free(b);
        
    }else{
        m = M;
    
    n1 = N - 1;
    for (j = 0, i = 0;i < n1;i++){
        if(i < j){
            tr = x[j];
            ti = y[j];
            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 = M_PI/n2;
        c = 1.0;
        s = 0.0;
        c1 = cos(e);
        s1 = -sign * sin(e);
        for (j = 0; j < n2; j++) { // [0..1)  [0..2)  [0..4)
            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 * 1.0);
            y[i] /= (N * 1.0);
        }
    }
        }
}

以上就是所有内容,下一篇文章就讲讲下傅立叶变换在开发中的应用。该歇会了。。。
博客原文:click here

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 216,240评论 6 498
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 92,328评论 3 392
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 162,182评论 0 353
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 58,121评论 1 292
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 67,135评论 6 388
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 51,093评论 1 295
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 40,013评论 3 417
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 38,854评论 0 273
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 45,295评论 1 310
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 37,513评论 2 332
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 39,678评论 1 348
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 35,398评论 5 343
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 40,989评论 3 325
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 31,636评论 0 22
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,801评论 1 268
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 47,657评论 2 368
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 44,558评论 2 352

推荐阅读更多精彩内容