fast fourier transform FFT radix-4 for N=1024 algorithm c++ source code implementation

//source:
//https://www.google.ch/patents/US6957241
//http://www.google.ch/patents/US20020083107
//https://www.beechwood.eu/fft-implementation-r2-dit-r4-dif-r8-dif/
//http://www.cmlab.csie.ntu.edu.tw/cml/dsp/training/coding/transform/fft.html
//http://dsp.stackexchange.com/questions/3481/radix-4-fft-implementation

//author marcin matysek (r)ewertyn.PL
//fast fourier transform FFT radix-4 for N=1024 algorithm c++ source code implementation

void fun_fourier_transform_FFT_radix_4_N_1024(int N,std::complex<double> tab[])
{
    const double pi=3.141592653589793238462;
    std::complex<double>  tab2[4096]={};    // tab2[]==N
    std::complex<double>  w1[1]={{1,0}};
    std::complex<double>  w2[1]={{1,0}};
    std::complex<double>  w3[1]={{1,0}};
    std::complex<double>  w4[1]={{1,0}};
    std::complex<double>  w5[1]={{1,0}};
    std::complex<double>  w6[1]={{1,0}};
    std::complex<double> tmp1,tmp2,tmp3,tmp4;
    double tmp5;


    w5[0].real()=0;
    w5[0].imag()=-1;
    w6[0].real()=0;
    w6[0].imag()=1;

    tmp5=2*pi/(N/1);

//stage 1
        for(int i=0;i<(N/4);i++)//256
        {
          w1[0].real()=cos(0*tmp5);
          w1[0].imag()=-sin(0*tmp5);
          w2[0].real()=cos(0*tmp5);
          w2[0].imag()=-sin(0*tmp5);
          w3[0].real()=cos(0*tmp5);
          w3[0].imag()=-sin(0*tmp5);
          w4[0].real()=cos(0*tmp5);
          w4[0].imag()=-sin(0*tmp5);

          tmp1=w1[0]*tab[i+0];
          tmp2=w2[0]*tab[i+N/4];
          tmp3=w3[0]*tab[i+N/2];
          tmp4=w4[0]*tab[i+3*N/4];
         //radix-4
          tab2[i]       =tmp1+tmp2+tmp3+tmp4;
          tab2[i+N/4]   =tmp1-tmp3+w5[0]*(tmp2-tmp4);
          tab2[i+N/2]   =tmp1-tmp2+tmp3-tmp4;
          tab2[i+3*N/4] =tmp1-tmp3+w6[0]*(tmp2-tmp4);
        }
/////////////////////////////////////////////////////////////////////////////////
    for(int m=0;m<1;m++)//stage 2
    {
        for(int i=0;i<4;i++)
        {
            for(int j=0;j<4*4*4;j++)
            {
                w1[0].real()= cos((1*0*i*(j+4*4*4*0)+1*1*i*j)*tmp5);
                w1[0].imag()=-sin((1*0*i*(j+4*4*4*0)+1*1*i*j)*tmp5);
                w2[0].real()= cos((1*1*i*(j+4*4*4*1)+1*0*i*j)*tmp5);
                w2[0].imag()=-sin((1*1*i*(j+4*4*4*1)+1*0*i*j)*tmp5);
                w3[0].real()= cos((1*2*i*(j+4*4*4*1)-1*1*i*j)*tmp5);
                w3[0].imag()=-sin((1*2*i*(j+4*4*4*1)-1*1*i*j)*tmp5);
                w4[0].real()= cos((1*3*i*(j+4*4*4*1)-1*2*i*j)*tmp5);
                w4[0].imag()=-sin((1*3*i*(j+4*4*4*1)-1*2*i*j)*tmp5);

                tmp1=w1[0]*tab2[0*N/(4*4)+N/(4)*i+N/(1)*m+j];
                tmp2=w2[0]*tab2[1*N/(4*4)+N/(4)*i+N/(1)*m+j];
                tmp3=w3[0]*tab2[2*N/(4*4)+N/(4)*i+N/(1)*m+j];
                tmp4=w4[0]*tab2[3*N/(4*4)+N/(4)*i+N/(1)*m+j];
                //radix-4
                tab[0*N/(4*4)+N/(4)*i+N/(1)*m+j]=tmp1+tmp2+       tmp3+tmp4;
                tab[1*N/(4*4)+N/(4)*i+N/(1)*m+j]=tmp1-tmp3+w5[0]*(tmp2-tmp4);
                tab[2*N/(4*4)+N/(4)*i+N/(1)*m+j]=tmp1-tmp2+       tmp3-tmp4;
                tab[3*N/(4*4)+N/(4)*i+N/(1)*m+j]=tmp1-tmp3+w6[0]*(tmp2-tmp4);
            }
        }
    }
////////////////////////////////////////////////////////////////////////////////
    for(int m=0;m<4;m++)//stage 3
    {
        for(int i=0;i<4;i++)
        {
            for(int j=0;j<4*4;j++)
            {
                w1[0].real()= cos((4*0*i*(j+4*4*0)+4*1*i*j)*tmp5);
                w1[0].imag()=-sin((4*0*i*(j+4*4*0)+4*1*i*j)*tmp5);
                w2[0].real()= cos((4*1*i*(j+4*4*1)+4*0*i*j)*tmp5);
                w2[0].imag()=-sin((4*1*i*(j+4*4*1)+4*0*i*j)*tmp5);
                w3[0].real()= cos((4*2*i*(j+4*4*1)-4*1*i*j)*tmp5);
                w3[0].imag()=-sin((4*2*i*(j+4*4*1)-4*1*i*j)*tmp5);
                w4[0].real()= cos((4*3*i*(j+4*4*1)-4*2*i*j)*tmp5);
                w4[0].imag()=-sin((4*3*i*(j+4*4*1)-4*2*i*j)*tmp5);

                tmp1=w1[0]*tab[0*N/(4*4*4)+N/(4*4)*i+N/(4)*m+j];
                tmp2=w2[0]*tab[1*N/(4*4*4)+N/(4*4)*i+N/(4)*m+j];
                tmp3=w3[0]*tab[2*N/(4*4*4)+N/(4*4)*i+N/(4)*m+j];
                tmp4=w4[0]*tab[3*N/(4*4*4)+N/(4*4)*i+N/(4)*m+j];
                //radix-4
                tab2[0*N/(4*4*4)+N/(4*4)*i+N/(4)*m+j]=tmp1+tmp2+       tmp3+tmp4;
                tab2[1*N/(4*4*4)+N/(4*4)*i+N/(4)*m+j]=tmp1-tmp3+w5[0]*(tmp2-tmp4);
                tab2[2*N/(4*4*4)+N/(4*4)*i+N/(4)*m+j]=tmp1-tmp2+       tmp3-tmp4;
                tab2[3*N/(4*4*4)+N/(4*4)*i+N/(4)*m+j]=tmp1-tmp3+w6[0]*(tmp2-tmp4);
            }
        }
    }
////////////////////////////////////////////////////////////////////////////////
    for(int m=0;m<4*4;m++)//stage 4
    {
        for(int i=0;i<4;i++)
        {
            for(int j=0;j<4;j++)
            {
                w1[0].real()= cos((4*4*0*i*(j+4*0)+4*4*1*i*j)*tmp5);
                w1[0].imag()=-sin((4*4*0*i*(j+4*0)+4*4*1*i*j)*tmp5);
                w2[0].real()= cos((4*4*1*i*(j+4*1)+4*4*0*i*j)*tmp5);
                w2[0].imag()=-sin((4*4*1*i*(j+4*1)+4*4*0*i*j)*tmp5);
                w3[0].real()= cos((4*4*2*i*(j+4*1)-4*4*1*i*j)*tmp5);
                w3[0].imag()=-sin((4*4*2*i*(j+4*1)-4*4*1*i*j)*tmp5);
                w4[0].real()= cos((4*4*3*i*(j+4*1)-4*4*2*i*j)*tmp5);
                w4[0].imag()=-sin((4*4*3*i*(j+4*1)-4*4*2*i*j)*tmp5);

                tmp1=w1[0]*tab2[0*N/(4*4*4*4)+N/(4*4*4)*i+N/(4*4)*m+j];
                tmp2=w2[0]*tab2[1*N/(4*4*4*4)+N/(4*4*4)*i+N/(4*4)*m+j];
                tmp3=w3[0]*tab2[2*N/(4*4*4*4)+N/(4*4*4)*i+N/(4*4)*m+j];
                tmp4=w4[0]*tab2[3*N/(4*4*4*4)+N/(4*4*4)*i+N/(4*4)*m+j];
                //radix-4
                tab[0*N/(4*4*4*4)+N/(4*4*4)*i+N/(4*4)*m+j]=tmp1+tmp2+       tmp3+tmp4;
                tab[1*N/(4*4*4*4)+N/(4*4*4)*i+N/(4*4)*m+j]=tmp1-tmp3+w5[0]*(tmp2-tmp4);
                tab[2*N/(4*4*4*4)+N/(4*4*4)*i+N/(4*4)*m+j]=tmp1-tmp2+       tmp3-tmp4;
                tab[3*N/(4*4*4*4)+N/(4*4*4)*i+N/(4*4)*m+j]=tmp1-tmp3+w6[0]*(tmp2-tmp4);
            }
        }
    }
///////////////////////////////////////////////////////////////////////////////
   for(int m=0;m<4*4*4;m++)//stage 5
    {
        for(int i=0;i<4;i++)
        {
            for(int j=0;j<1;j++)
            {
                w1[0].real()= cos((4*4*4*0*i*(j+1*0)+4*4*4*1*i*j)*tmp5);
                w1[0].imag()=-sin((4*4*4*0*i*(j+1*0)+4*4*4*1*i*j)*tmp5);
                w2[0].real()= cos((4*4*4*1*i*(j+1*1)+4*4*4*0*i*j)*tmp5);
                w2[0].imag()=-sin((4*4*4*1*i*(j+1*1)+4*4*4*0*i*j)*tmp5);
                w3[0].real()= cos((4*4*4*2*i*(j+1*1)-4*4*4*1*i*j)*tmp5);
                w3[0].imag()=-sin((4*4*4*2*i*(j+1*1)-4*4*4*1*i*j)*tmp5);
                w4[0].real()= cos((4*4*4*3*i*(j+1*1)-4*4*4*2*i*j)*tmp5);
                w4[0].imag()=-sin((4*4*4*3*i*(j+1*1)-4*4*4*2*i*j)*tmp5);

                tmp1=w1[0]*tab[0*N/(4*4*4*4*4)+N/(4*4*4*4)*i+N/(4*4*4)*m+j];
                tmp2=w2[0]*tab[1*N/(4*4*4*4*4)+N/(4*4*4*4)*i+N/(4*4*4)*m+j];
                tmp3=w3[0]*tab[2*N/(4*4*4*4*4)+N/(4*4*4*4)*i+N/(4*4*4)*m+j];
                tmp4=w4[0]*tab[3*N/(4*4*4*4*4)+N/(4*4*4*4)*i+N/(4*4*4)*m+j];
                //radix-4
                tab2[0*N/(4*4*4*4*4)+N/(4*4*4*4)*i+N/(4*4*4)*m+j]=tmp1+tmp2+       tmp3+tmp4;
                tab2[1*N/(4*4*4*4*4)+N/(4*4*4*4)*i+N/(4*4*4)*m+j]=tmp1-tmp3+w5[0]*(tmp2-tmp4);
                tab2[2*N/(4*4*4*4*4)+N/(4*4*4*4)*i+N/(4*4*4)*m+j]=tmp1-tmp2+       tmp3-tmp4;
                tab2[3*N/(4*4*4*4*4)+N/(4*4*4*4)*i+N/(4*4*4)*m+j]=tmp1-tmp3+w6[0]*(tmp2-tmp4);
            }
        }
    }
    for(int j=0;j<N;j++)
    {
      tab[j].real() =tab2[j].real();
      tab[j].imag() =tab2[j].imag();
    }
/*

    for(int j=0;j<N;j++)
    {
      tab[j].real() =tab2[j].real()*2/N;
      tab[j].imag() =tab2[j].imag()*2/N;
    }

*/

}
///////////////////////////////////////////////