1 //
2 // pll_example.c
3 //
4 // Demonstrates a basic phase-locked loop to track the phase of a
5 // complex sinusoid.
6 //
7 
8 #include <stdio.h>
9 #include <stdlib.h>
10 #include <complex.h>
11 #include <math.h>
12 #include <time.h>
13 
14 #include "liquid.h"
15 
16 // output to octave-friendly format
17 #define OUTPUT_FILENAME "pll_example.m"
18 
main()19 int main() {
20     // parameters
21     float phase_offset      = 0.8f;     // carrier phase offset
22     float frequency_offset  = 0.01f;    // carrier frequency offset
23     float wn                = 0.05f;    // pll bandwidth
24     float zeta              = 0.707f;   // pll damping factor
25     float K                 = 1000;     // pll loop gain
26     unsigned int n          = 256;      // number of samples
27     unsigned int d          = n/32;     // print every "d" lines
28 
29     //
30     float theta[n];         // input phase
31     float complex x[n];     // input sinusoid
32     float phi[n];           // output phase
33     float complex y[n];     // output sinusoid
34 
35     // generate iir loop filter object
36     iirfilt_rrrf H = iirfilt_rrrf_create_pll(wn, zeta, K);
37     iirfilt_rrrf_print(H);
38 
39     unsigned int i;
40 
41     // generate input
42     float t=phase_offset;
43     float dt = frequency_offset;
44     for (i=0; i<n; i++) {
45         theta[i] = t;
46         x[i] = cexpf(_Complex_I*theta[i]);
47 
48         t += dt;
49     }
50 
51     // run loop
52     float phi_hat=0.0f;
53     for (i=0; i<n; i++) {
54         y[i] = cexpf(_Complex_I*phi_hat);
55 
56         // compute error
57         float e = cargf(x[i]*conjf(y[i]));
58 
59         if ( (i%d)==0 )
60             printf("e(%3u) = %12.8f;\n", i, e);
61 
62         // filter error
63         iirfilt_rrrf_execute(H,e,&phi_hat);
64 
65         phi[i] = phi_hat;
66     }
67 
68     // destroy filter
69     iirfilt_rrrf_destroy(H);
70 
71     // open output file
72     FILE * fid = fopen(OUTPUT_FILENAME,"w");
73     fprintf(fid,"clear all;\n");
74     fprintf(fid,"n=%u;\n",n);
75 
76     for (i=0; i<n; i++) {
77         fprintf(fid,"theta(%3u) = %16.8e;\n", i+1, theta[i]);
78         fprintf(fid,"  phi(%3u) = %16.8e;\n", i+1, phi[i]);
79     }
80     fprintf(fid,"t=0:(n-1);\n");
81     fprintf(fid,"figure;\n");
82     fprintf(fid,"plot(t,theta,t,phi);\n");
83     fprintf(fid,"xlabel('sample index');\n");
84     fprintf(fid,"ylabel('phase');\n");
85 
86     fclose(fid);
87 
88     printf("output written to %s.\n", OUTPUT_FILENAME);
89 
90     printf("done.\n");
91     return 0;
92 }
93