fft - Decimation in time

Hello

I have this code of a fast fourier transform decimation in time(fft_DIT). I need to change into a fft-decimation in frequency. the difference between DIT and DIF is: -the order of the samples: -On DIT the input is bit-reversed order and the output is natural order; -On DIF the input is natural order and the output is bit-reversed order;

-the butterfy commutation: -On DIT Multiplication is done before additions; -On DIF Multiplication is done after additions;

I think these are the differences. My problem is understand this code and adapt it into a fft decimation in frequency..

I appreciate any help, thank you

if true
  % 
est=4;
N=2^est;
re=1/N*([0:N-1]);
im=zeros(1,N);
xnovo=re+j*im;
nchg=0; 
N2=N/2; flag=zeros(1,N);
for i=1:N-2,                         
 if(flag(i)==0)
    dv=N2; ic=1; final=0; pos=i;
    for k=1:est,                   
       if (pos/dv >= 1)
          pos=pos-dv;
          final=final+ic;
       end
       dv=dv/2; ic=ic*2;
    end
    if (i~=final)                  
       rev(1+nchg)=i;
       rev(N2+nchg)=final;
       flag(final)=1;
       nchg=nchg+1;
    end
 end   
end
rev(1:nchg);       
rev(N2:N2+nchg-1);  
for i=1:nchg,
 tmp=re(1+rev(i)); 
 re(1+rev(i))=re(1+rev(N2+i-1)); 
 re(1+rev(N2+i-1))=tmp;
 tmp=im(1+rev(i)); 
 im(1+rev(i))=im(1+rev(N2+i-1));
 im(1+rev(N2+i-1))=tmp;
end
grup=N;
butf=1;
phi=2*pi/N;
for i=0:est-1,                            
 grup=grup/2;
 for j=0:grup-1,                        
    for k=0:butf-1,                     
       ptr1=2*j*butf+k;
       ptr2=ptr1+butf;
       arg=phi*grup*k;
       C=cos(arg); S=sin(arg);
       retmp=C*re(1+ptr2)+S*im(1+ptr2);
       imtmp=C*im(1+ptr2)-S*re(1+ptr2);
       re(1+ptr2)=re(1+ptr1)-retmp;
       im(1+ptr2)=im(1+ptr1)-imtmp;
       re(1+ptr1)=re(1+ptr1)+retmp;
       im(1+ptr1)=im(1+ptr1)+imtmp;

      end
   end

 butf=butf*2;
end
for i=1:nchg,
   tmp=re(1+rev(i)); 
   re(1+rev(i))=re(1+rev(N2+i-1)); 
   re(1+rev(N2+i-1))=tmp;
   tmp=im(1+rev(i)); 
   im(1+rev(i))=im(1+rev(N2+i-1));
   im(1+rev(N2+i-1))=tmp;
end
y=(re.^2+im.^2)/N^2;
plot([0:N-1],y);
title('Densidade Espectral');
xlabel('Linha Espectral');
ylabel('Quadrado da Magnitude');

end

Answer 1

[Sunil]

if true

%

est=4;

N=2^est;

re=1/N*([0:N-1]);

im=zeros(1,N);

xnovo=re+j*im;

nchg=0;

N2=N/2; flag=zeros(1,N);

for i=1:N-2,

if(flag(i)==0)

dv=N2; ic=1; final=0; pos=i;

for k=1:est,

if (pos/dv >= 1)

pos=pos-dv;

final=final+ic;

end

dv=dv/2; ic=ic*2;

end

if (i~=final)

rev(1+nchg)=i;

rev(N2+nchg)=final;

flag(final)=1;

nchg=nchg+1;

end

end

end

rev(1:nchg);

rev(N2:N2+nchg-1);

for i=1:nchg,

tmp=re(1+rev(i));

re(1+rev(i))=re(1+rev(N2+i-1));

re(1+rev(N2+i-1))=tmp;

tmp=im(1+rev(i));

im(1+rev(i))=im(1+rev(N2+i-1));

im(1+rev(N2+i-1))=tmp;

end

grup=N;

butf=1;

phi=2*pi/N;

for i=0:est-1,

grup=grup/2;

for j=0:grup-1,

for k=0:butf-1,

ptr1=2*j*butf+k;

ptr2=ptr1+butf;

arg=phi*grup*k;

C=cos(arg); S=sin(arg);

retmp=C*re(1+ptr2)+S*im(1+ptr2);

imtmp=C*im(1+ptr2)-S*re(1+ptr2);

re(1+ptr2)=re(1+ptr1)-retmp;

im(1+ptr2)=im(1+ptr1)-imtmp;

re(1+ptr1)=re(1+ptr1)+retmp;

im(1+ptr1)=im(1+ptr1)+imtmp;

end

end

butf=butf*2;

end

for i=1:nchg,

tmp=re(1+rev(i));

re(1+rev(i))=re(1+rev(N2+i-1));

re(1+rev(N2+i-1))=tmp;

tmp=im(1+rev(i));

im(1+rev(i))=im(1+rev(N2+i-1));

im(1+rev(N2+i-1))=tmp;

end

y=(re.^2+im.^2)/N^2;

plot([0:N-1],y);

title('Densidade Espectral');

xlabel('Linha Espectral');

ylabel('Quadrado da Magnitude');

end