[586fe9]: fourier / ffracft.m Maximize Restore History

Download this file

ffracft.m    118 lines (96 with data), 2.9 kB

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
function frf=ffracft(f,a,varargin)
%FFRACFT Approximate fast fractional Fourier transform
% Usage: frf=ffracft(f,a)
% frf=ffracft(f,a,dim)
%
% `ffracft(f,a)` computes an approximation of the fractional Fourier
% transform of the signal *f* to the power *a*. If *f* is
% multi-dimensional, the transformation is applied along the first
% non-singleton dimension.
%
% `ffracft(f,a,dim)` does the same along dimension *dim*.
%
% `ffracft` takes the following flags at the end of the line of input
% arguments:
%
% 'origin' Rotate around the origin of the signal. This is the
% same action as the |dft|, but the signal will split in
% the middle, which may not be the correct action for
% data signals. This is the default.
%
% 'middle' Rotate around the middle of the signal. This will not
% break the signal in the middle, but the |dft| cannot be
% obtained in this way.
%
% Examples:
% ---------
%
% The following example shows a rotation of the |ltfatlogo| test
% signal:::
%
% sgram(ffracft(ltfatlogo,.3,'middle'),'lin','nf');
%
% See also: dfracft, hermbasis, pherm
%
% References: buma04
% AUTHOR: Christoph Wiesmeyr
% TESTING: ??
if nargin<2
error('%s: Too few input parameters.',upper(mfilename));
end;
definput.keyvals.p = 2;
definput.keyvals.dim = [];
definput.flags.center = {'origin','middle'};
[flags,keyvals,dim,p]=ltfatarghelper({'dim','p'},definput,varargin);
[f,L,Ls,W,dim,permutedsize,order]=assert_sigreshape_pre(f,[],dim,upper(mfilename));
% correct input
a=mod(a,4);
if flags.do_middle
f=fftshift(f);
end;
% special cases
switch(a)
case 0
frf=f;
case 1
frf=fft(f)/sqrt(L);
case 2
frf=flipud(f);
case 3
frf=fft(flipud(f));
otherwise
% reduce to interval 0.5 < a < 1.5
if (a>2.0), a = a-2; f = flipud(f); end
if (a>1.5), a = a-1; f = fft(f)/sqrt(L); end
if (a<0.5), a = a+1; f = ifft(f)*sqrt(L); end
% general setting
alpha = a*pi/2;
tana2 = tan(alpha/2);
sina = sin(alpha);
% oversample and zero pad f (sinc interpolation)
m=norm(f);
f=ifft(middlepad(fft(f),2*L))*sqrt(2);
f=middlepad(f,4*L);
% chirp multiplication
chrp = fftshift(exp(-i*pi/L*tana2/4*((-2*L):(2*L-1))'.^2));
f=f.*chrp;
% chirp convolution
c = pi/L/sina/4;
chrp2=fftshift(exp(i*c*((-2*L):(2*L-1))'.^2));
frf=(pconv(middlepad(chrp2,8*L),middlepad(f,8*L)));
frf(2*L+1:6*L)=[];
% chirp multiplication
frf=frf.*chrp;
% normalize and downsample
frf(L+1:3*L)=[];
ind=ceil(L/2);
ft=fft(frf);
ft(ind+1:ind+L)=[];
frf=ifft(ft);
frf = exp(-i*(1-a)*pi/4)*frf;
frf=normalize(frf)*m;
end;
if flags.do_middle
frf=ifftshift(frf);
end;
frf=assert_sigreshape_post(frf,dim,permutedsize,order);