diff options
| -rw-r--r-- | gr-filter/grc/CMakeLists.txt | 1 | ||||
| -rw-r--r-- | gr-filter/grc/filter_block_tree.xml | 1 | ||||
| -rw-r--r-- | gr-filter/grc/pfb_arb_resampler.xml | 61 | ||||
| -rw-r--r-- | gr-filter/grc/pfb_decimator.xml | 2 | ||||
| -rw-r--r-- | gr-filter/include/filter/CMakeLists.txt | 2 | ||||
| -rw-r--r-- | gr-filter/include/filter/pfb_arb_resampler_ccf.h | 138 | ||||
| -rw-r--r-- | gr-filter/include/filter/pfb_arb_resampler_fff.h | 139 | ||||
| -rw-r--r-- | gr-filter/include/filter/pfb_decimator_ccf.h | 2 | ||||
| -rw-r--r-- | gr-filter/lib/CMakeLists.txt | 2 | ||||
| -rw-r--r-- | gr-filter/lib/pfb_arb_resampler_ccf_impl.cc | 235 | ||||
| -rw-r--r-- | gr-filter/lib/pfb_arb_resampler_ccf_impl.h | 86 | ||||
| -rw-r--r-- | gr-filter/lib/pfb_arb_resampler_fff_impl.cc | 237 | ||||
| -rw-r--r-- | gr-filter/lib/pfb_arb_resampler_fff_impl.h | 85 | ||||
| -rwxr-xr-x | gr-filter/python/qa_pfb_arb_resampler.py | 95 | ||||
| -rw-r--r-- | gr-filter/swig/filter_swig.i | 6 |
15 files changed, 1090 insertions, 2 deletions
diff --git a/gr-filter/grc/CMakeLists.txt b/gr-filter/grc/CMakeLists.txt index 7fc129d8a..3110e0f1b 100644 --- a/gr-filter/grc/CMakeLists.txt +++ b/gr-filter/grc/CMakeLists.txt @@ -28,6 +28,7 @@ install(FILES hilbert_fc.xml iir_filter_ffd.xml interp_fir_filter_xxx.xml + pfb_arb_resampler.xml pfb_channelizer.xml pfb_decimator.xml pfb_interpolator.xml diff --git a/gr-filter/grc/filter_block_tree.xml b/gr-filter/grc/filter_block_tree.xml index 0ec9e5dd0..ecc97dd23 100644 --- a/gr-filter/grc/filter_block_tree.xml +++ b/gr-filter/grc/filter_block_tree.xml @@ -39,6 +39,7 @@ <block>hilbert_fc</block> <block>iir_filter_ffd</block> <block>interp_fir_filter_xxx</block> + <block>pfb_arb_resampler_xxx</block> <block>pfb_channelizer_ccf</block> <block>pfb_decimator_ccf</block> <block>pfb_interpolator_ccf</block> diff --git a/gr-filter/grc/pfb_arb_resampler.xml b/gr-filter/grc/pfb_arb_resampler.xml new file mode 100644 index 000000000..f3048000a --- /dev/null +++ b/gr-filter/grc/pfb_arb_resampler.xml @@ -0,0 +1,61 @@ +<?xml version="1.0"?> +<!-- +################################################### +##Polyphase Arbitrary Resampler +################################################### + --> +<block> + <name>Polyphase Arbitrary Resampler</name> + <key>pfb_arb_resampler_xxx</key> + <import>from gnuradio import filter</import> + <import>from gnuradio.filter import firdes</import> + <make>filter.pfb_arb_resampler_$(type)( + $rrate, + $taps, + $nfilts) + </make> + <callback>set_taps($taps)</callback> + <param> + <name>Type</name> + <key>type</key> + <type>enum</type> + <option> + <name>Complex->Complex (Real Taps)</name> + <key>ccf</key> + <opt>input:complex</opt> + <opt>output:complex</opt> + <opt>taps:real_vector</opt> + </option> + <option> + <name>Float->Float (Real Taps)</name> + <key>fff</key> + <opt>input:float</opt> + <opt>output:float</opt> + <opt>taps:real_vector</opt> + </option> + </param> + <param> + <name>Resampling Rate</name> + <key>rrate</key> + <type>real</type> + </param> + <param> + <name>Taps</name> + <key>taps</key> + <type>$type.taps</type> + </param> + <param> + <name>Number of Filters</name> + <key>nfilts</key> + <value>32</value> + <type>int</type> + </param> + <sink> + <name>in</name> + <type>complex</type> + </sink> + <source> + <name>out</name> + <type>complex</type> + </source> +</block> diff --git a/gr-filter/grc/pfb_decimator.xml b/gr-filter/grc/pfb_decimator.xml index 329d59c3d..b0540d3e2 100644 --- a/gr-filter/grc/pfb_decimator.xml +++ b/gr-filter/grc/pfb_decimator.xml @@ -9,7 +9,7 @@ <key>pfb_decimator_ccf</key> <import>from gnuradio import filter</import> <import>from gnuradio.filter import firdes</import> - <make>filter.pfb.decimator_ccf( + <make>filter.pfb_decimator_ccf( $decim, $taps, $channel) diff --git a/gr-filter/include/filter/CMakeLists.txt b/gr-filter/include/filter/CMakeLists.txt index 3d8e0c932..b0412aaba 100644 --- a/gr-filter/include/filter/CMakeLists.txt +++ b/gr-filter/include/filter/CMakeLists.txt @@ -100,6 +100,8 @@ install(FILES fractional_interpolator_ff.h hilbert_fc.h iir_filter_ffd.h + pfb_arb_resampler_ccf.h + pfb_arb_resampler_fff.h pfb_channelizer_ccf.h pfb_decimator_ccf.h pfb_interpolator_ccf.h diff --git a/gr-filter/include/filter/pfb_arb_resampler_ccf.h b/gr-filter/include/filter/pfb_arb_resampler_ccf.h new file mode 100644 index 000000000..1436674b2 --- /dev/null +++ b/gr-filter/include/filter/pfb_arb_resampler_ccf.h @@ -0,0 +1,138 @@ +/* -*- c++ -*- */ +/* + * Copyright 2009,2010,2012 Free Software Foundation, Inc. + * + * This file is part of GNU Radio + * + * GNU Radio is free software; you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 3, or (at your option) + * any later version. + * + * GNU Radio is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with GNU Radio; see the file COPYING. If not, write to + * the Free Software Foundation, Inc., 51 Franklin Street, + * Boston, MA 02110-1301, USA. + */ + + +#ifndef INCLUDED_PFB_ARB_RESAMPLER_CCF_H +#define INCLUDED_PFB_ARB_RESAMPLER_CCF_H + +#include <filter/api.h> +#include <gr_block.h> + +namespace gr { + namespace filter { + + /*! + * \class pfb_arb_resampler_ccf + * + * \brief Polyphase filterbank arbitrary resampler with + * gr_complex input, gr_complex output and float taps + * + * \ingroup filter_blk + * \ingroup pfb_blk + * + * This block takes in a signal stream and performs arbitrary + * resampling. The resampling rate can be any real number + * <EM>r</EM>. The resampling is done by constructing <EM>N</EM> + * filters where <EM>N</EM> is the interpolation rate. We then + * calculate <EM>D</EM> where <EM>D = floor(N/r)</EM>. + * + * Using <EM>N</EM> and <EM>D</EM>, we can perform rational + * resampling where <EM>N/D</EM> is a rational number close to the + * input rate <EM>r</EM> where we have <EM>N</EM> filters and we + * cycle through them as a polyphase filterbank with a stride of + * <EM>D</EM> so that <EM>i+1 = (i + D) % N</EM>. + * + * To get the arbitrary rate, we want to interpolate between two + * points. For each value out, we take an output from the current + * filter, <EM>i</EM>, and the next filter <EM>i+1</EM> and then + * linearly interpolate between the two based on the real + * resampling rate we want. + * + * The linear interpolation only provides us with an approximation + * to the real sampling rate specified. The error is a + * quantization error between the two filters we used as our + * interpolation points. To this end, the number of filters, + * <EM>N</EM>, used determines the quantization error; the larger + * <EM>N</EM>, the smaller the noise. You can design for a + * specified noise floor by setting the filter size (parameters + * <EM>filter_size</EM>). The size defaults to 32 filters, which + * is about as good as most implementations need. + * + * The trick with designing this filter is in how to specify the + * taps of the prototype filter. Like the PFB interpolator, the + * taps are specified using the interpolated filter rate. In this + * case, that rate is the input sample rate multiplied by the + * number of filters in the filterbank, which is also the + * interpolation rate. All other values should be relative to this + * rate. + * + * For example, for a 32-filter arbitrary resampler and using the + * GNU Radio's firdes utility to build the filter, we build a + * low-pass filter with a sampling rate of <EM>fs</EM>, a 3-dB + * bandwidth of <EM>BW</EM> and a transition bandwidth of + * <EM>TB</EM>. We can also specify the out-of-band attenuation to + * use, <EM>ATT</EM>, and the filter window function (a + * Blackman-harris window in this case). The first input is the + * gain of the filter, which we specify here as the interpolation + * rate (<EM>32</EM>). + * + * <B><EM>self._taps = gr.firdes.low_pass_2(32, 32*fs, BW, TB, + * attenuation_dB=ATT, window=gr.firdes.WIN_BLACKMAN_hARRIS)</EM></B> + * + * The theory behind this block can be found in Chapter 7.5 of + * the following book. + * + * <B><EM>f. harris, "Multirate Signal Processing for Communication + * Systems", Upper Saddle River, NJ: Prentice Hall, Inc. 2004.</EM></B> + */ + + class FILTER_API pfb_arb_resampler_ccf : virtual public gr_block + { + public: + // gr::filter::pfb_arb_resampler_ccf::sptr + typedef boost::shared_ptr<pfb_arb_resampler_ccf> sptr; + + /*! + * Build the polyphase filterbank arbitray resampler. + * \param rate (float) Specifies the resampling rate to use + * \param taps (vector/list of floats) The prototype filter to populate the filterbank. The taps + * should be generated at the filter_size sampling rate. + * \param filter_size (unsigned int) The number of filters in the filter bank. This is directly + * related to quantization noise introduced during the resampling. + * Defaults to 32 filters. + */ + static FILTER_API sptr make(float rate, + const std::vector<float> &taps, + unsigned int filter_size=32); + /*! + * Resets the filterbank's filter taps with the new prototype filter + * \param taps (vector/list of floats) The prototype filter to populate the filterbank. + */ + virtual void set_taps(const std::vector<float> &taps) = 0; + + /*! + * Return a vector<vector<>> of the filterbank taps + */ + virtual std::vector<std::vector<float> > taps() const = 0; + + /*! + * Print all of the filterbank taps to screen. + */ + virtual void print_taps() = 0; + + virtual void set_rate (float rate) = 0; + }; + + } /* namespace filter */ +} /* namespace gr */ + +#endif /* INCLUDED_PFB_ARB_RESAMPLER_CCF_H */ diff --git a/gr-filter/include/filter/pfb_arb_resampler_fff.h b/gr-filter/include/filter/pfb_arb_resampler_fff.h new file mode 100644 index 000000000..7449ea0cb --- /dev/null +++ b/gr-filter/include/filter/pfb_arb_resampler_fff.h @@ -0,0 +1,139 @@ +/* -*- c++ -*- */ +/* + * Copyright 2009-2012 Free Software Foundation, Inc. + * + * This file is part of GNU Radio + * + * GNU Radio is free software; you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 3, or (at your option) + * any later version. + * + * GNU Radio is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with GNU Radio; see the file COPYING. If not, write to + * the Free Software Foundation, Inc., 51 Franklin Street, + * Boston, MA 02110-1301, USA. + */ + + +#ifndef INCLUDED_PFB_ARB_RESAMPLER_FFF_H +#define INCLUDED_PFB_ARB_RESAMPLER_FFF_H + +#include <filter/api.h> +#include <gr_block.h> + +namespace gr { + namespace filter { + + /*! + * \class pfb_arb_resampler_fff + * + * \brief Polyphase filterbank arbitrary resampler with + * float input, float output and float taps + * + * \ingroup filter_blk + * \ingroup pfb_blk + * + * This block takes in a signal stream and performs arbitrary + * resampling. The resampling rate can be any real number + * <EM>r</EM>. The resampling is done by constructing <EM>N</EM> + * filters where <EM>N</EM> is the interpolation rate. We then + * calculate <EM>D</EM> where <EM>D = floor(N/r)</EM>. + * + * Using <EM>N</EM> and <EM>D</EM>, we can perform rational + * resampling where <EM>N/D</EM> is a rational number close to the + * input rate <EM>r</EM> where we have <EM>N</EM> filters and we + * cycle through them as a polyphase filterbank with a stride of + * <EM>D</EM> so that <EM>i+1 = (i + D) % N</EM>. + * + * To get the arbitrary rate, we want to interpolate between two + * points. For each value out, we take an output from the current + * filter, <EM>i</EM>, and the next filter <EM>i+1</EM> and then + * linearly interpolate between the two based on the real + * resampling rate we want. + * + * The linear interpolation only provides us with an approximation + * to the real sampling rate specified. The error is a + * quantization error between the two filters we used as our + * interpolation points. To this end, the number of filters, + * <EM>N</EM>, used determines the quantization error; the larger + * <EM>N</EM>, the smaller the noise. You can design for a + * specified noise floor by setting the filter size (parameters + * <EM>filter_size</EM>). The size defaults to 32 filters, which + * is about as good as most implementations need. + * + * The trick with designing this filter is in how to specify the + * taps of the prototype filter. Like the PFB interpolator, the + * taps are specified using the interpolated filter rate. In this + * case, that rate is the input sample rate multiplied by the + * number of filters in the filterbank, which is also the + * interpolation rate. All other values should be relative to this + * rate. + * + * For example, for a 32-filter arbitrary resampler and using the + * GNU Radio's firdes utility to build the filter, we build a + * low-pass filter with a sampling rate of <EM>fs</EM>, a 3-dB + * bandwidth of <EM>BW</EM> and a transition bandwidth of + * <EM>TB</EM>. We can also specify the out-of-band attenuation to + * use, <EM>ATT</EM>, and the filter window function (a + * Blackman-harris window in this case). The first input is the + * gain of the filter, which we specify here as the interpolation + * rate (<EM>32</EM>). + * + * <B><EM>self._taps = gr.firdes.low_pass_2(32, 32*fs, BW, TB, + * attenuation_dB=ATT, window=gr.firdes.WIN_BLACKMAN_hARRIS)</EM></B> + * + * The theory behind this block can be found in Chapter 7.5 of the + * following book. + * + * <B><EM>f. harris, "Multirate Signal Processing for Communication + * Systems", Upper Saddle River, NJ: Prentice Hall, Inc. 2004.</EM></B> + */ + + class FILTER_API pfb_arb_resampler_fff : virtual public gr_block + { + public: + // gr::filter::pfb_arb_resampler_fff::sptr + typedef boost::shared_ptr<pfb_arb_resampler_fff> sptr; + + /*! + * Build the polyphase filterbank arbitray resampler. + * \param rate (float) Specifies the resampling rate to use + * \param taps (vector/list of floats) The prototype filter to populate the filterbank. The taps + * should be generated at the filter_size sampling rate. + * \param filter_size (unsigned int) The number of filters in the filter bank. This is directly + * related to quantization noise introduced during the resampling. + * Defaults to 32 filters. + */ + static FILTER_API sptr make(float rate, + const std::vector<float> &taps, + unsigned int filter_size=32); + + /*! + * Resets the filterbank's filter taps with the new prototype filter + * \param taps (vector/list of floats) The prototype filter to populate the filterbank. + */ + virtual void set_taps(const std::vector<float> &taps) = 0; + + /*! + * Return a vector<vector<>> of the filterbank taps + */ + virtual std::vector<std::vector<float> > taps() const = 0; + + /*! + * Print all of the filterbank taps to screen. + */ + virtual void print_taps() = 0; + + virtual void set_rate (float rate) = 0; + }; + + } /* namespace filter */ +} /* namespace gr */ + +#endif /* INCLUDED_PFB_ARB_RESAMPLER_FFF_H */ diff --git a/gr-filter/include/filter/pfb_decimator_ccf.h b/gr-filter/include/filter/pfb_decimator_ccf.h index e0b8dce15..e41f16cd2 100644 --- a/gr-filter/include/filter/pfb_decimator_ccf.h +++ b/gr-filter/include/filter/pfb_decimator_ccf.h @@ -128,4 +128,4 @@ namespace gr { } /* namespace filter */ } /* namespace gr */ -#endif /* INCLUDED_FILTER_PFB_DECIMATOR_CCF_H */ +#endif /* INCLUDED_PFB_DECIMATOR_CCF_H */ diff --git a/gr-filter/lib/CMakeLists.txt b/gr-filter/lib/CMakeLists.txt index d1ddd2a92..8077baf9c 100644 --- a/gr-filter/lib/CMakeLists.txt +++ b/gr-filter/lib/CMakeLists.txt @@ -127,6 +127,8 @@ list(APPEND filter_sources fractional_interpolator_ff_impl.cc hilbert_fc_impl.cc iir_filter_ffd_impl.cc + pfb_arb_resampler_ccf_impl.cc + pfb_arb_resampler_fff_impl.cc pfb_channelizer_ccf_impl.cc pfb_decimator_ccf_impl.cc pfb_interpolator_ccf_impl.cc diff --git a/gr-filter/lib/pfb_arb_resampler_ccf_impl.cc b/gr-filter/lib/pfb_arb_resampler_ccf_impl.cc new file mode 100644 index 000000000..bb0906aa5 --- /dev/null +++ b/gr-filter/lib/pfb_arb_resampler_ccf_impl.cc @@ -0,0 +1,235 @@ +/* -*- c++ -*- */ +/* + * Copyright 2009,2010,2012 Free Software Foundation, Inc. + * + * This file is part of GNU Radio + * + * GNU Radio is free software; you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 3, or (at your option) + * any later version. + * + * GNU Radio is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with GNU Radio; see the file COPYING. If not, write to + * the Free Software Foundation, Inc., 51 Franklin Street, + * Boston, MA 02110-1301, USA. + */ + +#ifdef HAVE_CONFIG_H +#include "config.h" +#endif + +#include "pfb_arb_resampler_ccf_impl.h" +#include <gr_io_signature.h> +#include <cstdio> + +namespace gr { + namespace filter { + + pfb_arb_resampler_ccf::sptr + pfb_arb_resampler_ccf::make(float rate, + const std::vector<float> &taps, + unsigned int filter_size) + { + return gnuradio::get_initial_sptr + (new pfb_arb_resampler_ccf_impl(rate, taps, filter_size)); + } + + + pfb_arb_resampler_ccf_impl::pfb_arb_resampler_ccf_impl(float rate, + const std::vector<float> &taps, + unsigned int filter_size) + : gr_block("pfb_arb_resampler_ccf", + gr_make_io_signature(1, 1, sizeof(gr_complex)), + gr_make_io_signature(1, 1, sizeof(gr_complex))), + d_updated (false) + { + d_acc = 0; // start accumulator at 0 + + /* The number of filters is specified by the user as the filter + size; this is also the interpolation rate of the filter. We + use it and the rate provided to determine the decimation + rate. This acts as a rational resampler. The flt_rate is + calculated as the residual between the integer decimation + rate and the real decimation rate and will be used to + determine to interpolation point of the resampling process. + */ + d_int_rate = filter_size; + set_rate(rate); + + // Store the last filter between calls to work + d_last_filter = 0; + + d_start_index = 0; + + d_filters = std::vector<kernel::fir_filter_ccf*>(d_int_rate); + d_diff_filters = std::vector<kernel::fir_filter_ccf*>(d_int_rate); + + // Create an FIR filter for each channel and zero out the taps + std::vector<float> vtaps(0, d_int_rate); + for(unsigned int i = 0; i < d_int_rate; i++) { + d_filters[i] = new kernel::fir_filter_ccf(1, vtaps); + d_diff_filters[i] = new kernel::fir_filter_ccf(1, vtaps); + } + + // Now, actually set the filters' taps + std::vector<float> dtaps; + create_diff_taps(taps, dtaps); + create_taps(taps, d_taps, d_filters); + create_taps(dtaps, d_dtaps, d_diff_filters); + } + + pfb_arb_resampler_ccf_impl::~pfb_arb_resampler_ccf_impl() + { + for(unsigned int i = 0; i < d_int_rate; i++) { + delete d_filters[i]; + delete d_diff_filters[i]; + } + } + + void + pfb_arb_resampler_ccf_impl::create_taps(const std::vector<float> &newtaps, + std::vector< std::vector<float> > &ourtaps, + std::vector<kernel::fir_filter_ccf*> &ourfilter) + { + unsigned int ntaps = newtaps.size(); + d_taps_per_filter = (unsigned int)ceil((double)ntaps/(double)d_int_rate); + + // Create d_numchan vectors to store each channel's taps + ourtaps.resize(d_int_rate); + + // Make a vector of the taps plus fill it out with 0's to fill + // each polyphase filter with exactly d_taps_per_filter + std::vector<float> tmp_taps; + tmp_taps = newtaps; + while((float)(tmp_taps.size()) < d_int_rate*d_taps_per_filter) { + tmp_taps.push_back(0.0); + } + + // Partition the filter + for(unsigned int i = 0; i < d_int_rate; i++) { + // Each channel uses all d_taps_per_filter with 0's if not enough taps to fill out + ourtaps[d_int_rate-1-i] = std::vector<float>(d_taps_per_filter, 0); + for(unsigned int j = 0; j < d_taps_per_filter; j++) { + ourtaps[d_int_rate - 1 - i][j] = tmp_taps[i + j*d_int_rate]; + } + + // Build a filter for each channel and add it's taps to it + ourfilter[i]->set_taps(ourtaps[d_int_rate-1-i]); + } + + // Set the history to ensure enough input items for each filter + set_history (d_taps_per_filter + 1); + + d_updated = true; + } + + void + pfb_arb_resampler_ccf_impl::create_diff_taps(const std::vector<float> &newtaps, + std::vector<float> &difftaps) + { + // Calculate the differential taps (derivative filter) by taking the difference + // between two taps. Duplicate the last one to make both filters the same length. + float tap; + difftaps.clear(); + for(unsigned int i = 0; i < newtaps.size()-1; i++) { + tap = newtaps[i+1] - newtaps[i]; + difftaps.push_back(tap); + } + difftaps.push_back(tap); + } + + void + pfb_arb_resampler_ccf_impl::set_taps(const std::vector<float> &taps) + { + gruel::scoped_lock guard(d_mutex); + } + + std::vector<std::vector<float> > + pfb_arb_resampler_ccf_impl::taps() const + { + return d_taps; + } + + void + pfb_arb_resampler_ccf_impl::print_taps() + { + unsigned int i, j; + for(i = 0; i < d_int_rate; i++) { + printf("filter[%d]: [", i); + for(j = 0; j < d_taps_per_filter; j++) { + printf(" %.4e", d_taps[i][j]); + } + printf("]\n"); + } + } + + void + pfb_arb_resampler_ccf_impl::set_rate(float rate) + { + d_dec_rate = (unsigned int)floor(d_int_rate/rate); + d_flt_rate = (d_int_rate/rate) - d_dec_rate; + set_relative_rate(rate); + } + + int + pfb_arb_resampler_ccf_impl::general_work(int noutput_items, + gr_vector_int &ninput_items, + gr_vector_const_void_star &input_items, + gr_vector_void_star &output_items) + { + gr_complex *in = (gr_complex*)input_items[0]; + gr_complex *out = (gr_complex*)output_items[0]; + + if(d_updated) { + d_updated = false; + return 0; // history requirements may have changed. + } + + int i = 0, count = d_start_index; + unsigned int j; + gr_complex o0, o1; + + // Restore the last filter position + j = d_last_filter; + + // produce output as long as we can and there are enough input samples + int max_input = ninput_items[0] - (int)d_taps_per_filter; + while((i < noutput_items) && (count < max_input)) { + // start j by wrapping around mod the number of channels + while((j < d_int_rate) && (i < noutput_items)) { + // Take the current filter and derivative filter output + o0 = d_filters[j]->filter(&in[count]); + o1 = d_diff_filters[j]->filter(&in[count]); + + out[i] = o0 + o1*d_acc; // linearly interpolate between samples + i++; + + // Adjust accumulator and index into filterbank + d_acc += d_flt_rate; + j += d_dec_rate + (int)floor(d_acc); + d_acc = fmodf(d_acc, 1.0); + } + if(i < noutput_items) { // keep state for next entry + float ss = (int)(j / d_int_rate); // number of items to skip ahead by + count += ss; // we have fully consumed another input + j = j % d_int_rate; // roll filter around + } + } + + // Store the current filter position and start of next sample + d_last_filter = j; + d_start_index = std::max(0, count - ninput_items[0]); + + // consume all we've processed but no more than we can + consume_each(std::min(count, ninput_items[0])); + return i; + } + + } /* namespace filter */ +} /* namespace gr */ diff --git a/gr-filter/lib/pfb_arb_resampler_ccf_impl.h b/gr-filter/lib/pfb_arb_resampler_ccf_impl.h new file mode 100644 index 000000000..8e7e993cb --- /dev/null +++ b/gr-filter/lib/pfb_arb_resampler_ccf_impl.h @@ -0,0 +1,86 @@ +/* -*- c++ -*- */ +/* + * Copyright 2009,2010,2012 Free Software Foundation, Inc. + * + * This file is part of GNU Radio + * + * GNU Radio is free software; you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 3, or (at your option) + * any later version. + * + * GNU Radio is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with GNU Radio; see the file COPYING. If not, write to + * the Free Software Foundation, Inc., 51 Franklin Street, + * Boston, MA 02110-1301, USA. + */ + + +#ifndef INCLUDED_PFB_ARB_RESAMPLER_CCF_IMPL_H +#define INCLUDED_PFB_ARB_RESAMPLER_CCF_IMPL_H + +#include <filter/pfb_arb_resampler_ccf.h> +#include <filter/fir_filter.h> +#include <gruel/thread.h> + +namespace gr { + namespace filter { + + class FILTER_API pfb_arb_resampler_ccf_impl : public pfb_arb_resampler_ccf + { + private: + std::vector<kernel::fir_filter_ccf*> d_filters; + std::vector<kernel::fir_filter_ccf*> d_diff_filters; + std::vector< std::vector<float> > d_taps; + std::vector< std::vector<float> > d_dtaps; + unsigned int d_int_rate; // the number of filters (interpolation rate) + unsigned int d_dec_rate; // the stride through the filters (decimation rate) + float d_flt_rate; // residual rate for the linear interpolation + float d_acc; + unsigned int d_last_filter; + int d_start_index; + unsigned int d_taps_per_filter; + bool d_updated; + gruel::mutex d_mutex; // mutex to protect set/work access + + void create_diff_taps(const std::vector<float> &newtaps, + std::vector<float> &difftaps); + + /*! + * Resets the filterbank's filter taps with the new prototype filter + * \param newtaps (vector of floats) The prototype filter to populate the filterbank. + * The taps should be generated at the interpolated sampling rate. + * \param ourtaps (vector of floats) Reference to our internal member of holding the taps. + * \param ourfilter (vector of filters) Reference to our internal filter to set the taps for. + */ + void create_taps(const std::vector<float> &newtaps, + std::vector< std::vector<float> > &ourtaps, + std::vector<kernel::fir_filter_ccf*> &ourfilter); + + public: + pfb_arb_resampler_ccf_impl(float rate, + const std::vector<float> &taps, + unsigned int filter_size); + + ~pfb_arb_resampler_ccf_impl(); + + void set_taps(const std::vector<float> &taps); + std::vector<std::vector<float> > taps() const; + void print_taps(); + void set_rate(float rate); + + int general_work(int noutput_items, + gr_vector_int &ninput_items, + gr_vector_const_void_star &input_items, + gr_vector_void_star &output_items); + }; + + } /* namespace filter */ +} /* namespace gr */ + +#endif /* INCLUDED_PFB_ARB_RESAMPLER_CCF_IMPL_H */ diff --git a/gr-filter/lib/pfb_arb_resampler_fff_impl.cc b/gr-filter/lib/pfb_arb_resampler_fff_impl.cc new file mode 100644 index 000000000..79c19655a --- /dev/null +++ b/gr-filter/lib/pfb_arb_resampler_fff_impl.cc @@ -0,0 +1,237 @@ +/* -*- c++ -*- */ +/* + * Copyright 2009-2012 Free Software Foundation, Inc. + * + * This file is part of GNU Radio + * + * GNU Radio is free software; you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 3, or (at your option) + * any later version. + * + * GNU Radio is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with GNU Radio; see the file COPYING. If not, write to + * the Free Software Foundation, Inc., 51 Franklin Street, + * Boston, MA 02110-1301, USA. + */ + +#ifdef HAVE_CONFIG_H +#include "config.h" +#endif + +#include "pfb_arb_resampler_fff_impl.h" +#include <gr_io_signature.h> +#include <cstdio> + +namespace gr { + namespace filter { + + pfb_arb_resampler_fff::sptr + pfb_arb_resampler_fff::make(float rate, + const std::vector<float> &taps, + unsigned int filter_size) + { + return gnuradio::get_initial_sptr + (new pfb_arb_resampler_fff_impl(rate, taps, filter_size)); + } + + + pfb_arb_resampler_fff_impl::pfb_arb_resampler_fff_impl(float rate, + const std::vector<float> &taps, + unsigned int filter_size) + : gr_block("pfb_arb_resampler_fff", + gr_make_io_signature(1, 1, sizeof(float)), + gr_make_io_signature(1, 1, sizeof(float))), + d_updated(false) + { + d_acc = 0; // start accumulator at 0 + + /* The number of filters is specified by the user as the filter + size; this is also the interpolation rate of the filter. We + use it and the rate provided to determine the decimation + rate. This acts as a rational resampler. The flt_rate is + calculated as the residual between the integer decimation + rate and the real decimation rate and will be used to + determine to interpolation point of the resampling process. + */ + d_int_rate = filter_size; + set_rate(rate); + + // Store the last filter between calls to work + d_last_filter = 0; + + d_start_index = 0; + + d_filters = std::vector<kernel::fir_filter_fff*>(d_int_rate); + d_diff_filters = std::vector<kernel::fir_filter_fff*>(d_int_rate); + + // Create an FIR filter for each channel and zero out the taps + std::vector<float> vtaps(0, d_int_rate); + for(unsigned int i = 0; i < d_int_rate; i++) { + d_filters[i] = new kernel::fir_filter_fff(1, vtaps); + d_diff_filters[i] = new kernel::fir_filter_fff(1, vtaps); + } + + // Now, actually set the filters' taps + std::vector<float> dtaps; + create_diff_taps(taps, dtaps); + create_taps(taps, d_taps, d_filters); + create_taps(dtaps, d_dtaps, d_diff_filters); + } + + pfb_arb_resampler_fff_impl::~pfb_arb_resampler_fff_impl() + { + for(unsigned int i = 0; i < d_int_rate; i++) { + delete d_filters[i]; + delete d_diff_filters[i]; + } + } + + void + pfb_arb_resampler_fff_impl::create_taps(const std::vector<float> &newtaps, + std::vector< std::vector<float> > &ourtaps, + std::vector<kernel::fir_filter_fff*> &ourfilter) + { + unsigned int ntaps = newtaps.size(); + d_taps_per_filter = (unsigned int)ceil((double)ntaps/(double)d_int_rate); + + // Create d_numchan vectors to store each channel's taps + ourtaps.resize(d_int_rate); + + // Make a vector of the taps plus fill it out with 0's to fill + // each polyphase filter with exactly d_taps_per_filter + std::vector<float> tmp_taps; + tmp_taps = newtaps; + while((float)(tmp_taps.size()) < d_int_rate*d_taps_per_filter) { + tmp_taps.push_back(0.0); + } + + // Partition the filter + for(unsigned int i = 0; i < d_int_rate; i++) { + // Each channel uses all d_taps_per_filter with 0's if not enough taps to fill out + ourtaps[d_int_rate-1-i] = std::vector<float>(d_taps_per_filter, 0); + for(unsigned int j = 0; j < d_taps_per_filter; j++) { + ourtaps[d_int_rate - 1 - i][j] = tmp_taps[i + j*d_int_rate]; + } + + // Build a filter for each channel and add it's taps to it + ourfilter[i]->set_taps(ourtaps[d_int_rate-1-i]); + } + + // Set the history to ensure enough input items for each filter + set_history(d_taps_per_filter + 1); + + d_updated = true; + } + + void + pfb_arb_resampler_fff_impl::create_diff_taps(const std::vector<float> &newtaps, + std::vector<float> &difftaps) + { + // Calculate the differential taps (derivative filter) by taking the difference + // between two taps. Duplicate the last one to make both filters the same length. + float tap; + difftaps.clear(); + for(unsigned int i = 0; i < newtaps.size()-1; i++) { + tap = newtaps[i+1] - newtaps[i]; + difftaps.push_back(tap); + } + difftaps.push_back(tap); + } + + void + pfb_arb_resampler_fff_impl::set_taps(const std::vector<float> &taps) + { + gruel::scoped_lock guard(d_mutex); + } + + std::vector<std::vector<float> > + pfb_arb_resampler_fff_impl::taps() const + { + return d_taps; + } + + void + pfb_arb_resampler_fff_impl::print_taps() + { + unsigned int i, j; + for(i = 0; i < d_int_rate; i++) { + printf("filter[%d]: [", i); + for(j = 0; j < d_taps_per_filter; j++) { + printf(" %.4e", d_taps[i][j]); + } + printf("]\n"); + } + } + + void + pfb_arb_resampler_fff_impl::set_rate(float rate) + { + d_dec_rate = (unsigned int)floor(d_int_rate/rate); + d_flt_rate = (d_int_rate/rate) - d_dec_rate; + set_relative_rate(rate); + } + + int + pfb_arb_resampler_fff_impl::general_work(int noutput_items, + gr_vector_int &ninput_items, + gr_vector_const_void_star &input_items, + gr_vector_void_star &output_items) + { + gruel::scoped_lock guard(d_mutex); + + float *in = (float*)input_items[0]; + float *out = (float*)output_items[0]; + + if(d_updated) { + d_updated = false; + return 0; // history requirements may have changed. + } + + int i = 0, count = d_start_index; + unsigned int j; + float o0, o1; + + // Restore the last filter position + j = d_last_filter; + + // produce output as long as we can and there are enough input samples + int max_input = ninput_items[0] - (int)d_taps_per_filter; + while((i < noutput_items) && (count < max_input)) { + // start j by wrapping around mod the number of channels + while((j < d_int_rate) && (i < noutput_items)) { + // Take the current filter and derivative filter output + o0 = d_filters[j]->filter(&in[count]); + o1 = d_diff_filters[j]->filter(&in[count]); + + out[i] = o0 + o1*d_acc; // linearly interpolate between samples + i++; + + // Adjust accumulator and index into filterbank + d_acc += d_flt_rate; + j += d_dec_rate + (int)floor(d_acc); + d_acc = fmodf(d_acc, 1.0); + } + if(i < noutput_items) { // keep state for next entry + float ss = (int)(j / d_int_rate); // number of items to skip ahead by + count += ss; // we have fully consumed another input + j = j % d_int_rate; // roll filter around + } + } + + // Store the current filter position and start of next sample + d_last_filter = j; + d_start_index = std::max(0, count - ninput_items[0]); + + // consume all we've processed but no more than we can + consume_each(std::min(count, ninput_items[0])); + return i; + } + + } /* namespace filter */ +} /* namespace gr */ diff --git a/gr-filter/lib/pfb_arb_resampler_fff_impl.h b/gr-filter/lib/pfb_arb_resampler_fff_impl.h new file mode 100644 index 000000000..54e01375a --- /dev/null +++ b/gr-filter/lib/pfb_arb_resampler_fff_impl.h @@ -0,0 +1,85 @@ +/* -*- c++ -*- */ +/* + * Copyright 2009-2012 Free Software Foundation, Inc. + * + * This file is part of GNU Radio + * + * GNU Radio is free software; you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 3, or (at your option) + * any later version. + * + * GNU Radio is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with GNU Radio; see the file COPYING. If not, write to + * the Free Software Foundation, Inc., 51 Franklin Street, + * Boston, MA 02110-1301, USA. + */ + + +#ifndef INCLUDED_PFB_ARB_RESAMPLER_FFF_IMPL_H +#define INCLUDED_PFB_ARB_RESAMPLER_FFF_IMPL_H + +#include <filter/pfb_arb_resampler_fff.h> +#include <filter/fir_filter.h> +#include <gruel/thread.h> + +namespace gr { + namespace filter { + + class FILTER_API pfb_arb_resampler_fff_impl : public pfb_arb_resampler_fff + { + private: + std::vector<kernel::fir_filter_fff*> d_filters; + std::vector<kernel::fir_filter_fff*> d_diff_filters; + std::vector< std::vector<float> > d_taps; + std::vector< std::vector<float> > d_dtaps; + unsigned int d_int_rate; // the number of filters (interpolation rate) + unsigned int d_dec_rate; // the stride through the filters (decimation rate) + float d_flt_rate; // residual rate for the linear interpolation + float d_acc; + unsigned int d_last_filter; + int d_start_index; + unsigned int d_taps_per_filter; + bool d_updated; + gruel::mutex d_mutex; // mutex to protect set/work access + + void create_diff_taps(const std::vector<float> &newtaps, + std::vector<float> &difftaps); + + /*! + * Resets the filterbank's filter taps with the new prototype filter + * \param newtaps (vector of floats) The prototype filter to populate the filterbank. + * The taps should be generated at the interpolated sampling rate. + * \param ourtaps (vector of floats) Reference to our internal member of holding the taps. + * \param ourfilter (vector of filters) Reference to our internal filter to set the taps for. + */ + void create_taps(const std::vector<float> &newtaps, + std::vector< std::vector<float> > &ourtaps, + std::vector<kernel::fir_filter_fff*> &ourfilter); + public: + pfb_arb_resampler_fff_impl(float rate, + const std::vector<float> &taps, + unsigned int filter_size); + + ~pfb_arb_resampler_fff_impl(); + + void set_taps(const std::vector<float> &taps); + std::vector<std::vector<float> > taps() const; + void print_taps(); + void set_rate(float rate); + + int general_work(int noutput_items, + gr_vector_int &ninput_items, + gr_vector_const_void_star &input_items, + gr_vector_void_star &output_items); + }; + + } /* namespace filter */ +} /* namespace gr */ + +#endif /* INCLUDED_PFB_ARB_RESAMPLER_FFF_IMPL_H */ diff --git a/gr-filter/python/qa_pfb_arb_resampler.py b/gr-filter/python/qa_pfb_arb_resampler.py new file mode 100755 index 000000000..655b680f0 --- /dev/null +++ b/gr-filter/python/qa_pfb_arb_resampler.py @@ -0,0 +1,95 @@ +#!/usr/bin/env python +# +# Copyright 2012 Free Software Foundation, Inc. +# +# This file is part of GNU Radio +# +# GNU Radio is free software; you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation; either version 3, or (at your option) +# any later version. +# +# GNU Radio is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with GNU Radio; see the file COPYING. If not, write to +# the Free Software Foundation, Inc., 51 Franklin Street, +# Boston, MA 02110-1301, USA. +# + +from gnuradio import gr, gr_unittest +import filter_swig as filter +import math + +class test_pfb_arb_resampler(gr_unittest.TestCase): + + def setUp(self): + self.tb = gr.top_block() + + def tearDown(self): + self.tb = None + + def test_fff_000(self): + N = 1000 # number of samples to use + fs = 1000 # baseband sampling rate + rrate = 1.123 # resampling rate + + nfilts = 32 + taps = gr.firdes.low_pass_2(nfilts, nfilts*fs, fs/2, fs/10, + attenuation_dB=80, + window=gr.firdes.WIN_BLACKMAN_hARRIS) + + freq = 100 + signal = gr.sig_source_f(fs, gr.GR_SIN_WAVE, freq, 1) + head = gr.head(gr.sizeof_float, N) + pfb = filter.pfb_arb_resampler_fff(rrate, taps) + snk = gr.vector_sink_f() + + self.tb.connect(signal, head, pfb, snk) + self.tb.run() + + Ntest = 50 + L = len(snk.data()) + t = map(lambda x: float(x)/(fs*rrate), xrange(L)) + + phase = 0.53013 + expected_data = map(lambda x: math.sin(2.*math.pi*freq*x+phase), t) + + dst_data = snk.data() + self.assertFloatTuplesAlmostEqual(expected_data[-Ntest:], dst_data[-Ntest:], 3) + + def test_ccf_000(self): + N = 1000 # number of samples to use + fs = 1000 # baseband sampling rate + rrate = 1.123 # resampling rate + + nfilts = 32 + taps = gr.firdes.low_pass_2(nfilts, nfilts*fs, fs/2, fs/10, + attenuation_dB=80, + window=gr.firdes.WIN_BLACKMAN_hARRIS) + + freq = 100 + signal = gr.sig_source_c(fs, gr.GR_SIN_WAVE, freq, 1) + head = gr.head(gr.sizeof_gr_complex, N) + pfb = filter.pfb_arb_resampler_ccf(rrate, taps) + snk = gr.vector_sink_c() + + self.tb.connect(signal, head, pfb, snk) + self.tb.run() + + Ntest = 50 + L = len(snk.data()) + t = map(lambda x: float(x)/(fs*rrate), xrange(L)) + + phase = 0.53013 + expected_data = map(lambda x: math.cos(2.*math.pi*freq*x+phase) + \ + 1j*math.sin(2.*math.pi*freq*x+phase), t) + + dst_data = snk.data() + self.assertComplexTuplesAlmostEqual(expected_data[-Ntest:], dst_data[-Ntest:], 3) + +if __name__ == '__main__': + gr_unittest.run(test_pfb_arb_resampler, "test_pfb_arb_resampler.xml") diff --git a/gr-filter/swig/filter_swig.i b/gr-filter/swig/filter_swig.i index cc2e51d96..b72a67f15 100644 --- a/gr-filter/swig/filter_swig.i +++ b/gr-filter/swig/filter_swig.i @@ -59,6 +59,8 @@ #include "filter/interp_fir_filter_fff.h" #include "filter/interp_fir_filter_fsf.h" #include "filter/interp_fir_filter_scc.h" +#include "filter/pfb_arb_resampler_ccf.h" +#include "filter/pfb_arb_resampler_fff.h" #include "filter/pfb_channelizer_ccf.h" #include "filter/pfb_decimator_ccf.h" #include "filter/pfb_interpolator_ccf.h" @@ -97,6 +99,8 @@ %include "filter/interp_fir_filter_fff.h" %include "filter/interp_fir_filter_fsf.h" %include "filter/interp_fir_filter_scc.h" +%include "filter/pfb_arb_resampler_ccf.h" +%include "filter/pfb_arb_resampler_fff.h" %include "filter/pfb_channelizer_ccf.h" %include "filter/pfb_decimator_ccf.h" %include "filter/pfb_interpolator_ccf.h" @@ -132,6 +136,8 @@ GR_SWIG_BLOCK_MAGIC2(filter, interp_fir_filter_fcc); GR_SWIG_BLOCK_MAGIC2(filter, interp_fir_filter_fff); GR_SWIG_BLOCK_MAGIC2(filter, interp_fir_filter_fsf); GR_SWIG_BLOCK_MAGIC2(filter, interp_fir_filter_scc); +GR_SWIG_BLOCK_MAGIC2(filter, pfb_arb_resampler_ccf); +GR_SWIG_BLOCK_MAGIC2(filter, pfb_arb_resampler_fff); GR_SWIG_BLOCK_MAGIC2(filter, pfb_channelizer_ccf); GR_SWIG_BLOCK_MAGIC2(filter, pfb_decimator_ccf); GR_SWIG_BLOCK_MAGIC2(filter, pfb_interpolator_ccf); |