1 /** @addtogroup biquad
2  *  @{
3  */
4 /*
5   Copyright (C) 2016 D Levin (https://www.kfrlib.com)
6   This file is part of KFR
7 
8   KFR is free software: you can redistribute it and/or modify
9   it under the terms of the GNU General Public License as published by
10   the Free Software Foundation, either version 2 of the License, or
11   (at your option) any later version.
12 
13   KFR is distributed in the hope that it will be useful,
14   but WITHOUT ANY WARRANTY; without even the implied warranty of
15   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16   GNU General Public License for more details.
17 
18   You should have received a copy of the GNU General Public License
19   along with KFR.
20 
21   If GPL is not suitable for your project, you must purchase a commercial license to use KFR.
22   Buying a commercial license is mandatory as soon as you develop commercial activities without
23   disclosing the source code of your own applications.
24   See https://www.kfrlib.com for details.
25  */
26 #pragma once
27 
28 #include "biquad.hpp"
29 #include <cmath>
30 
31 namespace kfr
32 {
33 inline namespace CMT_ARCH_NAME
34 {
35 
36 /**
37  * @brief Calculates coefficients for the all-pass biquad filter
38  * @param frequency Normalized frequency (frequency_Hz / samplerate_Hz)
39  * @param Q Q factor
40  * @return Biquad filter coefficients
41  */
42 template <typename T = fbase>
biquad_allpass(identity<T> frequency,identity<T> Q)43 KFR_FUNCTION biquad_params<T> biquad_allpass(identity<T> frequency, identity<T> Q)
44 {
45     const T alpha = std::sin(frequency) / 2.0 * Q;
46     const T cs    = std::cos(frequency);
47 
48     const T b0 = 1.0 / (1.0 + alpha);
49     const T b1 = -2.0 * cs * b0;
50     const T b2 = (1.0 - alpha) * b0;
51     const T a0 = (1.0 - alpha) * b0;
52     const T a1 = -2.0 * cs * b0;
53     const T a2 = (1.0 + alpha) * b0;
54     return { b0, b1, b2, a0, a1, a2 };
55 }
56 
57 /**
58  * @brief Calculates coefficients for the low-pass biquad filter
59  * @param frequency Normalized frequency (frequency_Hz / samplerate_Hz)
60  * @param Q Q factor
61  * @return Biquad filter coefficients
62  */
63 template <typename T = fbase>
biquad_lowpass(identity<T> frequency,identity<T> Q)64 KFR_FUNCTION biquad_params<T> biquad_lowpass(identity<T> frequency, identity<T> Q)
65 {
66     const T K    = std::tan(c_pi<T, 1> * frequency);
67     const T K2   = K * K;
68     const T norm = 1 / (1 + K / Q + K2);
69     const T a0   = K2 * norm;
70     const T a1   = 2 * a0;
71     const T a2   = a0;
72     const T b1   = 2 * (K2 - 1) * norm;
73     const T b2   = (1 - K / Q + K2) * norm;
74     return { 1.0, b1, b2, a0, a1, a2 };
75 }
76 
77 /**
78  * @brief Calculates coefficients for the high-pass biquad filter
79  * @param frequency Normalized frequency (frequency_Hz / samplerate_Hz)
80  * @param Q Q factor
81  * @return Biquad filter coefficients
82  */
83 template <typename T = fbase>
biquad_highpass(identity<T> frequency,identity<T> Q)84 KFR_FUNCTION biquad_params<T> biquad_highpass(identity<T> frequency, identity<T> Q)
85 {
86     const T K    = std::tan(c_pi<T, 1> * frequency);
87     const T K2   = K * K;
88     const T norm = 1 / (1 + K / Q + K2);
89     const T a0   = 1 * norm;
90     const T a1   = -2 * a0;
91     const T a2   = a0;
92     const T b1   = 2 * (K2 - 1) * norm;
93     const T b2   = (1 - K / Q + K2) * norm;
94     return { 1.0, b1, b2, a0, a1, a2 };
95 }
96 
97 /**
98  * @brief Calculates coefficients for the band-pass biquad filter
99  * @param frequency Normalized frequency (frequency_Hz / samplerate_Hz)
100  * @param Q Q factor
101  * @return Biquad filter coefficients
102  */
103 template <typename T = fbase>
biquad_bandpass(identity<T> frequency,identity<T> Q)104 KFR_FUNCTION biquad_params<T> biquad_bandpass(identity<T> frequency, identity<T> Q)
105 {
106     const T K    = std::tan(c_pi<T, 1> * frequency);
107     const T K2   = K * K;
108     const T norm = 1 / (1 + K / Q + K2);
109     const T a0   = K / Q * norm;
110     const T a1   = 0;
111     const T a2   = -a0;
112     const T b1   = 2 * (K2 - 1) * norm;
113     const T b2   = (1 - K / Q + K2) * norm;
114     return { 1.0, b1, b2, a0, a1, a2 };
115 }
116 
117 /**
118  * @brief Calculates coefficients for the notch biquad filter
119  * @param frequency Normalized frequency (frequency_Hz / samplerate_Hz)
120  * @param Q Q factor
121  * @return Biquad filter coefficients
122  */
123 template <typename T = fbase>
biquad_notch(identity<T> frequency,identity<T> Q)124 KFR_FUNCTION biquad_params<T> biquad_notch(identity<T> frequency, identity<T> Q)
125 {
126     const T K    = std::tan(c_pi<T, 1> * frequency);
127     const T K2   = K * K;
128     const T norm = 1 / (1 + K / Q + K2);
129     const T a0   = (1 + K2) * norm;
130     const T a1   = 2 * (K2 - 1) * norm;
131     const T a2   = a0;
132     const T b1   = a1;
133     const T b2   = (1 - K / Q + K2) * norm;
134     return { 1.0, b1, b2, a0, a1, a2 };
135 }
136 
137 /**
138  * @brief Calculates coefficients for the peak biquad filter
139  * @param frequency Normalized frequency (frequency_Hz / samplerate_Hz)
140  * @param Q Q factor
141  * @param gain Gain in dB
142  * @return Biquad filter coefficients
143  */
144 template <typename T = fbase>
biquad_peak(identity<T> frequency,identity<T> Q,identity<T> gain)145 KFR_FUNCTION biquad_params<T> biquad_peak(identity<T> frequency, identity<T> Q, identity<T> gain)
146 {
147     biquad_params<T> result;
148     const T K  = std::tan(c_pi<T, 1> * frequency);
149     const T K2 = K * K;
150     const T V  = std::exp(std::abs(gain) * (1.0 / 20.0) * c_log_10<T>);
151 
152     if (gain >= 0)
153     { // boost
154         const T norm = 1 / (1 + 1 / Q * K + K2);
155         const T a0   = (1 + V / Q * K + K2) * norm;
156         const T a1   = 2 * (K2 - 1) * norm;
157         const T a2   = (1 - V / Q * K + K2) * norm;
158         const T b1   = a1;
159         const T b2   = (1 - 1 / Q * K + K2) * norm;
160         result       = { 1.0, b1, b2, a0, a1, a2 };
161     }
162     else
163     { // cut
164         const T norm = 1 / (1 + V / Q * K + K2);
165         const T a0   = (1 + 1 / Q * K + K2) * norm;
166         const T a1   = 2 * (K2 - 1) * norm;
167         const T a2   = (1 - 1 / Q * K + K2) * norm;
168         const T b1   = a1;
169         const T b2   = (1 - V / Q * K + K2) * norm;
170         result       = { 1.0, b1, b2, a0, a1, a2 };
171     }
172     return result;
173 }
174 
175 /**
176  * @brief Calculates coefficients for the low-shelf biquad filter
177  * @param frequency Normalized frequency (frequency_Hz / samplerate_Hz)
178  * @param gain Gain in dB
179  * @return Biquad filter coefficients
180  */
181 template <typename T = fbase>
biquad_lowshelf(identity<T> frequency,identity<T> gain)182 KFR_FUNCTION biquad_params<T> biquad_lowshelf(identity<T> frequency, identity<T> gain)
183 {
184     biquad_params<T> result;
185     const T K  = std::tan(c_pi<T, 1> * frequency);
186     const T K2 = K * K;
187     const T V  = std::exp(std::fabs(gain) * (1.0 / 20.0) * c_log_10<T>);
188 
189     if (gain >= 0)
190     { // boost
191         const T norm = 1 / (1 + c_sqrt_2<T> * K + K2);
192         const T a0   = (1 + std::sqrt(2 * V) * K + V * K2) * norm;
193         const T a1   = 2 * (V * K2 - 1) * norm;
194         const T a2   = (1 - std::sqrt(2 * V) * K + V * K2) * norm;
195         const T b1   = 2 * (K2 - 1) * norm;
196         const T b2   = (1 - c_sqrt_2<T> * K + K2) * norm;
197         result       = { 1.0, b1, b2, a0, a1, a2 };
198     }
199     else
200     { // cut
201         const T norm = 1 / (1 + std::sqrt(2 * V) * K + V * K2);
202         const T a0   = (1 + c_sqrt_2<T> * K + K2) * norm;
203         const T a1   = 2 * (K2 - 1) * norm;
204         const T a2   = (1 - c_sqrt_2<T> * K + K2) * norm;
205         const T b1   = 2 * (V * K2 - 1) * norm;
206         const T b2   = (1 - std::sqrt(2 * V) * K + V * K2) * norm;
207         result       = { 1.0, b1, b2, a0, a1, a2 };
208     }
209     return result;
210 }
211 
212 /**
213  * @brief Calculates coefficients for the high-shelf biquad filter
214  * @param frequency Normalized frequency (frequency_Hz / samplerate_Hz)
215  * @param gain Gain in dB
216  * @return Biquad filter coefficients
217  */
218 template <typename T = fbase>
biquad_highshelf(identity<T> frequency,identity<T> gain)219 KFR_FUNCTION biquad_params<T> biquad_highshelf(identity<T> frequency, identity<T> gain)
220 {
221     biquad_params<T> result;
222     const T K  = std::tan(c_pi<T, 1> * frequency);
223     const T K2 = K * K;
224     const T V  = std::exp(std::fabs(gain) * (1.0 / 20.0) * c_log_10<T>);
225 
226     if (gain >= 0)
227     { // boost
228         const T norm = 1 / (1 + c_sqrt_2<T> * K + K2);
229         const T a0   = (V + std::sqrt(2 * V) * K + K2) * norm;
230         const T a1   = 2 * (K2 - V) * norm;
231         const T a2   = (V - std::sqrt(2 * V) * K + K2) * norm;
232         const T b1   = 2 * (K2 - 1) * norm;
233         const T b2   = (1 - c_sqrt_2<T> * K + K2) * norm;
234         result       = { 1.0, b1, b2, a0, a1, a2 };
235     }
236     else
237     { // cut
238         const T norm = 1 / (V + std::sqrt(2 * V) * K + K2);
239         const T a0   = (1 + c_sqrt_2<T> * K + K2) * norm;
240         const T a1   = 2 * (K2 - 1) * norm;
241         const T a2   = (1 - c_sqrt_2<T> * K + K2) * norm;
242         const T b1   = 2 * (K2 - V) * norm;
243         const T b2   = (V - std::sqrt(2 * V) * K + K2) * norm;
244         result       = { 1.0, b1, b2, a0, a1, a2 };
245     }
246     return result;
247 }
248 } // namespace CMT_ARCH_NAME
249 } // namespace kfr
250