mirror of
https://github.com/mborgerding/kissfft.git
synced 2025-06-04 01:28:23 -04:00
'make snr' now tests 2d ffts
This commit is contained in:
parent
f3c4a9e9ca
commit
39d2974fe3
@ -4,6 +4,8 @@ endif
|
|||||||
ifeq "$(NUMFFTS)" ""
|
ifeq "$(NUMFFTS)" ""
|
||||||
NUMFFTS=10000
|
NUMFFTS=10000
|
||||||
endif
|
endif
|
||||||
|
NROWS=30
|
||||||
|
|
||||||
ifeq "$(DATATYPE)" ""
|
ifeq "$(DATATYPE)" ""
|
||||||
DATATYPE=float
|
DATATYPE=float
|
||||||
endif
|
endif
|
||||||
@ -51,8 +53,8 @@ mtime: all $(BENCHFFTW)
|
|||||||
done
|
done
|
||||||
|
|
||||||
snr: all
|
snr: all
|
||||||
@echo "### testing SNR for $(NFFT) point $(DATATYPE) FFTs"
|
@echo "testkiss( $(NFFT) , 1, '$(DATATYPE)' );" | octave -q
|
||||||
@echo "testkiss( $(NFFT) , '$(DATATYPE)' );" | octave -q
|
@echo "testkiss( $(NFFT) , $(NROWS), '$(DATATYPE)' );" | octave -q
|
||||||
|
|
||||||
test: snr time fftw
|
test: snr time fftw
|
||||||
|
|
||||||
|
@ -1,12 +1,19 @@
|
|||||||
function data = fileload( fname , prec , iscomplex )
|
function data = fileload( fname , prec , iscomplex ,nrows )
|
||||||
|
|
||||||
f = fopen(fname,"r", "native");
|
f = fopen(fname,"r", "native");
|
||||||
data = fread(f,Inf,prec)';
|
|
||||||
len=length(data);
|
data = fread(f,Inf,prec);
|
||||||
|
|
||||||
fclose(f);
|
fclose(f);
|
||||||
|
|
||||||
|
len = length(data);
|
||||||
if iscomplex,
|
if iscomplex,
|
||||||
data = (data(1:2:len) + j*data(2:2:len) );
|
data = (data(1:2:len) + j*data(2:2:len) );
|
||||||
|
len = len/2;
|
||||||
endif
|
endif
|
||||||
|
|
||||||
|
tdata = zeros(len/nrows,nrows);
|
||||||
|
tdata(:) = data;
|
||||||
|
data = tdata .';
|
||||||
|
|
||||||
endfunction
|
endfunction
|
||||||
|
@ -1,8 +1,12 @@
|
|||||||
function filesave( fname , prec , data )
|
function filesave( fname , prec , data )
|
||||||
|
|
||||||
f = fopen(fname,"w", "native");
|
# kiss fft stores data by rows
|
||||||
|
data = vec(data.');
|
||||||
|
|
||||||
len = length(data);
|
len = length(data);
|
||||||
|
|
||||||
|
f = fopen(fname,"w", "native");
|
||||||
|
|
||||||
if is_complex(data),
|
if is_complex(data),
|
||||||
flat=zeros( 1,2*len);
|
flat=zeros( 1,2*len);
|
||||||
flat(1:2:2*len) = real(data);
|
flat(1:2:2*len) = real(data);
|
||||||
|
@ -1,27 +1,35 @@
|
|||||||
function snr= testkiss( nfft , prec )
|
function snr= testkiss( nfft , nrows, prec )
|
||||||
|
printf('### testing SNR for %d x %d point %s FFTs\n' , nfft,nrows,prec);
|
||||||
|
|
||||||
if strcmp( prec ,'short')
|
if strcmp( prec ,'short')
|
||||||
scale_t2f=nfft;
|
scale_t2f=nfft*nrows;
|
||||||
scale_f2t=nfft;
|
scale_f2t=nfft*nrows;
|
||||||
else
|
else
|
||||||
scale_t2f=1;
|
scale_t2f=1;
|
||||||
scale_f2t=1/nfft;
|
scale_f2t=1/(nfft*nrows);
|
||||||
endif
|
endif
|
||||||
|
|
||||||
kfft= sprintf('./fftutil_%s',prec);
|
kfft= sprintf('./fftutil_%s',prec);
|
||||||
|
|
||||||
siglen = nfft;
|
sig=floor(32767*rand(nrows,nfft)) + j*floor(32767*rand(nrows,nfft));
|
||||||
sig=floor(32767*rand(1,siglen)) + j*floor(32767*rand(1,siglen));
|
|
||||||
|
|
||||||
filesave('time.dat',prec,sig);
|
filesave('time.dat',prec,sig);
|
||||||
|
if nrows > 1
|
||||||
|
cmd = sprintf('%s -r %d -n %d time.dat freq.dat',kfft,nrows,nfft);
|
||||||
|
else
|
||||||
cmd = sprintf('%s -n %d time.dat freq.dat',kfft,nfft);
|
cmd = sprintf('%s -n %d time.dat freq.dat',kfft,nfft);
|
||||||
|
endif
|
||||||
|
|
||||||
system(cmd);
|
system(cmd);
|
||||||
|
|
||||||
Fsigcomp=fileload('freq.dat',prec,1) * scale_t2f;
|
Fsigcomp=fileload('freq.dat',prec,1,nrows ) * scale_t2f;
|
||||||
Fsig=fft(sig);
|
|
||||||
|
if nrows == 1,
|
||||||
|
Fsig=fft(sig);
|
||||||
|
else
|
||||||
|
Fsig=fft2(sig);
|
||||||
|
endif
|
||||||
|
|
||||||
%x=linspace(0,2*pi*(nfft-1)/nfft,nfft);
|
|
||||||
%plot(x,abs(Fsig),'r',x,abs(Fsigcomp),'g')
|
|
||||||
diff = Fsig - Fsigcomp;
|
diff = Fsig - Fsigcomp;
|
||||||
noise_pow = sum( conj(diff).*diff );
|
noise_pow = sum( conj(diff).*diff );
|
||||||
sig_pow = sum( conj(Fsig).*Fsig );
|
sig_pow = sum( conj(Fsig).*Fsig );
|
||||||
@ -29,10 +37,14 @@ Fsig=fft(sig);
|
|||||||
avg_scale = mean( abs(Fsig) ./ abs(Fsigcomp) );
|
avg_scale = mean( abs(Fsig) ./ abs(Fsigcomp) );
|
||||||
var_scale = var( abs(Fsig) ./ abs(Fsigcomp) );
|
var_scale = var( abs(Fsig) ./ abs(Fsigcomp) );
|
||||||
|
|
||||||
|
if nrows > 1
|
||||||
|
cmd = sprintf('%s -r %d -i -n %d freq.dat time2.dat',kfft,nrows,nfft);
|
||||||
|
else
|
||||||
cmd = sprintf('%s -i -n %d freq.dat time2.dat',kfft,nfft);
|
cmd = sprintf('%s -i -n %d freq.dat time2.dat',kfft,nfft);
|
||||||
|
endif
|
||||||
|
|
||||||
system(cmd);
|
system(cmd);
|
||||||
sigcomp=fileload('time2.dat',prec,1) * scale_f2t;
|
sigcomp=fileload('time2.dat',prec,1,nrows) * scale_f2t;
|
||||||
%sig=ifft(Fsigcomp);
|
|
||||||
|
|
||||||
diff = sig - sigcomp;
|
diff = sig - sigcomp;
|
||||||
noise_pow = sum( conj(diff).*diff );
|
noise_pow = sum( conj(diff).*diff );
|
||||||
|
@ -4,6 +4,8 @@ endif
|
|||||||
ifeq "$(NUMFFTS)" ""
|
ifeq "$(NUMFFTS)" ""
|
||||||
NUMFFTS=10000
|
NUMFFTS=10000
|
||||||
endif
|
endif
|
||||||
|
NROWS=30
|
||||||
|
|
||||||
ifeq "$(DATATYPE)" ""
|
ifeq "$(DATATYPE)" ""
|
||||||
DATATYPE=float
|
DATATYPE=float
|
||||||
endif
|
endif
|
||||||
@ -51,8 +53,8 @@ mtime: all $(BENCHFFTW)
|
|||||||
done
|
done
|
||||||
|
|
||||||
snr: all
|
snr: all
|
||||||
@echo "### testing SNR for $(NFFT) point $(DATATYPE) FFTs"
|
@echo "testkiss( $(NFFT) , 1, '$(DATATYPE)' );" | octave -q
|
||||||
@echo "testkiss( $(NFFT) , '$(DATATYPE)' );" | octave -q
|
@echo "testkiss( $(NFFT) , $(NROWS), '$(DATATYPE)' );" | octave -q
|
||||||
|
|
||||||
test: snr time fftw
|
test: snr time fftw
|
||||||
|
|
||||||
|
@ -43,7 +43,11 @@ void fft_file(FILE * fin,FILE * fout,int nfft,int nrows,int isinverse,int useasc
|
|||||||
|
|
||||||
while ( fread( buf , sizeof(kiss_fft_cpx) * nfft * nrows ,1, fin ) > 0 ) {
|
while ( fread( buf , sizeof(kiss_fft_cpx) * nfft * nrows ,1, fin ) > 0 ) {
|
||||||
for (i=0;i<times;++i)
|
for (i=0;i<times;++i)
|
||||||
|
if (nrows!=1)
|
||||||
|
kiss_fft2d( st , buf ,bufout);
|
||||||
|
else
|
||||||
kiss_fft( st , buf ,bufout);
|
kiss_fft( st , buf ,bufout);
|
||||||
|
|
||||||
if (useascii) {
|
if (useascii) {
|
||||||
int i;
|
int i;
|
||||||
for (i=0;i<nfft*nrows;++i)
|
for (i=0;i<nfft*nrows;++i)
|
||||||
@ -68,12 +72,12 @@ int main(int argc,char ** argv)
|
|||||||
int nrows=1;
|
int nrows=1;
|
||||||
|
|
||||||
while (1) {
|
while (1) {
|
||||||
int c=getopt(argc,argv,"n:iax:c:");
|
int c=getopt(argc,argv,"n:iax:r:");
|
||||||
if (c==-1) break;
|
if (c==-1) break;
|
||||||
switch (c) {
|
switch (c) {
|
||||||
case 'a':useascii=1;break;
|
case 'a':useascii=1;break;
|
||||||
case 'n':nfft = atoi(optarg);break;
|
case 'n':nfft = atoi(optarg);break;
|
||||||
case 'c':nrows = atoi(optarg);break;
|
case 'r':nrows = atoi(optarg);break;
|
||||||
case 'i':isinverse=1;break;
|
case 'i':isinverse=1;break;
|
||||||
case 'x':times=atoi(optarg);break;
|
case 'x':times=atoi(optarg);break;
|
||||||
}
|
}
|
||||||
|
Loading…
x
Reference in New Issue
Block a user