1 //
2 // qdetector_example.c
3 //
4 // This example demonstrates the functionality of the qdetector object
5 // to detect an arbitrary signal in time in the presence of noise,
6 // carrier frequency/phase offsets, and fractional-sample timing
7 // offsets.
8 //
9
10 #include <stdio.h>
11 #include <stdlib.h>
12 #include <getopt.h>
13 #include <math.h>
14 #include <time.h>
15 #include "liquid.h"
16
17 #define OUTPUT_FILENAME "qdetector_cccf_example.m"
18
19 // print usage/help message
usage()20 void usage()
21 {
22 printf("qdetector_cccf_example\n");
23 printf("options:\n");
24 printf(" h : print usage/help\n");
25 printf(" n : number of sync symbols, default: 80\n");
26 printf(" k : samples/symbol, default: 2\n");
27 printf(" m : filter delay, default: 7 sybmols\n");
28 printf(" b : excess bandwidth factor, default: 0.3\n");
29 printf(" F : carrier frequency offset, default: -0.01\n");
30 printf(" T : fractional sample offset, default: 0\n");
31 printf(" S : SNR [dB], default: 20 dB\n");
32 printf(" t : detection threshold, default: 0.3\n");
33 printf(" r : carrier offset search range,default: 0.05\n");
34 }
35
main(int argc,char * argv[])36 int main(int argc, char*argv[])
37 {
38 // options
39 unsigned int sequence_len = 80; // number of sync symbols
40 unsigned int k = 2; // samples/symbol
41 unsigned int m = 7; // filter delay [symbols]
42 float beta = 0.3f; // excess bandwidth factor
43 int ftype = LIQUID_FIRFILT_ARKAISER;
44 float gamma = 10.0f; // channel gain
45 float tau = -0.3f; // fractional sample timing offset
46 float dphi = -0.01f; // carrier frequency offset
47 float phi = 0.5f; // carrier phase offset
48 float SNRdB = 20.0f; // signal-to-noise ratio [dB]
49 float threshold = 0.5f; // detection threshold
50 float range = 0.05f; // carrier offset search range [radians/sample]
51
52 int dopt;
53 while ((dopt = getopt(argc,argv,"hn:k:m:b:F:T:S:t:r:")) != EOF) {
54 switch (dopt) {
55 case 'h': usage(); return 0;
56 case 'n': sequence_len = atoi(optarg); break;
57 case 'k': k = atoi(optarg); break;
58 case 'm': m = atoi(optarg); break;
59 case 'b': beta = atof(optarg); break;
60 case 'F': dphi = atof(optarg); break;
61 case 'T': tau = atof(optarg); break;
62 case 'S': SNRdB = atof(optarg); break;
63 case 't': threshold = atof(optarg); break;
64 case 'r': range = atof(optarg); break;
65 default:
66 exit(1);
67 }
68 }
69
70 unsigned int i;
71
72 // validate input
73 if (tau < -0.5f || tau > 0.5f) {
74 fprintf(stderr,"error: %s, fractional sample offset must be in [-0.5,0.5]\n", argv[0]);
75 exit(1);
76 }
77
78 // generate synchronization sequence (QPSK symbols)
79 float complex sequence[sequence_len];
80 for (i=0; i<sequence_len; i++) {
81 sequence[i] = (rand() % 2 ? 1.0f : -1.0f) * M_SQRT1_2 +
82 (rand() % 2 ? 1.0f : -1.0f) * M_SQRT1_2 * _Complex_I;
83 }
84
85 //
86 float rxy = 0.0f;
87 float tau_hat = 0.0f;
88 float gamma_hat = 0.0f;
89 float dphi_hat = 0.0f;
90 float phi_hat = 0.0f;
91 int frame_detected = 0;
92
93 // create detector
94 qdetector_cccf q = qdetector_cccf_create_linear(sequence, sequence_len, ftype, k, m, beta);
95 qdetector_cccf_set_threshold(q, threshold);
96 qdetector_cccf_set_range (q, range);
97 qdetector_cccf_print(q);
98
99 //
100 unsigned int seq_len = qdetector_cccf_get_seq_len(q);
101 unsigned int buf_len = qdetector_cccf_get_buf_len(q);
102 unsigned int num_samples = 2*buf_len; // double buffer length to ensure detection
103 unsigned int num_symbols = buf_len;
104
105 // arrays
106 float complex y[num_samples]; // received signal
107 float complex syms_rx[num_symbols]; // recovered symbols
108
109 // get pointer to sequence and generate full sequence
110 float complex * v = (float complex*) qdetector_cccf_get_sequence(q);
111 unsigned int filter_delay = 15;
112 firfilt_crcf filter = firfilt_crcf_create_kaiser(2*filter_delay+1, 0.4f, 60.0f, -tau);
113 float nstd = powf(10.0f, -SNRdB/20.0f);
114 for (i=0; i<num_samples; i++) {
115 // add delay
116 firfilt_crcf_push(filter, i < seq_len ? v[i] : 0);
117 firfilt_crcf_execute(filter, &y[i]);
118
119 // channel gain
120 y[i] *= gamma;
121
122 // carrier offset
123 y[i] *= cexpf(_Complex_I*(dphi*i + phi));
124
125 // noise
126 y[i] += nstd*(randnf() + _Complex_I*randnf())*M_SQRT1_2;
127 }
128 firfilt_crcf_destroy(filter);
129
130 // run detection on sequence
131 for (i=0; i<num_samples; i++) {
132 v = qdetector_cccf_execute(q,y[i]);
133
134 if (v != NULL) {
135 printf("\nframe detected!\n");
136 frame_detected = 1;
137
138 // get statistics
139 rxy = qdetector_cccf_get_rxy(q);
140 tau_hat = qdetector_cccf_get_tau(q);
141 gamma_hat = qdetector_cccf_get_gamma(q);
142 dphi_hat = qdetector_cccf_get_dphi(q);
143 phi_hat = qdetector_cccf_get_phi(q);
144 break;
145 }
146 }
147
148 unsigned int num_syms_rx = 0; // output symbol counter
149 unsigned int counter = 0; // decimation counter
150 if (frame_detected) {
151 // recover symbols
152 firfilt_crcf mf = firfilt_crcf_create_rnyquist(ftype, k, m, beta, tau_hat);
153 firfilt_crcf_set_scale(mf, 1.0f / (float)(k*gamma_hat));
154 nco_crcf nco = nco_crcf_create(LIQUID_VCO);
155 nco_crcf_set_frequency(nco, dphi_hat);
156 nco_crcf_set_phase (nco, phi_hat);
157
158 for (i=0; i<buf_len; i++) {
159 //
160 float complex sample;
161 nco_crcf_mix_down(nco, v[i], &sample);
162 nco_crcf_step(nco);
163
164 // apply decimator
165 firfilt_crcf_push(mf, sample);
166 counter++;
167 if (counter == k-1)
168 firfilt_crcf_execute(mf, &syms_rx[num_syms_rx++]);
169 counter %= k;
170 }
171
172 nco_crcf_destroy(nco);
173 firfilt_crcf_destroy(mf);
174 }
175
176 // destroy objects
177 qdetector_cccf_destroy(q);
178
179 // print results
180 printf("\n");
181 printf("frame detected : %s\n", frame_detected ? "yes" : "no");
182 if (frame_detected) {
183 printf(" rxy : %8.3f\n", rxy);
184 printf(" gamma hat : %8.3f, actual=%8.3f (error=%8.3f)\n", gamma_hat, gamma, gamma_hat - gamma);
185 printf(" tau hat : %8.3f, actual=%8.3f (error=%8.3f) samples\n", tau_hat, tau, tau_hat - tau );
186 printf(" dphi hat : %8.5f, actual=%8.5f (error=%8.5f) rad/sample\n", dphi_hat, dphi, dphi_hat - dphi );
187 printf(" phi hat : %8.5f, actual=%8.5f (error=%8.5f) radians\n", phi_hat, phi, phi_hat - phi );
188 printf(" symbols rx : %u\n", num_syms_rx);
189 }
190 printf("\n");
191
192 //
193 // export results
194 //
195 FILE * fid = fopen(OUTPUT_FILENAME,"w");
196 fprintf(fid,"%% %s : auto-generated file\n", OUTPUT_FILENAME);
197 fprintf(fid,"clear all\n");
198 fprintf(fid,"close all\n");
199 fprintf(fid,"sequence_len= %u;\n", sequence_len);
200 fprintf(fid,"num_samples = %u;\n", num_samples);
201
202 fprintf(fid,"y = zeros(1,num_samples);\n");
203 for (i=0; i<num_samples; i++)
204 fprintf(fid,"y(%4u) = %12.8f + j*%12.8f;\n", i+1, crealf(y[i]), cimagf(y[i]));
205
206 fprintf(fid,"num_syms_rx = %u;\n", num_syms_rx);
207 fprintf(fid,"syms_rx = zeros(1,num_syms_rx);\n");
208 for (i=0; i<num_syms_rx; i++)
209 fprintf(fid,"syms_rx(%4u) = %12.8f + j*%12.8f;\n", i+1, crealf(syms_rx[i]), cimagf(syms_rx[i]));
210
211 fprintf(fid,"t=[0:(num_samples-1)];\n");
212 fprintf(fid,"figure;\n");
213 fprintf(fid,"subplot(4,1,1);\n");
214 fprintf(fid," plot(t,real(y), t,imag(y));\n");
215 fprintf(fid," grid on;\n");
216 fprintf(fid," xlabel('time');\n");
217 fprintf(fid," ylabel('received signal');\n");
218 fprintf(fid,"subplot(4,1,2:4);\n");
219 fprintf(fid," plot(real(syms_rx), imag(syms_rx), 'x');\n");
220 fprintf(fid," axis([-1 1 -1 1]*1.5);\n");
221 fprintf(fid," axis square;\n");
222 fprintf(fid," grid on;\n");
223 fprintf(fid," xlabel('real');\n");
224 fprintf(fid," ylabel('imag');\n");
225
226 fclose(fid);
227 printf("results written to '%s'\n", OUTPUT_FILENAME);
228
229 return 0;
230 }
231