shun_iwasawa a35b8f
#include "kiss_fft.h"
shun_iwasawa a35b8f
shun_iwasawa a35b8f
shun_iwasawa a35b8f
void check(kiss_fft_cpx  * in,kiss_fft_cpx  * out,int nfft,int isinverse)
shun_iwasawa a35b8f
{
shun_iwasawa a35b8f
    int bin,k;
shun_iwasawa a35b8f
    double errpow=0,sigpow=0;
shun_iwasawa a35b8f
    
shun_iwasawa a35b8f
    for (bin=0;bin
shun_iwasawa a35b8f
        double ansr = 0;
shun_iwasawa a35b8f
        double ansi = 0;
shun_iwasawa a35b8f
        double difr;
shun_iwasawa a35b8f
        double difi;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
        for (k=0;k
shun_iwasawa a35b8f
            double phase = -2*M_PI*bin*k/nfft;
shun_iwasawa a35b8f
            double re = cos(phase);
shun_iwasawa a35b8f
            double im = sin(phase);
shun_iwasawa a35b8f
            if (isinverse)
shun_iwasawa a35b8f
                im = -im;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
#ifdef FIXED_POINT
shun_iwasawa a35b8f
            re /= nfft;
shun_iwasawa a35b8f
            im /= nfft;
shun_iwasawa a35b8f
#endif            
shun_iwasawa a35b8f
shun_iwasawa a35b8f
            ansr += in[k].r * re - in[k].i * im;
shun_iwasawa a35b8f
            ansi += in[k].r * im + in[k].i * re;
shun_iwasawa a35b8f
        }
shun_iwasawa a35b8f
        difr = ansr - out[bin].r;
shun_iwasawa a35b8f
        difi = ansi - out[bin].i;
shun_iwasawa a35b8f
        errpow += difr*difr + difi*difi;
shun_iwasawa a35b8f
        sigpow += ansr*ansr+ansi*ansi;
shun_iwasawa a35b8f
    }
shun_iwasawa a35b8f
    printf("nfft=%d inverse=%d,snr = %f\n",nfft,isinverse,10*log10(sigpow/errpow) );
shun_iwasawa a35b8f
}
shun_iwasawa a35b8f
shun_iwasawa a35b8f
void test1d(int nfft,int isinverse)
shun_iwasawa a35b8f
{
shun_iwasawa a35b8f
    size_t buflen = sizeof(kiss_fft_cpx)*nfft;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    kiss_fft_cpx  * in = (kiss_fft_cpx*)malloc(buflen);
shun_iwasawa a35b8f
    kiss_fft_cpx  * out= (kiss_fft_cpx*)malloc(buflen);
shun_iwasawa a35b8f
    kiss_fft_cfg  cfg = kiss_fft_alloc(nfft,isinverse,0,0);
shun_iwasawa a35b8f
    int k;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    for (k=0;k
shun_iwasawa a35b8f
        in[k].r = (rand() % 65536) - 32768;
shun_iwasawa a35b8f
        in[k].i = (rand() % 65536) - 32768;
shun_iwasawa a35b8f
    }
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    kiss_fft(cfg,in,out);
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    check(in,out,nfft,isinverse);
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    free(in);
shun_iwasawa a35b8f
    free(out);
shun_iwasawa a35b8f
    free(cfg);
shun_iwasawa a35b8f
}
shun_iwasawa a35b8f
shun_iwasawa a35b8f
int main(int argc,char ** argv)
shun_iwasawa a35b8f
{
shun_iwasawa a35b8f
    if (argc>1) {
shun_iwasawa a35b8f
        int k;
shun_iwasawa a35b8f
        for (k=1;k
shun_iwasawa a35b8f
            test1d(atoi(argv[k]),0);
shun_iwasawa a35b8f
            test1d(atoi(argv[k]),1);
shun_iwasawa a35b8f
        }
shun_iwasawa a35b8f
    }else{
shun_iwasawa a35b8f
        test1d(32,0);
shun_iwasawa a35b8f
        test1d(32,1);
shun_iwasawa a35b8f
    }
shun_iwasawa a35b8f
    return 0;
shun_iwasawa a35b8f
}