Working on the Costas loop

This commit is contained in:
ha7ilm 2017-05-14 17:37:54 +02:00
parent 2450ed3439
commit 2a3c210250
4 changed files with 52 additions and 25 deletions

42
csdr.c
View file

@ -136,7 +136,7 @@ char usage[]=
" psk_modulator_u8_c <n_psk>\n"
" psk31_interpolate_sine_cc <interpolation>\n"
" duplicate_samples_ntimes_u8_u8 <sample_size_bytes> <ntimes>\n"
" bpsk_costas_loop_cc <loop_bandwidth> <damping_factor> <gain> [--dd | --decision_directed]\n"
" bpsk_costas_loop_cc <loop_bandwidth> <damping_factor> [--dd | --decision_directed] [--output_error | --output_dphase | --output_nco | --octave]\n"
" binary_slicer_f_u8\n"
" simple_agc_cc <rate> [reference [max_gain]]\n"
" firdes_peak_c <rate> <length> [window [--octave]]\n"
@ -2786,7 +2786,7 @@ int main(int argc, char *argv[])
}
}
if(!strcmp(argv[1],"bpsk_costas_loop_cc")) //<loop_bandwidth> <damping_factor> <gain> [--dd | --decision_directed]
if(!strcmp(argv[1],"bpsk_costas_loop_cc")) //<loop_bandwidth> <damping_factor> [--dd | --decision_directed] [--output_error | --output_dphase | --output_nco | --octave]
{
float loop_bandwidth;
if(argc<=2) return badsyntax("need required parameter (loop_bandwidth)");
@ -2796,25 +2796,47 @@ int main(int argc, char *argv[])
if(argc<=3) return badsyntax("need required parameter (damping_factor)");
sscanf(argv[3],"%f",&damping_factor);
float gain;
if(argc<=4) return badsyntax("need required parameter (gain)");
sscanf(argv[4],"%f",&gain);
int decision_directed = !!(argc>4 && (!strcmp(argv[4], "--dd") || !strcmp(argv[5], "--decision_directed")));
int decision_directed = !!(argc>5 && (!strcmp(argv[5], "--dd") || !strcmp(argv[5], "--decision_directed")));
if(decision_directed) { errhead(); fprintf(stderr, "decision directed mode\n"); }
int output_error = !!(argc>4+decision_directed && (!strcmp(argv[4+decision_directed], "--output_error")));
int output_dphase = !output_error & !!(argc>4+decision_directed && (!strcmp(argv[4+decision_directed], "--output_dphase")));
int output_nco = !output_dphase & !!(argc>4+decision_directed && (!strcmp(argv[4+decision_directed], "--output_nco")));
int octave = !output_nco & !!(argc>4+decision_directed && (!strcmp(argv[4+decision_directed], "--octave")));
if(decision_directed && !octave) { errhead(); fprintf(stderr, "decision directed mode\n"); }
bpsk_costas_loop_state_t state;
init_bpsk_costas_loop_cc(&state, decision_directed, damping_factor, loop_bandwidth, gain);
init_bpsk_costas_loop_cc(&state, decision_directed, damping_factor, loop_bandwidth);
if(!initialize_buffers()) return -2;
sendbufsize(the_bufsize);
float* buffer_output_error = (!(octave || output_error)) ? NULL : (float*)malloc(sizeof(float)*the_bufsize);
float* buffer_output_dphase = (!(octave || output_dphase)) ? NULL : (float*)malloc(sizeof(float)*the_bufsize);
complexf* buffer_output_nco = (!(octave || output_nco)) ? NULL : (complexf*)malloc(sizeof(complexf)*the_bufsize);
if(octave) fprintf(stderr, "error=[]; dphase=[]; nco=[];\n");
for(;;)
{
FEOF_CHECK;
FREAD_C;
bpsk_costas_loop_cc((complexf*)input_buffer, (complexf*)output_buffer, the_bufsize, &state);
FWRITE_C;
bpsk_costas_loop_cc((complexf*)input_buffer, (complexf*)output_buffer, the_bufsize,
buffer_output_error, buffer_output_dphase, buffer_output_nco,
&state);
if(output_error) fwrite(buffer_output_error, sizeof(float), the_bufsize, stdout);
else if(output_dphase) fwrite(buffer_output_dphase, sizeof(float), the_bufsize, stdout);
else if(output_nco) fwrite(buffer_output_nco, sizeof(complexf), the_bufsize, stdout);
else if(octave)
{
fprintf(stderr, "error=[error ");
for(int i=0; i<the_bufsize; i++) fprintf(stderr, "%f ", buffer_output_error[i]);
fprintf(stderr, "]; dphase=[dphase ");
for(int i=0; i<the_bufsize; i++) fprintf(stderr, "%f ", buffer_output_dphase[i]);
fprintf(stderr, "]; nco=[nco ");
for(int i=0; i<the_bufsize; i++) fprintf(stderr, "(%f)+(%f)*i", buffer_output_nco[i].i, buffer_output_nco[i].q);
fprintf(stderr, "];\n");
}
else { FWRITE_C; }
TRY_YIELD;
}
}

View file

