《76800点fft算法实现》由会员分享，可在线阅读，更多相关《76800点fft算法实现（9页珍藏版）》请在金锄头文库上搜索。
1、一、背景介绍：1.1线性相关：线性相关是讨论两信号之间的同步性(synchronism)或相似性(sim il ari ty) 或同相性(in-phase)或两信号的变化规律是否具有线性关系(linear rela tionship)或接近线性关系的程度。设有离散信号和，其线性相关函数为：G (翔)二为兀(并”5 +抑)打二一 S线性相关运算的简洁表示为：式中“ 表示线性相关运算符(correla tion opera tor)。1.2循环相关：循环相关是针对序列的循环移位的一种相关运算。有限长序列的循环移位是 指，也就是先让序列以N为周期进行周期延拓，然后再进行左移位，只朝一个方 向移位的原
2、因是：对周期序列向左移动一个位置，也就相当于向右移动了 N1 个位置。最后取(0, N-1)的N个值就得到了循环移位后的N个序列值。循环相关运算的简洁表示为：q(网二浙)0 7何式中表示循环相关运算符。1.3线性卷积：卷积运算是线性时不变系统分析的重要工具，很多滤波器的设计中都要用到 卷积运算。给出线性卷积运算的定义，设有离散信号x(n)和y (n),其线性卷积 为：+比答伽)= 0(- n)与线性相关运算不同的是：卷积运算时，y(n)要先反折得到y (-n)；m0 表示y (-n)序列右移，m0表示左移，不同的m得到不同的二.-值。其余与相关计算相同。线性卷积运算的简洁表示为：cxy(m)
3、= x(n)y(n)式中“*”表示线性卷积运算符。1.4循环卷积：循环卷积和循环相关类似，也是针对序列的循环移位的一种运算。有限长序 列的循环移位是指，也就是先让序列y(n)以N为周期进行周期 延拓，再进行反折，然后朝右移位，只朝一个方向移位的原因是：对周期序列向 右移动一个位置，也就相当于向左移动了 N1个位置。最后取(0，N-1)的N 个值就得到了循环移位后的N个序列值。设有序列和，其N点循环卷积为：对-1 勺(耐)=Y双初伽_程)凤何?i=0由于循环移位的关系最后得到的循环卷积的长度就是N点,m取0,1,2,，N-1。循环卷积的简洁表示为：爲(血)=x(n) X)式中和表示循环卷积运算符
4、。1.5 FFT算法实现(一)用FFT计算线性卷积用FFT可以实现两个序列的圆周卷积。在一定的条件下，可以使圆周卷积等 于线性卷积。一般情况，设两个序列的长度分别为N1和N2，要使圆周卷积等于 线性卷积的充要条件是FFT的长度N三N1 + N2。对于长度不足N的两个序列，分 别将他们补零延长到N。当两个序列中有一个序列比较长的时候，我们可以采用分段卷积的方法。有 两种方法：(1)重叠相加法。将长序列分成与短序列相仿的片段，分别用FFT对 它们作线性卷积，再将分段卷积各段重叠的部分相加构成总的卷积输出。(2)重 叠保留法。这种方法在长序列分段时，段与段之间保留有互相重叠的部分，在构 成总的卷积输
5、出时只需将各段线性卷积部分直接连接起来，省掉了输出段的直接 相加。(二)用FFT计算相关函数两个长为N的实离散时间序列x(n)与y(n)的互相关函数定义为：耳-1X-1厂住(m)=迟 jc(琢$(刑 + m)=三 jc(托=无(用)* Xm)r (n)的离散傅里叶变换为：R (k) = X*(k)Y(k),0 k = N + N -1 ,以便使用FFT，将x(n),y(n)补零至长为N。122) 用 FFT 计算 X(k),Y(k),0 k N-1。 R (k) = X*(k)Y(k),0 k N-1。xy3) 对R (k)作 IFFT；取后 N-1 项，得r (m) -N +1 m -1 ；
6、取前N 项，得xyxyr (m) 0 m N - 1 。xy仿真x1=1 2 3 4;x2=1 1 1;L1=length(x1)-1;L2=length(x2)-1;y1=conv(x1,x2);n1=0:1:L1+L2;subplot(211);stem(n1,y1) xlabel(n);ylabel(y1(n);title(直接计算线性卷积)N2=length(x1)+length(x2)-1; x21=x1,zeros(1,N2-length(x1);x22=x2,zeros(1,N2-length(x2);X21=fft(x21);X22=fft(x22); y2=ifft(X21.
7、*X22); n2=0:1:N2-1;subplot(212);stem(n2,y2); xlabel(n);ylabel(y2(n);title(由循环卷积计算线性卷积)111卩111111p100.511.522.533.544.55由循环卷积*算线性卷积II1C1111(j11po1仿真n=0:1:15;xan二exp(-(n-8)2/2);xbn=exp(-0.1*n).*sin(2*pi*0.0625*n); k=length(xbn);xan1=xan zeros(1,k-1);xbn1=xbn zeros(1,k-1);xak=fft(xan1);xbk=fft(xbn1);rm
8、=real(ifft(conj(xak).*xbk);rm1=rm(k+1:2*k-1) rm(1:k);m=(-k+1):(k-1);subplot(2,1,1); stem(m,rm1);xlabel(n);title(线性相关) subplot(2,1,2); xak=fft(xan); xbk=fft(xbn); rm=real(ifft(conj(xak).*xbk); stem(n,rm);xlabel(n);title(循环相关)76800 点 fftFFT 是数字信号分析和处理过程中最常用的重要变换之一。为便于在工程各 种需求中的广泛应用，各大FPGA生产商都已推出了 FFT的
9、IP核。此IP核应用 简便，减少了设计风险，缩短了开发周期。但是，由于IP核局限性很强，其速 度和点数等方面都已经不能满足目前通信或雷达信号侦听技术的快速发展。因此 利用FPGA丰富的逻辑资源，合理设计并行处理结构，实现快速大点数FFT运算 具有重要的实用价值。针对这一问题，借鉴图像处理中的二维信号处理方法，将一维的大点数FFT 算法转换成为适用于矩阵应用的二维FFT算法。在FPGA中设计并实现了并行处 理，在增加运算复杂度的前提下提高了运算的效率，解决了大点数FFT快速运算 的问题。2.1二维FFT的实现原理设序列x(n)的长度为N，且N为2的自然数次幕，其DFT为X(k)。贝I：X(k)
10、=T x(n)Wnk, k = 0,1, ，N -1,Nn=0若N = M xL，将输入的N个点首先写成L行M列的矩阵形式为：p0p1 pM 2pM-1PMPM +1 P2 M - 2P2 M -1PP PP(L-2)M(L-2)M+1(L-1) M - 2(L-1)M-1PP PP(L-1)M(L-1)M+1L M - 2LM-1假设：n = Mn + n , k = Lk + k ,1 0 1 0(n = 0,1，,Ll,n = 0,1，M 1),10(k =0,1,.,M 1,k =0,1,.,L1,);10将n, k代人X(k)的表达式，整理后得到:X( k ) =X(Lk +k )=
11、X(k ,K )=DFTx(n)=DFT x(Mn10X(k ,K ) = DFTx(n1 0 1 艺1 x (n , n )W (Mn 1 + nQ )(Lk10Nn=0 艺 1艺1 x(n ,n ) W “1 k010Ln 0 = 0 n1 = 0艺 1X (k ,n ) W nokW1 0 0 Nn0 = 0 艺 1X (k ,n )W nok1 =100Mn0 = 0由上式可以看出，方括号内是L点的FFT, 共M个。而最外层的求和项是 M点的FFT, 共L个。这样就把一个基于一维处理的FFT运算转换为基于二维 处理的FFT运算。1 0 1 0+ n ) = DFTx(n10, n )
12、=0,n0);+k 0) = W Nn0 k0Mn0 k1 =WMn0 k1=X (k , k201).二维FFT的具体实现步骤如下： 将N点数据表示成L行M列的矩阵； 对矩阵所有的列分别作L点FFT变换(列变换)，得到X (k ,n )；100 将X (k ,n )乘以中间旋转因子W nk得到X(k ,n )；1 0 0 N 1 0 0 将X(k,n)仍看作L行M列的矩阵，对矩阵所有的行分别作M点的FFT100变化(行变换)，得到X (k , k )；201 整序输出，将X (k ,k )整序为X (k ,k )，变换完成。2 0 1 2 1 02.2 mat lab 仿真 代码clcclo
13、se allclear allT=10e-6; B=50e6;K=B/T;%设定线性调频信号t=linspace (0,T/2,76800) ;%时间采样点76800yuanshi二exp(li *pi *K*t .”2);% 产生原始信号y二yuanshi,yuanshi,zeros(l,108544);% 补一倍数据，剩下的补零 ganrao=yuanshi(1,1:5000),zeros(1,16*16384-5000);% 干扰信号，取原始信 号的前5000点数据juzhen1=(reshape(y,16384,16);%将 y 变形为 16 行 16384 列的矩阵 a=fft(juzhen1,);%对矩阵的每列做 fftkO=O;nO=O;%设置旋转因子 while k016%旋转因子为wn while n016384wn(k0+1,n0+1)=exp(-1i*2*pi*k0*n0/(16*16384); n0=n0+1;endn0=0;k0=k0+1;endb=a. *wn;%每列做fft后乘以对应的旋转因子c=fft(b,2);%对每行做fft变换 jieguo=reshape (c, 1, 16*16384) ; %将得出的结果变形为1行16*16384列序列 ganraohou二jieguo.*fft(ganrao);%与干扰 fft 相乘juzhen
