Improved performance of FLL using a more robust error term.

2 files changed, 21 insertions, 14 deletions

diff --git a/gnuradio-core/src/lib/general/gr_fll_band_edge_cc.cc b/gnuradio-core/src/lib/general/gr_fll_band_edge_cc.cc
index 030e45ddf..7f2c468b7 100644
--- a/gnuradio-core/src/lib/general/gr_fll_band_edge_cc.cc
+++ b/gnuradio-core/src/lib/general/gr_fll_band_edge_cc.cc
@@ -53,7 +53,7 @@ gr_fll_band_edge_cc_sptr gr_make_fll_band_edge_cc (float samps_per_sym, float ro
 }
 
 
-static int ios[] = {sizeof(gr_complex), sizeof(float), sizeof(float), sizeof(float)};
+static int ios[] = {sizeof(gr_complex), sizeof(float), sizeof(float), sizeof(gr_complex)};
 static std::vector<int> iosig(ios, ios+sizeof(ios)/sizeof(int));
 gr_fll_band_edge_cc::gr_fll_band_edge_cc (float samps_per_sym, float rolloff,
 					  int filter_size, float alpha, float beta)
@@ -83,10 +83,11 @@ gr_fll_band_edge_cc::~gr_fll_band_edge_cc ()
 void
 gr_fll_band_edge_cc::set_alpha(float alpha) 
 { 
-  float eta = sqrt(2.0)/2.0;
-  float theta = alpha;
-  d_alpha = (4*eta*theta) / (1.0 + 2.0*eta*theta + theta*theta);
-  d_beta = (4*theta*theta) / (1.0 + 2.0*eta*theta + theta*theta);
+  //float eta = sqrt(2.0)/2.0;
+  //float theta = alpha;
+  //d_alpha = (4*eta*theta) / (1.0 + 2.0*eta*theta + theta*theta);
+  //d_beta = (4*theta*theta) / (1.0 + 2.0*eta*theta + theta*theta);
+  d_alpha = alpha;
 }
 
 void
@@ -160,11 +161,12 @@ gr_fll_band_edge_cc::work (int noutput_items,
   const gr_complex *in  = (const gr_complex *) input_items[0];
   gr_complex *out = (gr_complex *) output_items[0];
 
-  float *frq, *phs, *err;
+  float *frq, *phs;
+  gr_complex *err;
   if(output_items.size() > 2) {
     frq = (float *) output_items[1];
     phs = (float *) output_items[2];
-    err = (float *) output_items[3];
+    err = (gr_complex *) output_items[3];
   }
 
   if (d_updated) {
@@ -174,16 +176,17 @@ gr_fll_band_edge_cc::work (int noutput_items,
 
   int i;
   gr_complex nco_out;
-  float out_upper, out_lower;
+  gr_complex out_upper, out_lower;
   float error;
+  float avg_k = 0.1;
   for(i = 0; i < noutput_items; i++) {
     nco_out = gr_expj(d_phase);
     out[i] = in[i] * nco_out;
 
-    out_upper = norm(d_filter_upper->filter(&out[i]));
-    out_lower = norm(d_filter_lower->filter(&out[i]));
-    error = out_lower - out_upper;
-    d_error = 0.01*error + 0.99*d_error;  // average error
+    out_upper = (d_filter_upper->filter(&out[i]));
+    out_lower = (d_filter_lower->filter(&out[i]));
+    error = -real((out_upper + out_lower) * conj(out_upper - out_lower));
+    d_error = avg_k*error + avg_k*d_error;  // average error
 
     d_freq = d_freq + d_beta * d_error;
     d_phase = d_phase + d_freq + d_alpha * d_error;
diff --git a/gnuradio-core/src/lib/general/gr_fll_band_edge_cc.h b/gnuradio-core/src/lib/general/gr_fll_band_edge_cc.h
index 09baf7fde..178e18f3e 100644
--- a/gnuradio-core/src/lib/general/gr_fll_band_edge_cc.h
+++ b/gnuradio-core/src/lib/general/gr_fll_band_edge_cc.h
@@ -45,8 +45,12 @@ class gri_fft_complex;
  * (e.g., rolloff factor) of the modulated signal. The placement in frequency of the band-edges
  * is determined by the oversampling ratio (number of samples per symbol) and the excess bandwidth.
  * The size of the filters should be fairly large so as to average over a number of symbols.
- * The FLL works by calculating the power in both the upper and lower bands and comparing them. The
- * difference in power between the filters is proportional to the frequency offset.
+ *
+ * The FLL works by filtering the upper and lower band edges into x_u(t) and x_l(t), respectively.
+ * These are combined to form cc(t) = x_u(t) + x_l(t) and ss(t) = x_u(t) - x_l(t). Combining
+ * these to form the signal e(t) = Re{cc(t) \times ss(t)^*} (where ^* is the complex conjugate)
+ * provides an error signal at the DC term that is directly proportional to the carrier frequency.
+ * We then make a second-order loop using the error signal that is the running average of e(t).
  *
  * In theory, the band-edge filter is the derivative of the matched filter in frequency, 
  * (H_be(f) = \frac{H(f)}{df}. In practice, this comes down to a quarter sine wave at the point