git.gag.com Git - debian/gnuradio/blob - gnuradio-core/src/lib/filter/gr_pfb_arb_resampler_ccf.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_pfb_arb_resampler_ccf.h>
  28 #include <gr_fir_ccf.h>
  29 #include <gr_fir_util.h>
  30 #include <gr_io_signature.h>
  31 #include <cstdio>
  32
  33 gr_pfb_arb_resampler_ccf_sptr gr_make_pfb_arb_resampler_ccf (float rate,
  34                                                              const std::vector<float> &taps,
  35                                                              unsigned int filter_size)
  36 {
  37   return gr_pfb_arb_resampler_ccf_sptr (new gr_pfb_arb_resampler_ccf (rate, taps,
  38                                                                       filter_size));
  39 }
  40
  41
  42 gr_pfb_arb_resampler_ccf::gr_pfb_arb_resampler_ccf (float rate,
  43                                                     const std::vector<float> &taps,
  44                                                     unsigned int filter_size)
  45   : gr_block ("pfb_arb_resampler_ccf",
  46               gr_make_io_signature (1, 1, sizeof(gr_complex)),
  47               gr_make_io_signature (1, 1, sizeof(gr_complex))),
  48     d_updated (false)
  49 {
  50   /* The number of filters is specified by the user as the filter size;
  51      this is also the interpolation rate of the filter. We use it and the
  52      rate provided to determine the decimation rate. This acts as a
  53      rational resampler. The flt_rate is calculated as the residual
  54      between the integer decimation rate and the real decimation rate and
  55      will be used to determine to interpolation point of the resampling
  56      process.
  57   */
  58   d_int_rate = filter_size;
  59   d_dec_rate = (unsigned int)floor(d_int_rate/rate);
  60   d_flt_rate = (d_int_rate/rate) - d_dec_rate;
  61
  62   // The accumulator keeps track of overflow to increment the stride correctly.
  63   d_acc = 0;
  64
  65   // Store the last filter between calls to work
  66   d_last_filter = 0;
  67
  68   d_start_index = 0;
  69
  70   d_filters = std::vector<gr_fir_ccf*>(d_int_rate);
  71
  72   // Create an FIR filter for each channel and zero out the taps
  73   std::vector<float> vtaps(0, d_int_rate);
  74   for(unsigned int i = 0; i < d_int_rate; i++) {
  75     d_filters[i] = gr_fir_util::create_gr_fir_ccf(vtaps);
  76   }
  77
  78   // Now, actually set the filters' taps
  79   set_taps(taps);
  80 }
  81
  82 gr_pfb_arb_resampler_ccf::~gr_pfb_arb_resampler_ccf ()
  83 {
  84   for(unsigned int i = 0; i < d_int_rate; i++) {
  85     delete d_filters[i];
  86   }
  87 }
  88
  89 void
  90 gr_pfb_arb_resampler_ccf::set_taps (const std::vector<float> &taps)
  91 {
  92   unsigned int i,j;
  93
  94   unsigned int ntaps = taps.size();
  95   d_taps_per_filter = (unsigned int)ceil((double)ntaps/(double)d_int_rate);
  96
  97   // Create d_numchan vectors to store each channel's taps
  98   d_taps.resize(d_int_rate);
  99
 100   // Make a vector of the taps plus fill it out with 0's to fill
 101   // each polyphase filter with exactly d_taps_per_filter
 102   std::vector<float> tmp_taps;
 103   tmp_taps = taps;
 104   while((float)(tmp_taps.size()) < d_int_rate*d_taps_per_filter) {
 105     tmp_taps.push_back(0.0);
 106   }
 107
 108   // Partition the filter
 109   for(i = 0; i < d_int_rate; i++) {
 110     // Each channel uses all d_taps_per_filter with 0's if not enough taps to fill out
 111     d_taps[i] = std::vector<float>(d_taps_per_filter, 0);
 112     for(j = 0; j < d_taps_per_filter; j++) {
 113       d_taps[i][j] = tmp_taps[i + j*d_int_rate];  // add taps to channels in reverse order
 114     }
 115
 116     // Build a filter for each channel and add it's taps to it
 117     d_filters[i]->set_taps(d_taps[i]);
 118   }
 119
 120   // Set the history to ensure enough input items for each filter
 121   set_history (d_taps_per_filter);
 122
 123   d_updated = true;
 124 }
 125
 126 void
 127 gr_pfb_arb_resampler_ccf::print_taps()
 128 {
 129   unsigned int i, j;
 130   for(i = 0; i < d_int_rate; i++) {
 131     printf("filter[%d]: [", i);
 132     for(j = 0; j < d_taps_per_filter; j++) {
 133       printf(" %.4e", d_taps[i][j]);
 134     }
 135     printf("]\n");
 136   }
 137 }
 138
 139 int
 140 gr_pfb_arb_resampler_ccf::general_work (int noutput_items,
 141                                         gr_vector_int &ninput_items,
 142                                         gr_vector_const_void_star &input_items,
 143                                         gr_vector_void_star &output_items)
 144 {
 145   gr_complex *in = (gr_complex *) input_items[0];
 146   gr_complex *out = (gr_complex *) output_items[0];
 147
 148   if (d_updated) {
 149     d_updated = false;
 150     return 0;                // history requirements may have changed.
 151   }
 152
 153   int i = 0, j, count = d_start_index;
 154   gr_complex o0, o1;
 155
 156   // Restore the last filter position
 157   j = d_last_filter;
 158
 159   // produce output as long as we can and there are enough input samples
 160   while((i < noutput_items) && (count < ninput_items[0]-1)) {
 161
 162     // start j by wrapping around mod the number of channels
 163     while((j < d_int_rate) && (i < noutput_items)) {
 164       // Take the current filter output
 165       o0 = d_filters[j]->filter(&in[count]);
 166
 167       // Take the next filter output; wrap around to 0 if necessary
 168       if(j+1 == d_int_rate)
 169         // Use the sample of the next input item through the first filter
 170         o1 = d_filters[0]->filter(&in[count+1]);
 171       else {
 172         // Use the sample from the current input item through the nex filter
 173         o1 = d_filters[j+1]->filter(&in[count]);
 174       }
 175
 176       //out[i] = o0;                     // nearest-neighbor approach
 177       out[i] = o0 + (o1 - o0)*d_acc;     // linearly interpolate between samples
 178       i++;
 179
 180       // Accumulate the position in the stream for the interpolated point.
 181       // If it goes above 1, roll around to zero and increment the stride
 182       // length this time by the decimation rate plus 1 for the increment
 183       // due to the acculated position.
 184       d_acc += d_flt_rate;
 185       j += d_dec_rate + (int)floor(d_acc);
 186       d_acc = fmodf(d_acc, 1.0);
 187     }
 188     if(i < noutput_items) {              // keep state for next entry
 189       float ss = (int)(j / d_int_rate);  // number of items to skip ahead by
 190       count += ss;                       // we have fully consumed another input
 191       j = j % d_int_rate;                // roll filter around
 192     }
 193   }
 194
 195   // Store the current filter position and start of next sample
 196   d_last_filter = j;
 197   d_start_index = std::max(0, count - ninput_items[0]);
 198
 199   // consume all we've processed but no more than we can
 200   consume_each(std::min(count, ninput_items[0]));
 201   return i;
 202 }