Line data Source code
1 : /**
2 : * \file ButterworthFilter.hxx
3 : */
4 :
5 : #include "ButterworthFilter.h"
6 : #include <ATK/EQ/helpers.h>
7 :
8 : #include <boost/math/constants/constants.hpp>
9 : #include <boost/math/tools/polynomial.hpp>
10 :
11 : namespace ButterworthUtilities
12 : {
13 : template<typename DataType>
14 124 : void create_butterworth_analog_coefficients(int order, EQUtilities::ZPK<DataType>& zpk)
15 : {
16 124 : zpk.k = 1;
17 124 : zpk.z.clear(); // no zeros for this filter type
18 124 : zpk.p.clear();
19 312 : for(gsl::index i = -order+1; i < order; i += 2)
20 : {
21 188 : zpk.p.push_back(std::complex<DataType>(-std::cos(boost::math::constants::pi<DataType>() * i / (2 * order)), -std::sin(boost::math::constants::pi<DataType>() * i / (2 * order))));
22 : }
23 124 : }
24 :
25 : template<typename DataType, typename Container>
26 88 : void create_default_coeffs(size_t order, DataType Wn, Container& coefficients_in, Container& coefficients_out)
27 : {
28 176 : EQUtilities::ZPK<DataType> zpk;
29 :
30 88 : int fs = 2;
31 88 : create_butterworth_analog_coefficients(static_cast<int>(order), zpk);
32 88 : EQUtilities::populate_lp_coeffs(Wn, fs, order, zpk, coefficients_in, coefficients_out);
33 88 : }
34 :
35 : template<typename DataType, typename Container>
36 18 : void create_bp_coeffs(size_t order, DataType wc1, DataType wc2, Container& coefficients_in, Container& coefficients_out)
37 : {
38 36 : EQUtilities::ZPK<DataType> zpk;
39 :
40 18 : int fs = 2;
41 18 : create_butterworth_analog_coefficients(static_cast<int>(order/2), zpk);
42 18 : EQUtilities::populate_bp_coeffs(wc1, wc2, fs, order, zpk, coefficients_in, coefficients_out);
43 18 : }
44 :
45 : template<typename DataType, typename Container>
46 18 : void create_bs_coeffs(size_t order, DataType wc1, DataType wc2, Container& coefficients_in, Container& coefficients_out)
47 : {
48 36 : EQUtilities::ZPK<DataType> zpk;
49 :
50 18 : int fs = 2;
51 18 : create_butterworth_analog_coefficients(static_cast<int>(order/2), zpk);
52 18 : EQUtilities::populate_bs_coeffs(wc1, wc2, fs, order, zpk, coefficients_in, coefficients_out);
53 18 : }
54 : }
55 :
56 : namespace ATK
57 : {
58 : template <typename DataType>
59 20 : ButterworthLowPassCoefficients<DataType>::ButterworthLowPassCoefficients(gsl::index nb_channels)
60 20 : :Parent(nb_channels, nb_channels)
61 : {
62 20 : }
63 :
64 : template <typename DataType_>
65 18 : void ButterworthLowPassCoefficients<DataType_>::set_cut_frequency(CoeffDataType cut_frequency)
66 : {
67 18 : if(cut_frequency <= 0)
68 : {
69 1 : throw std::out_of_range("Frequency can't be negative");
70 : }
71 17 : this->cut_frequency = cut_frequency;
72 17 : setup();
73 17 : }
74 :
75 : template <typename DataType_>
76 1 : typename ButterworthLowPassCoefficients<DataType_>::CoeffDataType ButterworthLowPassCoefficients<DataType_>::get_cut_frequency() const
77 : {
78 1 : return cut_frequency;
79 : }
80 :
81 : template <typename DataType>
82 18 : void ButterworthLowPassCoefficients<DataType>::set_order(unsigned int order)
83 : {
84 18 : if(order == 0)
85 : {
86 1 : throw std::out_of_range("Order can't be null");
87 : }
88 17 : in_order = out_order = order;
89 17 : setup();
90 17 : }
91 :
92 : template <typename DataType>
93 1 : unsigned int ButterworthLowPassCoefficients<DataType>::get_order() const
94 : {
95 1 : return in_order;
96 : }
97 :
98 : template <typename DataType>
99 66 : void ButterworthLowPassCoefficients<DataType>::setup()
100 : {
101 66 : Parent::setup();
102 66 : coefficients_in.assign(in_order+1, 0);
103 66 : coefficients_out.assign(out_order, 0);
104 :
105 66 : ButterworthUtilities::create_default_coeffs(in_order, 2 * cut_frequency / input_sampling_rate, coefficients_in, coefficients_out);
106 66 : }
107 :
108 : template <typename DataType>
109 9 : ButterworthHighPassCoefficients<DataType>::ButterworthHighPassCoefficients(gsl::index nb_channels)
110 9 : :Parent(nb_channels, nb_channels)
111 : {
112 9 : }
113 :
114 : template <typename DataType_>
115 7 : void ButterworthHighPassCoefficients<DataType_>::set_cut_frequency(CoeffDataType cut_frequency)
116 : {
117 7 : if(cut_frequency <= 0)
118 : {
119 1 : throw std::out_of_range("Frequency can't be negative");
120 : }
121 6 : this->cut_frequency = cut_frequency;
122 6 : setup();
123 6 : }
124 :
125 : template <typename DataType_>
126 1 : typename ButterworthHighPassCoefficients<DataType_>::CoeffDataType ButterworthHighPassCoefficients<DataType_>::get_cut_frequency() const
127 : {
128 1 : return cut_frequency;
129 : }
130 :
131 : template <typename DataType>
132 7 : void ButterworthHighPassCoefficients<DataType>::set_order(unsigned int order)
133 : {
134 7 : if(order == 0)
135 : {
136 1 : throw std::out_of_range("Order can't be null");
137 : }
138 6 : in_order = out_order = order;
139 6 : setup();
140 6 : }
141 :
142 : template <typename DataType>
143 1 : unsigned int ButterworthHighPassCoefficients<DataType>::get_order() const
144 : {
145 1 : return in_order;
146 : }
147 :
148 : template <typename DataType>
149 22 : void ButterworthHighPassCoefficients<DataType>::setup()
150 : {
151 22 : Parent::setup();
152 22 : coefficients_in.assign(in_order+1, 0);
153 22 : coefficients_out.assign(out_order, 0);
154 :
155 22 : ButterworthUtilities::create_default_coeffs(in_order, (input_sampling_rate - 2 * cut_frequency) / input_sampling_rate, coefficients_in, coefficients_out);
156 50 : for(gsl::index i = in_order - 1; i >= 0; i -= 2)
157 : {
158 28 : coefficients_in[i] = - coefficients_in[i];
159 28 : coefficients_out[i] = - coefficients_out[i];
160 : }
161 22 : }
162 :
163 : template <typename DataType>
164 9 : ButterworthBandPassCoefficients<DataType>::ButterworthBandPassCoefficients(gsl::index nb_channels)
165 9 : :Parent(nb_channels, nb_channels)
166 : {
167 9 : }
168 :
169 : template <typename DataType_>
170 7 : void ButterworthBandPassCoefficients<DataType_>::set_cut_frequencies(std::pair<CoeffDataType, CoeffDataType> cut_frequencies)
171 : {
172 7 : if(cut_frequencies.first <= 0 || cut_frequencies.second <= 0)
173 : {
174 2 : throw std::out_of_range("Frequencies can't be negative");
175 : }
176 5 : this->cut_frequencies = cut_frequencies;
177 5 : setup();
178 5 : }
179 :
180 : template <typename DataType_>
181 3 : void ButterworthBandPassCoefficients<DataType_>::set_cut_frequencies(CoeffDataType f0, CoeffDataType f1)
182 : {
183 3 : set_cut_frequencies(std::make_pair(f0, f1));
184 1 : }
185 :
186 : template <typename DataType_>
187 2 : std::pair<typename ButterworthBandPassCoefficients<DataType_>::CoeffDataType, typename ButterworthBandPassCoefficients<DataType_>::CoeffDataType> ButterworthBandPassCoefficients<DataType_>::get_cut_frequencies() const
188 : {
189 2 : return cut_frequencies;
190 : }
191 :
192 : template <typename DataType>
193 6 : void ButterworthBandPassCoefficients<DataType>::set_order(unsigned int order)
194 : {
195 6 : if(order == 0)
196 : {
197 1 : throw std::out_of_range("Order can't be null");
198 : }
199 5 : in_order = out_order = 2 * order;
200 5 : setup();
201 5 : }
202 :
203 : template <typename DataType>
204 1 : unsigned int ButterworthBandPassCoefficients<DataType>::get_order() const
205 : {
206 1 : return in_order / 2;
207 : }
208 :
209 : template <typename DataType>
210 18 : void ButterworthBandPassCoefficients<DataType>::setup()
211 : {
212 18 : Parent::setup();
213 18 : coefficients_in.assign(in_order+1, 0);
214 18 : coefficients_out.assign(out_order, 0);
215 :
216 18 : ButterworthUtilities::create_bp_coeffs(in_order, 2 * cut_frequencies.first / input_sampling_rate, 2 * cut_frequencies.second / input_sampling_rate, coefficients_in, coefficients_out);
217 18 : }
218 :
219 : template <typename DataType>
220 9 : ButterworthBandStopCoefficients<DataType>::ButterworthBandStopCoefficients(gsl::index nb_channels)
221 9 : :Parent(nb_channels, nb_channels)
222 : {
223 9 : }
224 :
225 : template <typename DataType_>
226 7 : void ButterworthBandStopCoefficients<DataType_>::set_cut_frequencies(std::pair<CoeffDataType, CoeffDataType> cut_frequencies)
227 : {
228 7 : if(cut_frequencies.first <= 0 || cut_frequencies.second <= 0)
229 : {
230 2 : throw std::out_of_range("Frequencies can't be negative");
231 : }
232 5 : this->cut_frequencies = cut_frequencies;
233 5 : setup();
234 5 : }
235 :
236 : template <typename DataType_>
237 3 : void ButterworthBandStopCoefficients<DataType_>::set_cut_frequencies(CoeffDataType f0, CoeffDataType f1)
238 : {
239 3 : set_cut_frequencies(std::make_pair(f0, f1));
240 1 : }
241 :
242 : template <typename DataType_>
243 2 : std::pair<typename ButterworthBandStopCoefficients<DataType_>::CoeffDataType, typename ButterworthBandStopCoefficients<DataType_>::CoeffDataType> ButterworthBandStopCoefficients<DataType_>::get_cut_frequencies() const
244 : {
245 2 : return cut_frequencies;
246 : }
247 :
248 : template <typename DataType>
249 6 : void ButterworthBandStopCoefficients<DataType>::set_order(unsigned int order)
250 : {
251 6 : if(order == 0)
252 : {
253 1 : throw std::out_of_range("Order can't be null");
254 : }
255 5 : in_order = out_order = 2 * order;
256 5 : setup();
257 5 : }
258 :
259 : template <typename DataType>
260 1 : unsigned int ButterworthBandStopCoefficients<DataType>::get_order() const
261 : {
262 1 : return in_order / 2;
263 : }
264 :
265 : template <typename DataType>
266 18 : void ButterworthBandStopCoefficients<DataType>::setup()
267 : {
268 18 : Parent::setup();
269 18 : coefficients_in.assign(in_order+1, 0);
270 18 : coefficients_out.assign(out_order, 0);
271 :
272 18 : ButterworthUtilities::create_bs_coeffs(in_order, 2 * cut_frequencies.first / input_sampling_rate, 2 * cut_frequencies.second / input_sampling_rate, coefficients_in, coefficients_out);
273 18 : }
274 : }
|