FFT介绍及python源码编写

本文深入探讨了FFT(快速傅立叶变换)与DFT(离散傅立叶变换)的关系,解析了FFT的来源、原理及其在Python中的实现。详细介绍了DFT的背景、频域采样定理,并通过递归方法阐述了FFT如何减少计算复杂度,最后给出了Python代码示例。

大纲

  1. FFT的来源DFT
  2. FFT与DFT的关系
  3. FFT的python实现

一、FFT的来源DFT
要了解DFT,就必须要先搞懂DFT,FFT可以看作是为了计算方便而简化之后的DFT,而要了解DFT就需要了解它和DTFT和DFS的关系,关于DTFT和DFS的知识网络上已有许多介绍,不再赘述。
DFT即离散傅立叶变换,它的产生是为了解决DTFT(离散时间傅立叶变换)在频域上连续的问题。众所周知,计算机所能处理的只有离散、有限长的序列,这样就首先排除了我们在信号与系统中所熟知的FS、FT、DFS的应用,但着眼于DFT时我们又发现虽然它在时域上满足离散、有限长的条件,但他在频域上却是连续且周期的,所以我们应当对其在频域上的值进行截断和离散化。这就为我们带来了两个问题:

1.频域采样要满足什么条件,才能使得时域图像不失真?
2.频域采样、截断之后,怎样重构时域信号?

  1. 对于频域采样的问题,很容易就联想到奈奎斯特采样定理,二者的原理其实也是一样的,我们通过将带限信号与时域冲激信号相乘进行采样,再对所得频域图像正变换得到周期化(时域离散对应频域周期)之后的频域图像,判断频域图像不混叠的条件从而得到采样定理。频域采样定理同理,我们通过将一个时限信号的频域序列与频域冲激信号相乘进行频域采样,再由卷积定理反变换得到时域图像,发现时域序列进行了周期延拓,判断时域序列不混叠的条件从而得到频域采样定理。

频域取样脉冲S(e^jw)(每2pi内采样N个点,采样频率ws):
在这里插入图片描述
取样后信号频域值:
在这里插入图片描述
对频域取样脉冲进行傅立叶级数展开可得:
在这里插入图片描述

又通过DFT时移性质和δ ( t )的DFT为1可求得频域脉冲信号离散时间傅立叶反变换:
在这里插入图片描述
再由卷积定理得频域取样后时域值:
在这里插入图片描述
可以发现频域的采样对应的就是时域信号的周期延拓。对于实际运用来说,我们要进行处理的序列都是因果序列(长度为M+1)如下:
在这里插入图片描述

进行周期延拓之后如果N<M+1就会发生混叠的状况如下
在这里插入图片描述
至此得出了频域采样定理:进行频域采样之后能保证理想重构时域序列的条件是频域采样点数N大于时域序列长度M,由此也就推导出了DFT的正向变换式:
在这里插入图片描述
2. 对于还原时域信号的问题,我们要从时域频域对偶特性和DFS入手,当我们对时限信号的频域图像进行采样时,在第一问中就已经知道实际上对时域信号进行了一个周期延拓的操作,即一个域离散、另一个域周期。由此对一个周期序列,我们已知IDFS的公式:
在这里插入图片描述

以及DFS的公式:
在这里插入图片描述
从而很容易发现这样一个周期序列做DFS之后的值与原信号做DFT之后的值之差一个1/N,所以我们可以用DFT之后的值计算IDFS乘以系数1/N,再对计算结果加矩形窗只取0~M-1个值就得到了原序列。得到IDFT公式如下
在这里插入图片描述
插入:此段可以跳过,为记录时域脉冲FT的一点疑惑

已知冲激函数dalta(t)的傅立叶变换为常数1,那么时域脉冲函数sigma[delta(t-nT),n=(-inf,inf)]由线性性质及时移性质求得FT应该为sigma[exp(-jwnT),n=(-inf,inf)]; 但是如果按照时域脉冲函数的傅立叶级数表示sigma[1/Texp(jnWt),n=(-inf,inf)], 再由1的傅立叶变换为2pisigma(w)和频移特性就求得时域脉冲函数FT为2pi/Tsigma[delta(w-kW),k=(-inf,inf)]. 两个结果在形式上并不相似,那如何让二者统一起来呢?参见下述思路
在这里插入图片描述

二、FFT与DFT的关系
首先由将DFT中时间序列取值分为奇、偶两部分得:
在这里插入图片描述
而后又有虚指数函数的周期性得:
在这里插入图片描述
从而可以将DFT后一半频率的取值与前一半频率的取值联系起来,这样就只需要计算前一半频率对应DFT值,减少了一半的计算量:
在这里插入图片描述
而我们可以发现其实对于每一个F偶和F奇而言,它们都可以视作是对一个新序列(原始序列的奇采样或者偶采样,长度为原来的1/2,u的长度也为原来的1/2)的DFT,他们本身也可以进行上述操作,这样我们就可以递归上述操作,直到最后序列长度为2,此时就可直接使用最原始的公式求之,以长度为16的序列为例(F,Wr的下标表示这一序列应该有多少个值):
在这里插入图片描述
这就是FFT的基本思路了,不断的递归到最后需要计算的实际只有两个频率的值,这样就十分便于计算了。
如果需要计算的是二维DFT,实际上就是按维度先后进行的两个一维DFT:
在这里插入图片描述
至此,完成了FFT的所有介绍。
三、FFT的python实现

"""
作于2020.08.04 尝试写出二维FFT算法
"""
import numpy as np
import cv2 as cv


def fft_2d(img):
    """
    完成了图片的填,分维度计算FFT
    :param img:输入图片,应该是单通道的灰度图
    :return: FFT之后返回值,复数矩阵
    """
    shape=img.shape
    fft_res = np.zeros(shape,dtype=complex)  # 返回值
    N=2
    while(N<shape[0]):
        N = N<<1  # N二进制左移一位,相当于不断乘2
    num1 = N-shape[0]  # 因为FFT要不断除以2递归,所以要求序列的长度应该是2的幂,用此方法计算出需要补多少个零
    N=2
    while(N<shape[1]):
        N = N<<1
    num2 = N-shape[1]
    # 先对图片第一维(行)进行填充后做一维FFT,再取结果中需要的部分
    fft_res = fft_1d(np.pad(img,((num1//2,num1-num1//2),(0,0)),'edge'),0)[num1//2:num1//2+shape[0],:]
    # 对第二维运算
    fft_res=fft_1d(np.pad(fft_res,((num2//2,num2-num2//2),(0,0)),'edge'),1)[:,num2//2:num2//2+shape[1]]
    return fft_res


def fft_1d(img,axis):
    """
    构造旋转因子矩阵 FFT准备工作
    :param img: 待处理的序列
   
评论 5
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值