1 //
2 // symsync_crcf_kaiser_example.c
3 //
4 // This is a simplified example of the symync family of objects to show how
5 // symbol timing can be recovered after the matched filter output.
6 //
7
8 #include <stdio.h>
9 #include <stdlib.h>
10 #include <string.h>
11 #include <math.h>
12 #include <getopt.h>
13 #include <time.h>
14 #include <assert.h>
15
16 #include "liquid.h"
17
18 #define OUTPUT_FILENAME "symsync_crcf_kaiser_example.m"
19
20 // print usage/help message
21 void usage(void);
usage(void)22 void usage(void)
23 {
24 printf("symsync_crcf_example [options]\n");
25 printf(" -h : print help\n");
26 printf(" -k <samp/sym> : filter samples/symbol, default: 2\n");
27 printf(" -m <delay> : filter delay (symbols), default: 3\n");
28 printf(" -b <beta> : filter excess bandwidth, default: 0.5\n");
29 printf(" -B <n> : filter polyphase banks, default: 32\n");
30 printf(" -s <snr> : signal-to-noise ratio, default: 30 dB\n");
31 printf(" -w <bw> : timing pll bandwidth, default: 0.02\n");
32 printf(" -n <n> : number of symbols, default: 400\n");
33 printf(" -t <tau> : timing offset, default: -0.2\n");
34 }
35
main(int argc,char * argv[])36 int main(int argc, char*argv[])
37 {
38 srand(time(NULL));
39
40 // options
41 unsigned int k = 2; // samples/symbol (input)
42 unsigned int m = 3; // filter delay (symbols)
43 float beta = 0.5f; // filter excess bandwidth factor
44 unsigned int num_filters = 32; // number of filters in the bank
45 float SNRdB = 30.0f; // signal-to-noise ratio
46 float bt = 0.02f; // loop filter bandwidth
47 unsigned int num_symbols = 400; // number of data symbols
48 float tau = -0.20f; // fractional symbol offset
49
50 // Nyquist filter type
51 liquid_firfilt_type ftype = LIQUID_FIRFILT_KAISER;
52
53 int dopt;
54 while ((dopt = getopt(argc,argv,"uhk:m:b:B:s:w:n:t:")) != EOF) {
55 switch (dopt) {
56 case 'h': usage(); return 0;
57 case 'k': k = atoi(optarg); break;
58 case 'm': m = atoi(optarg); break;
59 case 'b': beta = atof(optarg); break;
60 case 'B': num_filters = atoi(optarg); break;
61 case 's': SNRdB = atof(optarg); break;
62 case 'w': bt = atof(optarg); break;
63 case 'n': num_symbols = atoi(optarg); break;
64 case 't': tau = atof(optarg); break;
65 default:
66 exit(1);
67 }
68 }
69
70 unsigned int i;
71
72 // validate input
73 if (k < 2) {
74 fprintf(stderr,"error: k (samples/symbol) must be at least 2\n");
75 exit(1);
76 } else if (m < 1) {
77 fprintf(stderr,"error: m (filter delay) must be greater than 0\n");
78 exit(1);
79 } else if (beta <= 0.0f || beta > 1.0f) {
80 fprintf(stderr,"error: beta (excess bandwidth factor) must be in (0,1]\n");
81 exit(1);
82 } else if (num_filters == 0) {
83 fprintf(stderr,"error: number of polyphase filters must be greater than 0\n");
84 exit(1);
85 } else if (bt <= 0.0f) {
86 fprintf(stderr,"error: timing PLL bandwidth must be greater than 0\n");
87 exit(1);
88 } else if (num_symbols == 0) {
89 fprintf(stderr,"error: number of symbols must be greater than 0\n");
90 exit(1);
91 } else if (tau < -1.0f || tau > 1.0f) {
92 fprintf(stderr,"error: timing phase offset must be in [-1,1]\n");
93 exit(1);
94 }
95
96 // derived values
97 unsigned int num_samples = k*num_symbols;
98 float complex x[num_samples]; // interpolated samples
99 float complex y[num_samples]; // received signal (with noise)
100 float tau_hat[num_samples]; // instantaneous timing offset estimate
101 float complex sym_out[num_symbols + 64];// synchronized symbols
102
103 // create sequence of Nyquist-interpolated QPSK symbols
104 firinterp_crcf interp = firinterp_crcf_create_prototype(ftype,k,m,beta,tau);
105 for (i=0; i<num_symbols; i++) {
106 // generate random QPSK symbol
107 float complex s = ( rand() % 2 ? M_SQRT1_2 : -M_SQRT1_2 ) +
108 ( rand() % 2 ? M_SQRT1_2 : -M_SQRT1_2 ) * _Complex_I;
109
110 // interpolate symbol
111 firinterp_crcf_execute(interp, s, &x[i*k]);
112 }
113 firinterp_crcf_destroy(interp);
114
115
116 // add noise
117 float nstd = powf(10.0f, -SNRdB/20.0f);
118 for (i=0; i<num_samples; i++)
119 y[i] = x[i] + nstd*(randnf() + _Complex_I*randnf());
120
121
122 // create and run symbol synchronizer
123 symsync_crcf decim = symsync_crcf_create_kaiser(k, m, beta, num_filters);
124 symsync_crcf_set_lf_bw(decim,bt); // set loop filter bandwidth
125
126 // NOTE: we could just synchronize entire block (see following line);
127 // however we would like to save the instantaneous timing offset
128 // estimate for plotting purposes
129 //symsync_crcf_execute(d, y, num_samples, sym_out, &num_symbols_sync);
130
131 unsigned int num_symbols_sync = 0;
132 unsigned int num_written=0;
133 for (i=0; i<num_samples; i++) {
134 // save instantaneous timing offset estimate
135 tau_hat[i] = symsync_crcf_get_tau(decim);
136
137 // execute one sample at a time
138 symsync_crcf_execute(decim, &y[i], 1, &sym_out[num_symbols_sync], &num_written);
139
140 // increment number of symbols synchronized
141 num_symbols_sync += num_written;
142 }
143 symsync_crcf_destroy(decim);
144
145 // print last several symbols to screen
146 printf("output symbols:\n");
147 for (i=num_symbols_sync-10; i<num_symbols_sync; i++)
148 printf(" sym_out(%2u) = %8.4f + j*%8.4f;\n", i+1, crealf(sym_out[i]), cimagf(sym_out[i]));
149
150
151 //
152 // export output file
153 //
154 FILE* fid = fopen(OUTPUT_FILENAME,"w");
155 fprintf(fid,"%% %s, auto-generated file\n\n", OUTPUT_FILENAME);
156 fprintf(fid,"close all;\nclear all;\n\n");
157
158 fprintf(fid,"k=%u;\n",k);
159 fprintf(fid,"m=%u;\n",m);
160 fprintf(fid,"beta=%12.8f;\n",beta);
161 fprintf(fid,"num_filters=%u;\n",num_filters);
162 fprintf(fid,"num_symbols=%u;\n",num_symbols);
163
164 for (i=0; i<num_samples; i++)
165 fprintf(fid,"x(%3u) = %12.8f + j*%12.8f;\n", i+1, crealf(x[i]), cimagf(x[i]));
166
167 for (i=0; i<num_samples; i++)
168 fprintf(fid,"y(%3u) = %12.8f + j*%12.8f;\n", i+1, crealf(y[i]), cimagf(y[i]));
169
170 for (i=0; i<num_samples; i++)
171 fprintf(fid,"tau_hat(%3u) = %12.8f;\n", i+1, tau_hat[i]);
172
173 for (i=0; i<num_symbols_sync; i++)
174 fprintf(fid,"sym_out(%3u) = %12.8f + j*%12.8f;\n", i+1, crealf(sym_out[i]), cimagf(sym_out[i]));
175
176 fprintf(fid,"i0 = 1:round( 0.5*num_symbols );\n");
177 fprintf(fid,"i1 = round( 0.5*num_symbols ):num_symbols;\n");
178 fprintf(fid,"figure;\n");
179 fprintf(fid,"hold on;\n");
180 fprintf(fid,"plot(real(sym_out(i0)),imag(sym_out(i0)),'x','MarkerSize',4,'Color',[0.6 0.6 0.6]);\n");
181 fprintf(fid,"plot(real(sym_out(i1)),imag(sym_out(i1)),'o','MarkerSize',4,'Color',[0 0.25 0.5]);\n");
182 fprintf(fid,"hold off;\n");
183 fprintf(fid,"axis square;\n");
184 fprintf(fid,"grid on;\n");
185 fprintf(fid,"axis([-1 1 -1 1]*1.6);\n");
186 fprintf(fid,"xlabel('In-phase');\n");
187 fprintf(fid,"ylabel('Quadrature');\n");
188 fprintf(fid,"legend(['first 50%%'],['last 50%%'],1);\n");
189
190 fprintf(fid,"figure;\n");
191 fprintf(fid,"tt = 0:(length(tau_hat)-1);\n");
192 fprintf(fid,"b = floor(num_filters*tau_hat + 0.5);\n");
193 fprintf(fid,"stairs(tt,tau_hat*num_filters);\n");
194 fprintf(fid,"hold on;\n");
195 fprintf(fid,"plot(tt,b,'-k','Color',[0 0 0]);\n");
196 fprintf(fid,"hold off;\n");
197 fprintf(fid,"xlabel('time');\n");
198 fprintf(fid,"ylabel('filterbank index');\n");
199 fprintf(fid,"grid on;\n");
200 fprintf(fid,"axis([0 length(tau_hat) -1 num_filters]);\n");
201
202 fclose(fid);
203 printf("results written to %s.\n", OUTPUT_FILENAME);
204
205 // clean it up
206 printf("done.\n");
207 return 0;
208 }
209