Add fft_fc: FFT for real input signal

This commit is contained in:
Tatu Peltola 2016-10-23 19:27:37 +03:00
parent 38d567d96e
commit 5a46f15d1b
4 changed files with 100 additions and 0 deletions

91
csdr.c
View file

@ -1314,6 +1314,97 @@ int main(int argc, char *argv[])
TRY_YIELD;
}
}
if(!strcmp(argv[1],"fft_fc"))
{
/* For real FFT, the parameter is the number of output complex bins
instead of the actual FFT size.
Thus, the number of input samples used for each FFT is twice the given parameter
and for this reason, out_of_every_n_samples is also doubled
to get correct amount of overlap.
This is not very neat but makes it easier to replace fft_cc by fft_fc
in some applications. */
if(argc<=3) return badsyntax("need required parameters (fft_out_size, out_of_every_n_samples)");
int fft_in_size=0, fft_out_size=0;
sscanf(argv[2],"%d",&fft_out_size);
if(log2n(fft_out_size)==-1) return badsyntax("fft_out_size should be power of 2");
fft_in_size = 2*fft_out_size;
int every_n_samples;
sscanf(argv[3],"%d",&every_n_samples);
every_n_samples *= 2;
int benchmark=0;
int octave=0;
window_t window = WINDOW_DEFAULT;
if(argc>=5)
{
window=firdes_get_window_from_string(argv[4]);
}
if(argc>=6)
{
benchmark|=!strcmp("--benchmark",argv[5]);
octave|=!strcmp("--octave",argv[5]);
}
if(argc>=7)
{
benchmark|=!strcmp("--benchmark",argv[6]);
octave|=!strcmp("--octave",argv[6]);
}
if(!initialize_buffers()) return -2;
sendbufsize(fft_out_size);
//make FFT plan
float* input=(float*)fft_malloc(sizeof(float)*fft_in_size);
float* windowed=(float*)fft_malloc(sizeof(float)*fft_in_size);
complexf* output=(complexf*)fft_malloc(sizeof(complexf)*fft_out_size);
if(benchmark) fprintf(stderr,"fft_cc: benchmarking...");
FFT_PLAN_T* plan=make_fft_r2c(fft_in_size, windowed, output, benchmark);
if(benchmark) fprintf(stderr," done\n");
//if(octave) printf("setenv(\"GNUTERM\",\"X11 noraise\");y=zeros(1,%d);semilogy(y,\"ydatasource\",\"y\");\n",fft_size); // TODO
float *windowt;
windowt = precalculate_window(fft_in_size, window);
for(;;)
{
FEOF_CHECK;
if(every_n_samples>fft_in_size)
{
fread(input, sizeof(float), fft_in_size, stdin);
//skipping samples before next FFT (but fseek doesn't work for pipes)
for(int seek_remain=every_n_samples-fft_in_size;seek_remain>0;seek_remain-=the_bufsize)
{
fread(temp_f, sizeof(complexf), MIN_M(the_bufsize,seek_remain), stdin);
}
}
else
{
//overlapped FFT
for(int i=0;i<fft_in_size-every_n_samples;i++) input[i]=input[i+every_n_samples];
fread(input+fft_in_size-every_n_samples, sizeof(float), every_n_samples, stdin);
}
//apply_window_c(input,windowed,fft_size,window);
apply_precalculated_window_f(input,windowed,fft_in_size,windowt);
fft_execute(plan);
if(octave)
{
#if 0
// TODO
printf("fftdata=[");
//we have to swap the two parts of the array to get a valid spectrum
for(int i=fft_size/2;i<fft_size;i++) printf("(%g)+(%g)*i ",iof(output,i),qof(output,i));
for(int i=0;i<fft_size/2;i++) printf("(%g)+(%g)*i ",iof(output,i),qof(output,i));
printf(
"];\n"
"y=abs(fftdata);\n"
"refreshdata;\n"
);
#endif
}
else fwrite(output, sizeof(complexf), fft_out_size, stdout);
TRY_YIELD;
}
}
#define LOGPOWERCF_BUFSIZE 64
if(!strcmp(argv[1],"logpower_cf"))
{

View file

@ -22,6 +22,7 @@ struct fft_plan_s
#include "libcsdr.h"
FFT_PLAN_T* make_fft_c2c(int size, complexf* input, complexf* output, int forward, int benchmark);
FFT_PLAN_T* make_fft_r2c(int size, float* input, complexf* output, int benchmark);
void fft_execute(FFT_PLAN_T* plan);
void fft_destroy(FFT_PLAN_T* plan);

View file

@ -952,6 +952,13 @@ void apply_precalculated_window_c(complexf* input, complexf* output, int size, f
}
}
void apply_precalculated_window_f(float* input, float* output, int size, float *windowt)
{
for(int i=0;i<size;i++) //@apply_precalculated_window_f
{
output[i] = input[i] * windowt[i];
}
}
void apply_window_f(float* input, float* output, int size, window_t window)
{

View file

@ -138,6 +138,7 @@ void rational_resampler_get_lowpass_f(float* output, int output_size, int interp
float *precalculate_window(int size, window_t window);
void apply_window_c(complexf* input, complexf* output, int size, window_t window);
void apply_precalculated_window_c(complexf* input, complexf* output, int size, float *windowt);
void apply_precalculated_window_f(float* input, float* output, int size, float *windowt);
void apply_window_f(float* input, float* output, int size, window_t window);
void logpower_cf(complexf* input, float* output, int size, float add_db);
void accumulate_power_cf(complexf* input, float* output, int size);