分数阶傅里叶变换离散算法Ozaktas

离散分数阶Fourier变换(DFRFT)算法

FRFT这篇文献发表于：作者：1一、分数阶Fourier变换的定义二、分数阶与其他时频分析工具（Wigner-Ville分布）的关系三、离散分数阶傅立叶变换的计算2一、分数阶Fourier变换的定义3二、分数阶傅里叶变换与Wigner-Ville分布首先，看一下Wigner-Ville分布是傅里叶变换4经过一系列变换后变为由以上可得，等式的右边是的Wigner-Ville分布，左边是的Wigner-Ville分布也就是说的Wigner-Ville分布，是由的Wigner-Ville分布旋转а角得到。5所以分数阶Fourier变换有一个重要的性质，分数阶Fourier变换是角度为α的时频面旋转.这个性质建立起分数阶Fourier变换与时频分布间的直接联系,并且为分数阶Fourier域理解为一种统一的时频变换域奠定了理论基础,同时也为分数阶Fourier变换在信号处理领域中的应用提供了有利条件。6tωuvαα7三、离散分数阶傅立叶变换的计算目前DFRFT的四种离散化算法8在这篇文献中，第二种，采用分解的方法。91.第一种分解方法可以把以上改写为10假定p∈[-1,1],经过量纲归一化的信号x(t）的分数阶傅里叶变换，可以分解为以下三个步骤：（1）用chirp信号调制信号f(x):（2）调制信号与另一个chirp信号卷积:（3）用chirp信号调制卷积后的信号：式1式2式311具体细节：第一步：将函数，与线性调频函数相乘(式1)。注意，g(x)的频率带宽与时间带宽乘积可以是，f(x)的相应带宽乘积的两倍，所以要求g(x)的采样间隔为1／(2Δx)。如果，()样本值的采样间隔是1／Δx，那么就需要对这些样本值进行插值，然后再与线性调频函数的离散采样值相乘，以得到所希望的g(x)的采样。第二步：将g(x)与一线性调频函数作卷积式(式(2))。注意，由于g(x)是带限信号，所以线性调频函数也可以用其带限形式代替而不会有任何影响。12132、第二种分解方法为了简化计算，人们提出更加有效的分解计算方法。假定x(t)的wigner-ville分布限定在以原点为中心，直径为Δx的圆内。若令，则与chirp信号乘积后的信号在频域具有带宽Δx。可以用Shannon插值表示14简要介绍一下Shannon插值

Shannon定理到设信号，如果存在，使

，，则称是B频率截断的的，这时，只要采样间隔按间隔进行采样就不会损失信息，而且，可按如下公式构造原信号

上式Shannon插值公式。

利用采样序列15163、MATLAB程序functionFaf=frft(f,a)%ThefastFractionalFourierTransform%input:f=samplesofthesignal%a=fractionalpower%output:Faf=fastFractionalFouriertransformerror(nargchk(2,2,nargin));f=f(:);N=length(f);shft=rem((0:N-1)+fix(N/2),N)+1;sN=sqrt(N);a=mod(a,4);17%dospecialcasesif(a==0),Faf=f;return;end;if(a==2),Faf=flipud(f);return;end;if(a==1),Faf(shft,1)=fft(f(shft))/sN;return;endif(a==3),Faf(shft,1)=ifft(f(shft))*sN;return;end%reducetointerval0.5<a<1.5if(a>2.0),a=a-2;f=flipud(f);endif(a>1.5),a=a-1;f(shft,1)=fft(f(shft))/sN;endif(a<0.5),a=a+1;f(shft,1)=ifft(f(shft))*sN;end%thegeneralcasefor0.5<a<1.5alpha=a*pi/2;tana2=tan(alpha/2);sina=sin(alpha);f=[zeros(N-1,1);interp(f);zeros(N-1,1)];18%chirppremultiplicationchrp=exp(-i*pi/N*tana2/4*(-2*N+2:2*N-2)'.^2);f=chrp.*f;%chirpconvolutionc=pi/N/sina/4;Faf=fconv(exp(i*c*(-(4*N-4):4*N-4)'.^2),f);Faf=Faf(4*N-3:8*N-7)*sqrt(c/pi);%chirppostmultiplicationFaf=chrp.*Faf;%normalizingconstantFaf=exp(-i*(1-a)*pi/4)*Faf(N:2:end-N+1);19functionxint=interp(x)%sincinterpolationN=length(x);y=zeros(2*N-1,1);y(1:2:2*N-1)=x;xint=fconv(y(1:2*N-1),sinc([-(2*N-3):(2*N-3)]'/2));xint=xint(2*N-2:end-2*N+3);functi