GNU Radio 3.7.0 C++ API
|
00001 /* -*- c++ -*- */ 00002 /* 00003 * Copyright 2002,2012 Free Software Foundation, Inc. 00004 * 00005 * This file is part of GNU Radio 00006 * 00007 * GNU Radio is free software; you can redistribute it and/or modify 00008 * it under the terms of the GNU General Public License as published by 00009 * the Free Software Foundation; either version 3, or (at your option) 00010 * any later version. 00011 * 00012 * GNU Radio is distributed in the hope that it will be useful, 00013 * but WITHOUT ANY WARRANTY; without even the implied warranty of 00014 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00015 * GNU General Public License for more details. 00016 * 00017 * You should have received a copy of the GNU General Public License 00018 * along with GNU Radio; see the file COPYING. If not, write to 00019 * the Free Software Foundation, Inc., 51 Franklin Street, 00020 * Boston, MA 02110-1301, USA. 00021 */ 00022 00023 #ifndef INCLUDED_IIR_FILTER_H 00024 #define INCLUDED_IIR_FILTER_H 00025 00026 #include <gnuradio/filter/api.h> 00027 #include <vector> 00028 #include <stdexcept> 00029 00030 namespace gr { 00031 namespace filter { 00032 namespace kernel { 00033 00034 /*! 00035 * \brief base class template for Infinite Impulse Response filter (IIR) 00036 */ 00037 template<class i_type, class o_type, class tap_type> 00038 class iir_filter 00039 { 00040 public: 00041 /*! 00042 * \brief Construct an IIR with the given taps. 00043 * 00044 * This filter uses the Direct Form I implementation, where 00045 * \p fftaps contains the feed-forward taps, and \p fbtaps the feedback ones. 00046 * 00047 * \p fftaps and \p fbtaps must have equal numbers of taps 00048 * 00049 * \p oldstyle: The old style of the IIR filter uses feedback 00050 * taps that are negative of what most definitions use (scipy 00051 * and Matlab among them). This parameter keeps using the old 00052 * GNU Radio style and is set to TRUE by default. When taps 00053 * generated from scipy, Matlab, or gr_filter_design, use the 00054 * new style by setting this to FALSE. 00055 * 00056 * The input and output satisfy a difference equation of the form 00057 00058 \f[ 00059 y[n] \pm \sum_{k=1}^{M} a_k y[n-k] = \sum_{k=0}^{N} b_k x[n-k] 00060 \f] 00061 00062 * with the corresponding rational system function 00063 00064 \f[ 00065 H(z) = \frac{\sum_{k=0}^{N} b_k z^{-k}}{1 \pm \sum_{k=1}^{M} a_k z^{-k}} 00066 \f] 00067 */ 00068 iir_filter(const std::vector<tap_type>& fftaps, 00069 const std::vector<tap_type>& fbtaps, 00070 bool oldstyle=true) 00071 throw (std::invalid_argument) 00072 { 00073 d_oldstyle = oldstyle; 00074 set_taps(fftaps, fbtaps); 00075 } 00076 00077 iir_filter() : d_latest_n(0),d_latest_m(0) { } 00078 00079 ~iir_filter() {} 00080 00081 /*! 00082 * \brief compute a single output value. 00083 * \returns the filtered input value. 00084 */ 00085 o_type filter(const i_type input); 00086 00087 /*! 00088 * \brief compute an array of N output values. 00089 * \p input must have N valid entries. 00090 */ 00091 void filter_n(o_type output[], const i_type input[], long n); 00092 00093 /*! 00094 * \return number of taps in filter. 00095 */ 00096 unsigned ntaps_ff() const { return d_fftaps.size(); } 00097 unsigned ntaps_fb() const { return d_fbtaps.size(); } 00098 00099 /*! 00100 * \brief install new taps. 00101 */ 00102 void set_taps(const std::vector<tap_type> &fftaps, 00103 const std::vector<tap_type> &fbtaps) 00104 throw (std::invalid_argument) 00105 { 00106 d_latest_n = 0; 00107 d_latest_m = 0; 00108 d_fftaps = fftaps; 00109 00110 if(d_oldstyle) { 00111 d_fbtaps = fbtaps; 00112 } 00113 else { 00114 // New style negates taps a[1:N-1] to fit with most IIR 00115 // tap generating programs. 00116 d_fbtaps.resize(fbtaps.size()); 00117 d_fbtaps[0] = fbtaps[0]; 00118 for(size_t i = 1; i < fbtaps.size(); i++) { 00119 d_fbtaps[i] = -fbtaps[i]; 00120 } 00121 } 00122 00123 int n = fftaps.size(); 00124 int m = fbtaps.size(); 00125 d_prev_input.clear(); 00126 d_prev_output.clear(); 00127 d_prev_input.resize(2 * n, 0); 00128 d_prev_output.resize(2 * m, 0); 00129 } 00130 00131 protected: 00132 bool d_oldstyle; 00133 std::vector<tap_type> d_fftaps; 00134 std::vector<tap_type> d_fbtaps; 00135 int d_latest_n; 00136 int d_latest_m; 00137 std::vector<tap_type> d_prev_output; 00138 std::vector<i_type> d_prev_input; 00139 }; 00140 00141 // 00142 // general case. We may want to specialize this 00143 // 00144 template<class i_type, class o_type, class tap_type> 00145 o_type 00146 iir_filter<i_type, o_type, tap_type>::filter(const i_type input) 00147 { 00148 tap_type acc; 00149 unsigned i = 0; 00150 unsigned n = ntaps_ff(); 00151 unsigned m = ntaps_fb(); 00152 00153 if(n == 0) 00154 return (o_type)0; 00155 00156 int latest_n = d_latest_n; 00157 int latest_m = d_latest_m; 00158 00159 acc = d_fftaps[0] * input; 00160 for(i = 1; i < n; i ++) 00161 acc += (d_fftaps[i] * d_prev_input[latest_n + i]); 00162 for(i = 1; i < m; i ++) 00163 acc += (d_fbtaps[i] * d_prev_output[latest_m + i]); 00164 00165 // store the values twice to avoid having to handle wrap-around in the loop 00166 d_prev_output[latest_m] = acc; 00167 d_prev_output[latest_m+m] = acc; 00168 d_prev_input[latest_n] = input; 00169 d_prev_input[latest_n+n] = input; 00170 00171 latest_n--; 00172 latest_m--; 00173 if(latest_n < 0) 00174 latest_n += n; 00175 if(latest_m < 0) 00176 latest_m += m; 00177 00178 d_latest_m = latest_m; 00179 d_latest_n = latest_n; 00180 return (o_type)acc; 00181 } 00182 00183 template<class i_type, class o_type, class tap_type> 00184 void 00185 iir_filter<i_type, o_type, tap_type>::filter_n(o_type output[], 00186 const i_type input[], 00187 long n) 00188 { 00189 for(int i = 0; i < n; i++) 00190 output[i] = filter(input[i]); 00191 } 00192 00193 } /* namespace kernel */ 00194 } /* namespace filter */ 00195 } /* namespace gr */ 00196 00197 #endif /* INCLUDED_IIR_FILTER_H */ 00198