//===============================================================//
void fft(float *xr, float *xi, int nr, float T)
{
	int   i, j, k, l, n, n2, nr1, i1, j1, k2, k1n2;
	float c, s, p, tr, ti, trc, tic, ars;

    n=(int)(pow(2,nr));
    n2=n/2;
    nr1=nr-1;

    k=0;
    for(l=1; l<=nr; l++)
	{
loop1:  for(j=1; j<=n2; j++)
		{
	    k2=k/(int)(pow(2,nr1));
	    p=(float)(IBIT(k2,nr));
	    ars=(float)(6.2831852*p/(float)(n));
	    c=(float)(cos(ars));
	    s=(float)(sin(ars));
            k1n2=k+n2;
            tr=xr[k1n2]*c+xi[k1n2]*s*T;
            ti=xi[k1n2]*c-xr[k1n2]*s*T;
            xr[k1n2]=xr[k]-tr;
            xi[k1n2]=xi[k]-ti;
            xr[k]=xr[k]+tr;
            xi[k]=xi[k]+ti;
            k++;		
		}
		k+=n2;
        if(k<n)
		{  goto loop1;  }
		else
		{
			k=0;
            nr1=nr1-1;
            n2 /=2;
		}
	}

    for(j=1; j<=n; j++)
	{
		i=IBIT(j-1,nr)+1;
		if(i>j)
		{
			j1=j-1;
			i1=i-1;
			trc   =xr[j1];
            tic   =xi[j1];
            xr[j1]=xr[i1];
            xi[j1]=xi[i1];
            xr[i1]=trc;
            xi[i1]=tic;
		}
	}

    if(T<0.0)
    {
		for(j=0; j<=n; j++)
		{
			xr[j]/=n;
		    xi[j]/=n;
		}
    }
}

int IBIT(int j, int m)
{
	int i, it, j1, j2;
    it=0;
    j1=j;
    for(i=1; i<=m; i++)
	{
		j2=j1/2;
		it=it*2+(j1-2*j2);
		j1=j2;
	}
    return it;
}