git.gag.com Git - debian/gnuradio/blob - gr-atsc/src/lib/GrAtscBitTimingLoop.cc

   1 /* -*- c++ -*- */
   2 /*
   3  * Copyright 2002 Free Software Foundation, Inc.
   4  *
   5  * This file is part of GNU Radio
   6  *
   7  * GNU Radio is free software; you can redistribute it and/or modify
   8  * it under the terms of the GNU General Public License as published by
   9  * the Free Software Foundation; either version 3, or (at your option)
  10  * any later version.
  11  *
  12  * GNU Radio is distributed in the hope that it will be useful,
  13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  15  * GNU General Public License for more details.
  16  *
  17  * You should have received a copy of the GNU General Public License
  18  * along with GNU Radio; see the file COPYING.  If not, write to
  19  * the Free Software Foundation, Inc., 51 Franklin Street,
  20  * Boston, MA 02110-1301, USA.
  21  */
  22
  23 #include <cmath>
  24 #include <GrAtscBitTimingLoop.h>
  25 #include "fpll_btloop_coupling.h"
  26 #include <algorithm>
  27 #include <atsc_consts.h>
  28 #include <stdio.h>
  29 #include <assert.h>
  30
  31 using std::abs;
  32
  33 static const int     DEC = 2;   // nominal decimation factor
  34
  35 /*
  36  * I strongly suggest that you not mess with these...
  37  */
  38 static const double   DEFAULT_TIMING_RATE = 2.19e-4 / FPLL_BTLOOP_COUPLING_CONST;
  39 static const double   DEFAULT_LOOP_TAP    = 0.05;
  40
  41
  42 GrAtscBitTimingLoop::GrAtscBitTimingLoop ()
  43   : VrDecimatingSigProc<float,float> (1, DEC),
  44     next_input(0), w (1.0), mu (0.5), last_right(0),
  45     debug_no_update (false)
  46 {
  47   d_timing_rate = DEFAULT_TIMING_RATE;
  48   loop.set_taps (DEFAULT_LOOP_TAP);
  49
  50   history = 1500;       // spare input samples in case we need them.
  51
  52 #ifdef _BT_DIAG_OUTPUT_
  53   fp_loop = fopen ("loop.out", "w");
  54   if (fp_loop == 0){
  55     perror ("loop.out");
  56     exit (1);
  57   }
  58
  59   fp_ps = fopen ("ps.out", "w");
  60   if (fp_ps == 0){
  61     perror ("ps.out");
  62     exit (1);
  63   }
  64 #endif
  65 }
  66
  67 //
  68 // We are nominally a 2x decimator, but our actual rate varies slightly
  69 // depending on the difference between the transmitter and receiver
  70 // sampling clocks.  Hence, we need to compute our input ranges
  71 // explictly.
  72
  73 int
  74 GrAtscBitTimingLoop::forecast(VrSampleRange output,
  75                               VrSampleRange inputs[]) {
  76   /* dec:1 ratio with history */
  77   for(unsigned int i=0;i<numberInputs;i++) {
  78     inputs[i].index=next_input;
  79     inputs[i].size=output.size*decimation + history-1;
  80   }
  81   return 0;
  82 }
  83
  84 inline double
  85 GrAtscBitTimingLoop::filter_error (double e)
  86 {
  87   static const double limit = 50 * FPLL_BTLOOP_COUPLING_CONST;
  88
  89   // first limit
  90
  91   if (e > limit)
  92     e = limit;
  93   else if (e < -limit)
  94     e = -limit;
  95
  96   return loop.filter (e);
  97 }
  98
  99 int
 100 GrAtscBitTimingLoop::work (VrSampleRange output, void *ao[],
 101                            VrSampleRange inputs[], void *ai[])
 102 {
 103   iType  *in = ((iType **)ai)[0];
 104   oType  *out = ((oType **)ao)[0];
 105
 106   // Force in-order computation of output stream.
 107   // This is required because of our slightly variable decimation factor
 108   sync (output.index);
 109
 110
 111   // We are tasked with producing output.size output samples.
 112   // We will consume approximately 2 * output.size input samples.
 113
 114
 115   unsigned int  ii = 0;         // input index
 116   unsigned int  k;              // output index
 117
 118   // We look at a window of 3 samples that we call left (oldest),
 119   // middle, right (newest).  Each time through the loop, the previous
 120   // right becomes the new left, and the new samples are middle and
 121   // right.
 122   //
 123   // The basic game plan is to drive the average difference between
 124   // right and left to zero.  Given that all transitions are
 125   // equiprobable (the data is white) and that the composite matched
 126   // filter is symmetric (raised cosine) it turns out that in the
 127   // average, if we drive that difference to zero, (implying that the
 128   // average slope at the middle point is zero), we'll be sampling
 129   // middle at the maximum or minimum point in the pulse.
 130
 131   iType left;
 132   iType middle;
 133   iType right = last_right;
 134
 135   for (k = 0; k < output.size; k++){
 136
 137     left = right;
 138     middle = produce_sample (in, ii);
 139     right = produce_sample (in, ii);
 140
 141     // assert (ii < inputs[0].size);
 142     if (!(ii < inputs[0].size)){
 143       fprintf (stderr, "ii < inputs[0].size\n");
 144       fprintf (stderr, "ii = %d, inputs[0].size = %lu, k = %d, output.size = %lu\n",
 145                ii, inputs[0].size, k, output.size);
 146       assert (0);
 147     }
 148
 149
 150     out[k] = middle;    // produce our output
 151
 152     double timing_error = -middle * ((double) right - left);
 153
 154     // update_timing_control_word
 155
 156     double filtered_timing_error = filter_error (timing_error);
 157
 158     if (!debug_no_update){
 159       mu += filtered_timing_error * d_timing_rate;
 160     }
 161
 162 #ifdef _BT_DIAG_OUTPUT_
 163     float       iodata[8];
 164     iodata[0] = left;
 165     iodata[1] = middle;
 166     iodata[2] = right;
 167     iodata[3] = timing_error;
 168     iodata[4] = filtered_timing_error;
 169     iodata[5] = mu;
 170     iodata[6] = w;
 171     iodata[7] = 0;
 172     if (fwrite (iodata, sizeof (iodata), 1, fp_loop) != 1){
 173       perror ("fwrite: loop");
 174       exit (1);
 175     }
 176 #endif
 177
 178   }
 179
 180   last_right = right;
 181   next_input += ii;     // update next_input so forecast can get us what we need
 182   return output.size;
 183 }
 184
 185 /*!
 186  * Produce samples equally spaced in time that are referenced
 187  * to the transmitter's sample clock, not ours.
 188  *
 189  * See pp 523-527 of "Digital Communication Receivers", Meyr,
 190  * Moeneclaey and Fechtel, Wiley, 1998.
 191  */
 192
 193 GrAtscBitTimingLoop::iType
 194 GrAtscBitTimingLoop::produce_sample (const iType *in, unsigned int &index)
 195 {
 196   // update mu and index as function of control word, w
 197
 198   double sum = mu + w;
 199   double f = floor (sum);
 200   int incr = (int) f;           // mostly 1, rarely 0 or 2
 201   mu = sum - f;
 202
 203   assert (0 <= incr && incr <= 2);
 204   assert (0.0 <= mu && mu <= 1.0);
 205
 206   index += incr;
 207
 208   iType n = intr.interpolate (&in[index], mu);
 209
 210 #if defined(_BT_DIAG_OUTPUT_) && 0
 211   float iodata[4];
 212   iodata[0] = incr;
 213   iodata[1] = mu;
 214   iodata[2] = w;
 215   iodata[3] = 0;
 216   if (fwrite (iodata, sizeof (iodata), 1, fp_ps) != 1){
 217     perror ("fwrite: ps");
 218     exit (1);
 219   }
 220 #endif
 221
 222   return n;
 223 }