@ -270,7 +270,7 @@
</param>
<param>
<key>value</key>
<value>2**13</value>
<value>2**12</value>
</param>
</block>
<block>
@ -693,7 +693,7 @@
</param>
<param>
<key>commandline</key>
<value>csdr bpsk_costas_loop_cc $(csdr =2*pi/100) 0.707 5 2&gt;/tmp/cout</value>
<value>csdr bpsk_costas_loop_cc 0.01 0.707 --dd --octave 2&gt;/home/pcfl/Asztal/szakdoga/dipterv1/costaslog</value>
</param>
<param>
<key>comment</key>

View file

@ -2076,21 +2076,23 @@ char* timing_recovery_get_string_from_algorithm(timing_recovery_algorithm_t algo
return "INVALID";
}
void init_bpsk_costas_loop_cc(bpsk_costas_loop_state_t* s, int decision_directed, float damping_factor, float bandwidth, float gain)
void init_bpsk_costas_loop_cc(bpsk_costas_loop_state_t* s, int decision_directed, float damping_factor, float bandwidth)
{
float bandwidth_omega = 2*M_PI*bandwidth;
s->alpha = (damping_factor*2*bandwidth_omega)/gain;
float sampling_rate = 1; //the bandwidth is normalized to the sampling rate
s->beta = (bandwidth_omega*bandwidth_omega)/(sampling_rate*gain);
//based on: http://gnuradio.squarespace.com/blog/2011/8/13/control-loop-gain-values.html
float bandwidth_omega = 2*M_PI*bandwidth; //so that the bandwidth should be around 0.01 by default (2pi/100), and the damping_factor should be default 0.707
float denomiator = 1+2*damping_factor*bandwidth_omega+bandwidth_omega*bandwidth_omega;
s->alpha = (4*damping_factor*bandwidth_omega)/denomiator;
s->beta = (4*bandwidth_omega*bandwidth_omega)/denomiator;
s->iir_temp = s->dphase = s->nco_phase = 0;
}
void bpsk_costas_loop_cc(complexf* input, complexf* output, int input_size, bpsk_costas_loop_state_t* s)
void bpsk_costas_loop_cc(complexf* input, complexf* output, int input_size, float* output_error, float* output_dphase, complexf* output_nco, bpsk_costas_loop_state_t* s)
{
for(int i=0;i<input_size;i++)
{
complexf nco_sample;
e_powj(&nco_sample, s->nco_phase);
e_powj(&nco_sample, -s->nco_phase);
if(output_nco) output_nco[i]=nco_sample;
cmult(&output[i], &input[i], &nco_sample);
float error = 0;
if(s->decision_directed)
@ -2104,10 +2106,14 @@ void bpsk_costas_loop_cc(complexf* input, complexf* output, int input_size, bpsk
while(error>PI) error -= 2*PI;
}
}
else error = iof(output,i)*qof(output,i);
else error = PI*iof(output,i)*qof(output,i);
if(output_error) output_error[i]=error;
s->dphase = error * s->alpha + s->iir_temp;
s->iir_temp += error * s->beta;
//fprintf(stderr, " error = %f, dphase = %f, nco_phase = %f\n", error, s->dphase, s->nco_phase);
s->iir_temp += error * s->beta; //iir_temp could be named current_freq. See Tom Rondeau's article for better understanding.
if(s->dphase>PI) s->dphase=PI;
if(s->dphase<-PI) s->dphase=-PI;
if(output_dphase) output_dphase[i]=s->dphase;
//fprintf(stderr, " error = %f; dphase = %f; nco_phase = %f;\n", error, s->dphase, s->nco_phase);
//step NCO
s->nco_phase += s->dphase;
@ -2171,7 +2177,6 @@ void bpsk_costas_loop_c1mc(complexf* input, complexf* output, int input_size, bp
cmult(&output[i], &input[i], &vco_sample);
}
}
#endif
void simple_agc_cc(complexf* input, complexf* output, int input_size, float rate, float reference, float max_gain, float* current_gain)

View file

@ -375,7 +375,8 @@ typedef struct bpsk_costas_loop_state_s
} bpsk_costas_loop_state_t;
void plain_interpolate_cc(complexf* input, complexf* output, int input_size, int interpolation);
void bpsk_costas_loop_cc(complexf* input, complexf* output, int input_size, bpsk_costas_loop_state_t* state);
void bpsk_costas_loop_cc(complexf* input, complexf* output, int input_size, float* output_error, float* output_dphase, complexf* output_nco, bpsk_costas_loop_state_t* s);
void init_bpsk_costas_loop_cc(bpsk_costas_loop_state_t* s, int decision_directed, float damping_factor, float bandwidth);
void simple_agc_cc(complexf* input, complexf* output, int input_size, float rate, float reference, float max_gain, float* current_gain);
void firdes_add_resonator_c(complexf* output, int length, float rate, window_t window, int add, int normalize);
@ -405,5 +406,4 @@ int apply_real_fir_cc(complexf* input, complexf* output, int input_size, float*
void generic_slicer_f_u8(float* input, unsigned char* output, int input_size, int n_symbols);
void plain_interpolate_cc(complexf* input, complexf* output, int input_size, int interpolation);;
void normalize_fir_f(float* input, float* output, int length);
void init_bpsk_costas_loop_cc(bpsk_costas_loop_state_t* s, int decision_directed, float damping_factor, float bandwidth, float gain);
float* add_const_cc(complexf* input, complexf* output, int input_size, complexf x);