git.gag.com Git - debian/gnuradio/blob - gnuradio-core/src/lib/filter/gr_fft_filter_ccc.cc

   1 /* -*- c++ -*- */
   2 /*
   3  * Copyright 2005 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 /*
  24  * WARNING: This file is automatically generated by generate_gr_fft_filter_XXX.py
  25  * Any changes made to this file will be overwritten.
  26  */
  27
  28 #ifdef HAVE_CONFIG_H
  29 #include "config.h"
  30 #endif
  31
  32 #include <gr_fft_filter_ccc.h>
  33 #include <gr_io_signature.h>
  34 #include <gri_fft.h>
  35 #include <math.h>
  36 #include <assert.h>
  37 #include <stdexcept>
  38 #include <gr_firdes.h>
  39
  40 #include <iostream>
  41
  42 gr_fft_filter_ccc_sptr gr_make_fft_filter_ccc (int decimation, const std::vector<gr_complex> &taps)
  43 {
  44   return gr_fft_filter_ccc_sptr (new gr_fft_filter_ccc (decimation, taps));
  45 }
  46
  47
  48 gr_fft_filter_ccc::gr_fft_filter_ccc (int decimation, const std::vector<gr_complex> &taps)
  49   : gr_sync_decimator ("fft_filter_ccc",
  50                        gr_make_io_signature (1, 1, sizeof (gr_complex)),
  51                        gr_make_io_signature (1, 1, sizeof (gr_complex)),
  52                        decimation),
  53     d_fftsize(-1), d_fwdfft(0), d_invfft(0), d_updated(false)
  54 {
  55   // if (decimation != 1)
  56   //    throw std::invalid_argument("gr_fft_filter_ccc: decimation must be 1");
  57
  58   set_history(1);
  59   actual_set_taps(taps);
  60 }
  61
  62 gr_fft_filter_ccc::~gr_fft_filter_ccc ()
  63 {
  64   delete d_fwdfft;
  65   delete d_invfft;
  66 }
  67
  68 static void
  69 print_vector_complex(const std::string label, const std::vector<gr_complex> &x)
  70 {
  71   std::cout << label;
  72   for (unsigned i = 0; i < x.size(); i++)
  73     std::cout << x[i] << " ";
  74   std::cout << "\n";
  75 }
  76
  77 void
  78 gr_fft_filter_ccc::set_taps (const std::vector<gr_complex> &taps)
  79 {
  80   d_new_taps = taps;
  81   d_updated = true;
  82 }
  83
  84 /*
  85  * determines d_ntaps, d_nsamples, d_fftsize, d_xformed_taps
  86  */
  87 void
  88 gr_fft_filter_ccc::actual_set_taps (const std::vector<gr_complex> &taps)
  89 {
  90   int i = 0;
  91   compute_sizes(taps.size());
  92
  93   d_tail.resize(tailsize());
  94   for (i = 0; i < tailsize(); i++)
  95     d_tail[i] = 0;
  96
  97   gr_complex *in = d_fwdfft->get_inbuf();
  98   gr_complex *out = d_fwdfft->get_outbuf();
  99
 100   float scale = 1.0 / d_fftsize;
 101
 102   // Compute forward xform of taps.
 103   // Copy taps into first ntaps slots, then pad with zeros
 104   for (i = 0; i < d_ntaps; i++)
 105     in[i] = taps[i] * scale;
 106
 107   for (; i < d_fftsize; i++)
 108     in[i] = 0;
 109
 110   d_fwdfft->execute();          // do the xform
 111
 112   // now copy output to d_xformed_taps
 113   for (i = 0; i < d_fftsize; i++)
 114     d_xformed_taps[i] = out[i];
 115
 116   //print_vector_complex("transformed taps:", d_xformed_taps);
 117 }
 118
 119 // determine and set d_ntaps, d_nsamples, d_fftsize
 120
 121 void
 122 gr_fft_filter_ccc::compute_sizes(int ntaps)
 123 {
 124   int old_fftsize = d_fftsize;
 125   d_ntaps = ntaps;
 126   d_fftsize = (int) (2 * pow(2.0, ceil(log(ntaps) / log(2))));
 127   d_nsamples = d_fftsize - d_ntaps + 1;
 128
 129   if (0)
 130     fprintf(stderr, "gr_fft_filter: ntaps = %d, fftsize = %d, nsamples = %d\n",
 131             d_ntaps, d_fftsize, d_nsamples);
 132
 133   assert(d_fftsize == d_ntaps + d_nsamples -1 );
 134
 135   if (d_fftsize != old_fftsize){        // compute new plans
 136     delete d_fwdfft;
 137     delete d_invfft;
 138     d_fwdfft = new gri_fft_complex(d_fftsize, true);
 139     d_invfft = new gri_fft_complex(d_fftsize, false);
 140     d_xformed_taps.resize(d_fftsize);
 141   }
 142
 143   set_output_multiple(d_nsamples);
 144 }
 145
 146 int
 147 gr_fft_filter_ccc::work (int noutput_items,
 148                          gr_vector_const_void_star &input_items,
 149                          gr_vector_void_star &output_items)
 150 {
 151   gr_complex *in = (gr_complex *) input_items[0];
 152   gr_complex *out = (gr_complex *) output_items[0];
 153
 154   if (d_updated){
 155     actual_set_taps(d_new_taps);
 156     d_updated = false;
 157     return 0;                           // output multiple may have changed
 158   }
 159
 160   assert(noutput_items % d_nsamples == 0);
 161
 162   int dec_ctr = 0;
 163   int j = 0;
 164   int ninput_items = noutput_items * decimation();
 165
 166   for (int i = 0; i < ninput_items; i += d_nsamples){
 167
 168     memcpy(d_fwdfft->get_inbuf(), &in[i], d_nsamples * sizeof(gr_complex));
 169
 170     for (j = d_nsamples; j < d_fftsize; j++)
 171       d_fwdfft->get_inbuf()[j] = 0;
 172
 173     d_fwdfft->execute();        // compute fwd xform
 174
 175     gr_complex *a = d_fwdfft->get_outbuf();
 176     gr_complex *b = &d_xformed_taps[0];
 177     gr_complex *c = d_invfft->get_inbuf();
 178
 179     for (j = 0; j < d_fftsize; j++)     // filter in the freq domain
 180       c[j] = a[j] * b[j];
 181
 182     d_invfft->execute();        // compute inv xform
 183
 184     // add in the overlapping tail
 185
 186     for (j = 0; j < tailsize(); j++)
 187       d_invfft->get_outbuf()[j] += d_tail[j];
 188
 189     // copy nsamples to output
 190
 191     //memcpy(out, d_invfft->get_outbuf(), d_nsamples * sizeof(gr_complex));
 192     //out += d_nsamples;
 193
 194     j = dec_ctr;
 195     while (j < d_nsamples) {
 196       *out++ = d_invfft->get_outbuf()[j];
 197       j += decimation();
 198     }
 199     dec_ctr = (j - d_nsamples);
 200
 201     // stash the tail
 202     memcpy(&d_tail[0], d_invfft->get_outbuf() + d_nsamples,
 203            tailsize() * sizeof(gr_complex));
 204   }
 205
 206   assert((out - (gr_complex *) output_items[0]) == noutput_items);
 207   assert(dec_ctr == 0);
 208
 209   return noutput_items;
 210 }