1 /* Calf DSP Library 2 * Basic one-pole one-zero filters based on bilinear transform. 3 * Copyright (C) 2001-2007 Krzysztof Foltman 4 * 5 * This program is free software; you can redistribute it and/or 6 * modify it under the terms of the GNU Lesser General Public 7 * License as published by the Free Software Foundation; either 8 * version 2 of the License, or (at your option) any later version. 9 * 10 * This program is distributed in the hope that it will be useful, 11 * but WITHOUT ANY WARRANTY; without even the implied warranty of 12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 13 * Lesser General Public License for more details. 14 * 15 * You should have received a copy of the GNU Lesser General 16 * Public License along with this program; if not, write to the 17 * Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, 18 * Boston, MA 02111-1307, USA. 19 */ 20 #ifndef __CALF_ONEPOLE_H 21 #define __CALF_ONEPOLE_H 22 23 #include "primitives.h" 24 25 namespace dsp { 26 27 /** 28 * one-pole filter, for floating point values 29 * coefficient calculation is based on bilinear transform, and the code itself is based on my very old OneSignal lib 30 * lp and hp are *somewhat* tested, allpass is not tested at all 31 * don't use this for integers because it won't work 32 */ 33 template<class T = float, class Coeff = float> 34 class onepole 35 { 36 public: 37 typedef std::complex<double> cfloat; 38 39 T x1, y1; 40 Coeff a0, a1, b1; 41 onepole()42 onepole() 43 { 44 reset(); 45 } 46 47 /// Set coefficients for a lowpass filter set_lp(float fc,float sr)48 void set_lp(float fc, float sr) 49 { 50 // x x 51 // x+1 x-1 52 Coeff x = tan (M_PI * fc / (2 * sr)); 53 Coeff q = 1/(1+x); 54 a0 = a1 = x*q; 55 b1 = (x-1)*q; 56 } 57 58 /// Set coefficients for an allpass filter set_ap(float fc,float sr)59 void set_ap(float fc, float sr) 60 { 61 // x-1 x+1 62 // x+1 x-1 63 Coeff x = tan (M_PI * fc / (2 * sr)); 64 Coeff q = 1/(1+x); 65 b1 = a0 = (x-1)*q; 66 a1 = 1; 67 } 68 69 /// Set coefficients for an allpass filter, using omega instead of fc and sr 70 /// omega = (PI / 2) * fc / sr set_ap_w(float w)71 void set_ap_w(float w) 72 { 73 // x-1 x+1 74 // x+1 x-1 75 Coeff x = tan (w); 76 Coeff q = 1/(1+x); 77 b1 = a0 = (x-1)*q; 78 a1 = 1; 79 } 80 81 /// Set coefficients for a highpass filter set_hp(float fc,float sr)82 void set_hp(float fc, float sr) 83 { 84 // x -x 85 // x+1 x-1 86 Coeff x = tan (M_PI * fc / (2 * sr)); 87 Coeff q = 1/(1+x); 88 a0 = q; 89 a1 = -a0; 90 b1 = (x-1)*q; 91 } 92 93 /// Process one sample process(T in)94 inline T process(T in) 95 { 96 T out = in * a0 + x1 * a1 - y1 * b1; 97 x1 = in; 98 y1 = out; 99 return out; 100 } 101 102 /// Process one sample, assuming it's a lowpass filter (optimized special case) process_lp(T in)103 inline T process_lp(T in) 104 { 105 T out = (in + x1) * a0 - y1 * b1; 106 x1 = in; 107 y1 = out; 108 return out; 109 } 110 111 /// Process one sample, assuming it's a highpass filter (optimized special case) process_hp(T in)112 inline T process_hp(T in) 113 { 114 T out = (in - x1) * a0 - y1 * b1; 115 x1 = in; 116 y1 = out; 117 return out; 118 } 119 120 /// Process one sample, assuming it's an allpass filter (optimized special case) process_ap(T in)121 inline T process_ap(T in) 122 { 123 T out = (in - y1) * a0 + x1; 124 x1 = in; 125 y1 = out; 126 return out; 127 } 128 129 /// Process one sample using external state variables process_ap(T in,float & x1,float & y1)130 inline T process_ap(T in, float &x1, float &y1) 131 { 132 T out = (in - y1) * a0 + x1; 133 x1 = in; 134 y1 = out; 135 return out; 136 } 137 138 /// Process one sample using external state variables, including filter coeff process_ap(T in,float & x1,float & y1,float a0)139 inline T process_ap(T in, float &x1, float &y1, float a0) 140 { 141 T out = (in - y1) * a0 + x1; 142 x1 = in; 143 y1 = out; 144 return out; 145 } 146 empty()147 inline bool empty() const { 148 return y1 == 0; 149 } 150 sanitize()151 inline void sanitize() 152 { 153 dsp::sanitize(x1); 154 dsp::sanitize(y1); 155 } 156 reset()157 inline void reset() 158 { 159 dsp::zero(x1); 160 dsp::zero(y1); 161 } 162 163 template<class U> copy_coeffs(const onepole<U> & src)164 inline void copy_coeffs(const onepole<U> &src) 165 { 166 a0 = src.a0; 167 a1 = src.a1; 168 b1 = src.b1; 169 } 170 171 /// Return the filter's gain at frequency freq 172 /// @param freq Frequency to look up 173 /// @param sr Filter sample rate (used to convert frequency to angular frequency) freq_gain(float freq,float sr)174 float freq_gain(float freq, float sr) const 175 { 176 freq *= 2.0 * M_PI / sr; 177 cfloat z = 1.0 / exp(cfloat(0.0, freq)); 178 179 return std::abs(h_z(z)); 180 } 181 182 /// Return H(z) the filter's gain at frequency freq 183 /// @param z Z variable (e^jw) h_z(const cfloat & z)184 cfloat h_z(const cfloat &z) const 185 { 186 return (cfloat(a0) + double(a1) * z) / (cfloat(1.0) + double(b1) * z); 187 } 188 }; 189 190 }; 191 192 #endif 193