matlab bitrevorder,求MATLAB写的FFT源程序！

function  X = radix2(x,inverse)

x = [1 1 1 1 0 0 0 0 ];

inverse = 1;

switch inverse

case 1

xcheck=x;

p=nextpow2(length(x));

x=[x zeros(1,(2^p)-length(x))];

N=length(x);

M=N/2;                                                                     % Index Helping for Controlling range of indices of the Butterfly at Each Stage (Shrinking)

for stage=1:log2(N);                                                       % No of times decimation has to occur

for index=0:(N/(2^(stage-1))):(N-1);                              % Adjusting the index variations of the butterfly in each stage

for n=0:M-1;                                                  % Within a stage, for a given index (single) as reference, develop a local Butterfly Index

a=x(n+index+1)+ x(n+index+M+1);

AA = exp((-j*(2*pi)/N)*(2^(stage-1))*(n));

a = fi(AA);

b=(x(n+index+1)- x(n+index+M+1)).* AA;

x(n+1+index)=a;

x(n+M+1+index)=b;                                    % In place Computation

end;

end;

M=M/2;                                                                     % Used for creating Butterfly Pairs (INDEXing the wings) (Shrinkage)

end;

X=bitrevorder(x);                                                          % Bit reversing X[k] to obtain X[k] 0

% Ycheck=fft(xcheck,N)

% Cross Check the answer using inbuilt FFT

case 2

xcheck=x;

x=x';

p=nextpow2(length(x));

x=[x zeros(1,(2^p)-length(x))];

N=length(x);

M=N/2;                                                                     % Index Helping for Controlling range of indices of the Butterfly at Each Stage (Shrinking)

for stage=1:log2(N);                                                       % No of times decimation has to occur

for index=0:(N/(2^(stage-1))):(N-1);                              % Adjusting the index variations of the butterfly in each stage

for n=0:M-1;                                                  % Within a stage, for a given index (single) as reference, develop a local Butterfly Index

a=x(n+index+1)+ x(n+index+M+1);

b=(x(n+index+1)- x(n+index+M+1)).*exp((-j*(2*pi)/N)*(2^(stage-1))*(n));

x(n+1+index)=a;

x(n+M+1+index)=b;                                    % In place Computation

end;

end;

M=M/2;                                                                     % Used for creating Butterfly Pairs (INDEXing the wings) (Shrinkage)

end;

X=bitrevorder(x) ;                                                         % Bit reversing X[k] to obtain X[k] 0

X = X'/N

end

% Ycheck=ifft(xcheck,N)

% Cross Check the answer using inbuilt IFFT

end