Merge pull request #21 from tejeez/master

Support for receivers providing real signal
This commit is contained in:
András Retzler 2017-05-03 18:33:11 +02:00 committed by GitHub
commit 03724f3a7a
6 changed files with 174 additions and 0 deletions

137
csdr.c
View file

@ -836,6 +836,55 @@ int main(int argc, char *argv[])
return 0;
}
if(!strcmp(argv[1],"shift_addition_fc"))
{
bigbufs=1;
float starting_phase=0;
float rate;
int fd;
if(fd=init_fifo(argc,argv))
{
while(!read_fifo_ctl(fd,"%g\n",&rate)) usleep(10000);
}
else
{
if(argc<=2) return badsyntax("need required parameter (rate)");
sscanf(argv[2],"%g",&rate);
}
if(!sendbufsize(initialize_buffers())) return -2;
for(;;)
{
shift_addition_data_t data=shift_addition_init(rate);
fprintf(stderr,"shift_addition_fc: reinitialized to %g\n",rate);
int remain, current_size;
float* ibufptr;
float* obufptr;
for(;;)
{
FEOF_CHECK;
if(!FREAD_R) break;
remain=the_bufsize;
ibufptr=input_buffer;
obufptr=output_buffer;
while(remain)
{
current_size=(remain>1024)?1024:remain;
starting_phase=shift_addition_fc(ibufptr, (complexf*)obufptr, current_size, data, starting_phase);
ibufptr+=current_size;
obufptr+=current_size*2;
remain-=current_size;
}
FWRITE_C;
if(read_fifo_ctl(fd,"%g\n",&rate)) break;
TRY_YIELD;
}
}
return 0;
}
if(!strcmp(argv[1],"shift_addition_cc_test"))
{
if(argc<=2) return badsyntax("need required parameter (rate)");
@ -1495,6 +1544,94 @@ 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.
Number of input samples used for each FFT is twice the given parameter.
This 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);
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

@ -1246,6 +1246,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

@ -140,6 +140,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);

View file

@ -51,6 +51,33 @@ float shift_addition_cc(complexf *input, complexf* output, int input_size, shift
return starting_phase;
}
float shift_addition_fc(float *input, complexf* output, int input_size, shift_addition_data_t d, float starting_phase)
{
//The original idea was taken from wdsp:
//http://svn.tapr.org/repos_sdr_hpsdr/trunk/W5WC/PowerSDR_HPSDR_mRX_PS/Source/wdsp/shift.c
//However, this method introduces noise (from floating point rounding errors), which increases until the end of the buffer.
//fprintf(stderr, "cosd=%g sind=%g\n", d.cosdelta, d.sindelta);
float cosphi=cos(starting_phase);
float sinphi=sin(starting_phase);
float cosphi_last, sinphi_last;
for(int i=0;i<input_size;i++) //@shift_addition_cc: work
{
iof(output,i)=cosphi*input[i];
qof(output,i)=sinphi*input[i];
//using the trigonometric addition formulas
//cos(phi+delta)=cos(phi)cos(delta)-sin(phi)*sin(delta)
cosphi_last=cosphi;
sinphi_last=sinphi;
cosphi=cosphi_last*d.cosdelta-sinphi_last*d.sindelta;
sinphi=sinphi_last*d.cosdelta+cosphi_last*d.sindelta;
}
starting_phase+=d.rate*PI*input_size;
while(starting_phase>PI) starting_phase-=2*PI; //@shift_addition_cc: normalize starting_phase
while(starting_phase<-PI) starting_phase+=2*PI;
return starting_phase;
}
shift_addition_data_t shift_addition_init(float rate)
{
rate*=2;

View file

@ -31,6 +31,7 @@ typedef struct shift_addition_data_s
} shift_addition_data_t;
shift_addition_data_t shift_addition_init(float rate);
float shift_addition_cc(complexf *input, complexf* output, int input_size, shift_addition_data_t d, float starting_phase);
float shift_addition_fc(float *input, complexf* output, int input_size, shift_addition_data_t d, float starting_phase);
void shift_addition_cc_test(shift_addition_data_t d);
float agc_ff(float* input, float* output, int input_size, float reference, float attack_rate, float decay_rate, float max_gain, short hang_time, short attack_wait_time, float gain_filter_alpha, float last_gain);