Merge branch 'dfsg-orig'
[debian/gnuradio] / gnuradio-core / src / lib / general / gr_fll_band_edge_cc.cc
1 /* -*- c++ -*- */
2 /*
3  * Copyright 2009 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 #ifdef HAVE_CONFIG_H
24 #include "config.h"
25 #endif
26
27 #include <gr_fll_band_edge_cc.h>
28 #include <gr_fir_ccc.h>
29 #include <gr_fir_util.h>
30 #include <gri_fft.h>
31 #include <gr_io_signature.h>
32 #include <gr_expj.h>
33 #include <gr_math.h>
34 #include <cstdio>
35
36 #define M_TWOPI (2*M_PI)
37
38 float sinc(float x)
39 {
40   if(x == 0)
41     return 1;
42   else
43     return sin(M_PI*x)/(M_PI*x);
44 }
45   
46
47
48 gr_fll_band_edge_cc_sptr gr_make_fll_band_edge_cc (float samps_per_sym, float rolloff,
49                                                    int filter_size, float gain_alpha, float gain_beta)
50 {
51   return gr_fll_band_edge_cc_sptr (new gr_fll_band_edge_cc (samps_per_sym, rolloff,
52                                                             filter_size, gain_alpha, gain_beta));
53 }
54
55
56 static int ios[] = {sizeof(gr_complex), sizeof(float), sizeof(float), sizeof(gr_complex)};
57 static std::vector<int> iosig(ios, ios+sizeof(ios)/sizeof(int));
58 gr_fll_band_edge_cc::gr_fll_band_edge_cc (float samps_per_sym, float rolloff,
59                                           int filter_size, float alpha, float beta)
60   : gr_sync_block ("fll_band_edge_cc",
61                    gr_make_io_signature (1, 1, sizeof(gr_complex)),
62                    gr_make_io_signaturev (1, 4, iosig)),
63     d_alpha(alpha), d_beta(beta), d_updated (false)
64 {
65   // base this on the number of samples per symbol
66   d_max_freq =  M_TWOPI * (2.0/samps_per_sym);
67   d_min_freq = -M_TWOPI * (2.0/samps_per_sym);
68
69   d_freq = 0;
70   d_phase = 0;
71
72   set_alpha(alpha);
73
74   design_filter(samps_per_sym, rolloff, filter_size);
75 }
76
77 gr_fll_band_edge_cc::~gr_fll_band_edge_cc ()
78 {
79   delete d_filter_lower;
80   delete d_filter_upper;
81 }
82
83 void
84 gr_fll_band_edge_cc::set_alpha(float alpha) 
85
86   //float eta = sqrt(2.0)/2.0;
87   //float theta = alpha;
88   //d_alpha = (4*eta*theta) / (1.0 + 2.0*eta*theta + theta*theta);
89   //d_beta = (4*theta*theta) / (1.0 + 2.0*eta*theta + theta*theta);
90   d_alpha = alpha;
91 }
92
93 void
94 gr_fll_band_edge_cc::design_filter(float samps_per_sym, float rolloff, int filter_size)
95 {
96   int M = rint(filter_size / samps_per_sym);
97   float power = 0;
98   std::vector<float> bb_taps;
99   for(int i = 0; i < filter_size; i++) {
100     float k = -M + i*2.0/samps_per_sym;
101     float tap = sinc(rolloff*k - 0.5) + sinc(rolloff*k + 0.5);
102     power += tap;
103
104     bb_taps.push_back(tap);
105   }
106
107   int N = (bb_taps.size() - 1.0)/2.0;
108   std::vector<gr_complex> taps_lower;
109   std::vector<gr_complex> taps_upper;
110   for(unsigned int i = 0; i < bb_taps.size(); i++) {
111     float tap = bb_taps[i] / power;
112
113     float k = (-N + (int)i)/(2.0*samps_per_sym);
114     
115     gr_complex t1 = tap * gr_expj(-2*M_PI*(1+rolloff)*k);
116     gr_complex t2 = tap * gr_expj(2*M_PI*(1+rolloff)*k);
117
118     taps_lower.push_back(t1);
119     taps_upper.push_back(t2);
120   }
121
122   std::vector<gr_complex> vtaps(0, taps_lower.size());
123   d_filter_upper = gr_fir_util::create_gr_fir_ccc(vtaps);
124   d_filter_lower = gr_fir_util::create_gr_fir_ccc(vtaps);
125
126   d_filter_lower->set_taps(taps_lower);
127   d_filter_upper->set_taps(taps_upper);
128
129   d_updated = true;
130
131   // Set the history to ensure enough input items for each filter
132   set_history(filter_size+1);
133
134 }
135
136 void
137 gr_fll_band_edge_cc::print_taps()
138 {
139   unsigned int i;
140   std::vector<gr_complex> taps_upper = d_filter_upper->get_taps();
141   std::vector<gr_complex> taps_lower = d_filter_lower->get_taps();
142
143   printf("Upper Band-edge: [");
144   for(i = 0; i < taps_upper.size(); i++) {
145     printf(" %.4e + %.4ej,", taps_upper[i].real(), taps_upper[i].imag());
146   }
147   printf("]\n\n");
148
149   printf("Lower Band-edge: [");
150   for(i = 0; i < taps_lower.size(); i++) {
151     printf(" %.4e + %.4ej,", taps_lower[i].real(), taps_lower[i].imag());
152   }
153   printf("]\n\n");
154 }
155
156 int
157 gr_fll_band_edge_cc::work (int noutput_items,
158                            gr_vector_const_void_star &input_items,
159                            gr_vector_void_star &output_items)
160 {
161   const gr_complex *in  = (const gr_complex *) input_items[0];
162   gr_complex *out = (gr_complex *) output_items[0];
163
164   float *frq, *phs;
165   gr_complex *err;
166   if(output_items.size() > 2) {
167     frq = (float *) output_items[1];
168     phs = (float *) output_items[2];
169     err = (gr_complex *) output_items[3];
170   }
171
172   if (d_updated) {
173     d_updated = false;
174     return 0;                // history requirements may have changed.
175   }
176
177   int i;
178   gr_complex nco_out;
179   gr_complex out_upper, out_lower;
180   float error;
181   float avg_k = 0.1;
182   for(i = 0; i < noutput_items; i++) {
183     nco_out = gr_expj(d_phase);
184     out[i] = in[i] * nco_out;
185
186     out_upper = (d_filter_upper->filter(&out[i]));
187     out_lower = (d_filter_lower->filter(&out[i]));
188     error = -real((out_upper + out_lower) * conj(out_upper - out_lower));
189     d_error = avg_k*error + avg_k*d_error;  // average error
190
191     d_freq = d_freq + d_beta * d_error;
192     d_phase = d_phase + d_freq + d_alpha * d_error;
193
194     if(d_phase > M_PI)
195       d_phase -= M_TWOPI;
196     else if(d_phase < -M_PI)
197       d_phase += M_TWOPI;
198
199     if (d_freq > d_max_freq)
200       d_freq = d_max_freq;
201     else if (d_freq < d_min_freq)
202       d_freq = d_min_freq;
203
204     if(output_items.size() > 2) {
205       frq[i] = d_freq;
206       phs[i] = d_phase;
207       err[i] = d_error;
208     }
209   }
210
211
212   return noutput_items;
213 }