1 /* FluidSynth - A Software Synthesizer
2 *
3 * Copyright (C) 2003 Peter Hanappe and others.
4 *
5 * This library is free software; you can redistribute it and/or
6 * modify it under the terms of the GNU Library General Public License
7 * as published by the Free Software Foundation; either version 2 of
8 * the License, or (at your option) any later version.
9 *
10 * This library is distributed in the hope that it will be useful, but
11 * WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 * Library General Public License for more details.
14 *
15 * You should have received a copy of the GNU Library General Public
16 * License along with this library; if not, write to the Free
17 * Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
18 * 02110-1301, USA
19 */
20
21 #include "fluid_iir_filter.h"
22 #include "fluid_sys.h"
23 #include "fluid_conv.h"
24
25 /**
26 * Applies a lowpass filter with variable cutoff frequency and quality factor.
27 * Also modifies filter state accordingly.
28 * @param iir_filter Filter parameter
29 * @param dsp_buf Pointer to the synthesized audio data
30 * @param count Count of samples in dsp_buf
31 */
32 /*
33 * Variable description:
34 * - dsp_a1, dsp_a2, dsp_b0, dsp_b1, dsp_b2: Filter coefficients
35 *
36 * A couple of variables are used internally, their results are discarded:
37 * - dsp_i: Index through the output buffer
38 * - dsp_phase_fractional: The fractional part of dsp_phase
39 * - dsp_coeff: A table of four coefficients, depending on the fractional phase.
40 * Used to interpolate between samples.
41 * - dsp_process_buffer: Holds the processed signal between stages
42 * - dsp_centernode: delay line for the IIR filter
43 * - dsp_hist1: same
44 * - dsp_hist2: same
45 */
46 void
fluid_iir_filter_apply(fluid_iir_filter_t * iir_filter,fluid_real_t * dsp_buf,int count)47 fluid_iir_filter_apply(fluid_iir_filter_t* iir_filter,
48 fluid_real_t *dsp_buf, int count)
49 {
50 /* IIR filter sample history */
51 fluid_real_t dsp_hist1 = iir_filter->hist1;
52 fluid_real_t dsp_hist2 = iir_filter->hist2;
53
54 /* IIR filter coefficients */
55 fluid_real_t dsp_a1 = iir_filter->a1;
56 fluid_real_t dsp_a2 = iir_filter->a2;
57 fluid_real_t dsp_b02 = iir_filter->b02;
58 fluid_real_t dsp_b1 = iir_filter->b1;
59 int dsp_filter_coeff_incr_count = iir_filter->filter_coeff_incr_count;
60
61 fluid_real_t dsp_centernode;
62 int dsp_i;
63
64 /* filter (implement the voice filter according to SoundFont standard) */
65
66 /* Check for denormal number (too close to zero). */
67 if (fabs (dsp_hist1) < 1e-20) dsp_hist1 = 0.0f; /* FIXME JMG - Is this even needed? */
68
69 /* Two versions of the filter loop. One, while the filter is
70 * changing towards its new setting. The other, if the filter
71 * doesn't change.
72 */
73
74 if (dsp_filter_coeff_incr_count > 0)
75 {
76 fluid_real_t dsp_a1_incr = iir_filter->a1_incr;
77 fluid_real_t dsp_a2_incr = iir_filter->a2_incr;
78 fluid_real_t dsp_b02_incr = iir_filter->b02_incr;
79 fluid_real_t dsp_b1_incr = iir_filter->b1_incr;
80
81
82 /* Increment is added to each filter coefficient filter_coeff_incr_count times. */
83 for (dsp_i = 0; dsp_i < count; dsp_i++)
84 {
85 /* The filter is implemented in Direct-II form. */
86 dsp_centernode = dsp_buf[dsp_i] - dsp_a1 * dsp_hist1 - dsp_a2 * dsp_hist2;
87 dsp_buf[dsp_i] = dsp_b02 * (dsp_centernode + dsp_hist2) + dsp_b1 * dsp_hist1;
88 dsp_hist2 = dsp_hist1;
89 dsp_hist1 = dsp_centernode;
90
91 if (dsp_filter_coeff_incr_count-- > 0)
92 {
93 fluid_real_t old_b02 = dsp_b02;
94 dsp_a1 += dsp_a1_incr;
95 dsp_a2 += dsp_a2_incr;
96 dsp_b02 += dsp_b02_incr;
97 dsp_b1 += dsp_b1_incr;
98
99 /* Compensate history to avoid the filter going havoc with large frequency changes */
100 if (iir_filter->compensate_incr && fabs(dsp_b02) > 0.001) {
101 fluid_real_t compensate = old_b02 / dsp_b02;
102 dsp_centernode *= compensate;
103 dsp_hist1 *= compensate;
104 dsp_hist2 *= compensate;
105 }
106 }
107 } /* for dsp_i */
108 }
109 else /* The filter parameters are constant. This is duplicated to save time. */
110 {
111 for (dsp_i = 0; dsp_i < count; dsp_i++)
112 { /* The filter is implemented in Direct-II form. */
113 dsp_centernode = dsp_buf[dsp_i] - dsp_a1 * dsp_hist1 - dsp_a2 * dsp_hist2;
114 dsp_buf[dsp_i] = dsp_b02 * (dsp_centernode + dsp_hist2) + dsp_b1 * dsp_hist1;
115 dsp_hist2 = dsp_hist1;
116 dsp_hist1 = dsp_centernode;
117 }
118 }
119
120 iir_filter->hist1 = dsp_hist1;
121 iir_filter->hist2 = dsp_hist2;
122 iir_filter->a1 = dsp_a1;
123 iir_filter->a2 = dsp_a2;
124 iir_filter->b02 = dsp_b02;
125 iir_filter->b1 = dsp_b1;
126 iir_filter->filter_coeff_incr_count = dsp_filter_coeff_incr_count;
127
128 fluid_check_fpe ("voice_filter");
129 }
130
131
132 void
fluid_iir_filter_reset(fluid_iir_filter_t * iir_filter)133 fluid_iir_filter_reset(fluid_iir_filter_t* iir_filter)
134 {
135 iir_filter->hist1 = 0;
136 iir_filter->hist2 = 0;
137 iir_filter->last_fres = -1.;
138 iir_filter->filter_startup = 1;
139 }
140
141 void
fluid_iir_filter_set_fres(fluid_iir_filter_t * iir_filter,fluid_real_t fres)142 fluid_iir_filter_set_fres(fluid_iir_filter_t* iir_filter,
143 fluid_real_t fres)
144 {
145 iir_filter->fres = fres;
146 iir_filter->last_fres = -1.;
147 }
148
149
150 void
fluid_iir_filter_set_q_dB(fluid_iir_filter_t * iir_filter,fluid_real_t q_dB)151 fluid_iir_filter_set_q_dB(fluid_iir_filter_t* iir_filter,
152 fluid_real_t q_dB)
153 {
154 /* The 'sound font' Q is defined in dB. The filter needs a linear
155 q. Convert. */
156 iir_filter->q_lin = (fluid_real_t) (pow(10.0f, q_dB / 20.0f));
157
158 /* SF 2.01 page 59:
159 *
160 * The SoundFont specs ask for a gain reduction equal to half the
161 * height of the resonance peak (Q). For example, for a 10 dB
162 * resonance peak, the gain is reduced by 5 dB. This is done by
163 * multiplying the total gain with sqrt(1/Q). `Sqrt' divides dB
164 * by 2 (100 lin = 40 dB, 10 lin = 20 dB, 3.16 lin = 10 dB etc)
165 * The gain is later factored into the 'b' coefficients
166 * (numerator of the filter equation). This gain factor depends
167 * only on Q, so this is the right place to calculate it.
168 */
169 iir_filter->filter_gain = (fluid_real_t) (1.0 / sqrt(iir_filter->q_lin));
170
171 /* The synthesis loop will have to recalculate the filter coefficients. */
172 iir_filter->last_fres = -1.;
173
174 }
175
176
177 static inline void
fluid_iir_filter_calculate_coefficients(fluid_iir_filter_t * iir_filter,int transition_samples,fluid_real_t output_rate)178 fluid_iir_filter_calculate_coefficients(fluid_iir_filter_t* iir_filter,
179 int transition_samples,
180 fluid_real_t output_rate)
181 {
182
183 /*
184 * Those equations from Robert Bristow-Johnson's `Cookbook
185 * formulae for audio EQ biquad filter coefficients', obtained
186 * from Harmony-central.com / Computer / Programming. They are
187 * the result of the bilinear transform on an analogue filter
188 * prototype. To quote, `BLT frequency warping has been taken
189 * into account for both significant frequency relocation and for
190 * bandwidth readjustment'. */
191
192 fluid_real_t omega = (fluid_real_t) (2.0 * M_PI *
193 (iir_filter->last_fres / ((float) output_rate)));
194 fluid_real_t sin_coeff = (fluid_real_t) sin(omega);
195 fluid_real_t cos_coeff = (fluid_real_t) cos(omega);
196 fluid_real_t alpha_coeff = sin_coeff / (2.0f * iir_filter->q_lin);
197 fluid_real_t a0_inv = 1.0f / (1.0f + alpha_coeff);
198
199 /* Calculate the filter coefficients. All coefficients are
200 * normalized by a0. Think of `a1' as `a1/a0'.
201 *
202 * Here a couple of multiplications are saved by reusing common expressions.
203 * The original equations should be:
204 * iir_filter->b0=(1.-cos_coeff)*a0_inv*0.5*iir_filter->filter_gain;
205 * iir_filter->b1=(1.-cos_coeff)*a0_inv*iir_filter->filter_gain;
206 * iir_filter->b2=(1.-cos_coeff)*a0_inv*0.5*iir_filter->filter_gain; */
207
208 fluid_real_t a1_temp = -2.0f * cos_coeff * a0_inv;
209 fluid_real_t a2_temp = (1.0f - alpha_coeff) * a0_inv;
210 fluid_real_t b1_temp = (1.0f - cos_coeff) * a0_inv * iir_filter->filter_gain;
211 /* both b0 -and- b2 */
212 fluid_real_t b02_temp = b1_temp * 0.5f;
213
214 iir_filter->compensate_incr = 0;
215
216 if (iir_filter->filter_startup || (transition_samples == 0))
217 {
218 /* The filter is calculated, because the voice was started up.
219 * In this case set the filter coefficients without delay.
220 */
221 iir_filter->a1 = a1_temp;
222 iir_filter->a2 = a2_temp;
223 iir_filter->b02 = b02_temp;
224 iir_filter->b1 = b1_temp;
225 iir_filter->filter_coeff_incr_count = 0;
226 iir_filter->filter_startup = 0;
227 // printf("Setting initial filter coefficients.\n");
228 }
229 else
230 {
231
232 /* The filter frequency is changed. Calculate an increment
233 * factor, so that the new setting is reached after one buffer
234 * length. x_incr is added to the current value FLUID_BUFSIZE
235 * times. The length is arbitrarily chosen. Longer than one
236 * buffer will sacrifice some performance, though. Note: If
237 * the filter is still too 'grainy', then increase this number
238 * at will.
239 */
240
241 iir_filter->a1_incr = (a1_temp - iir_filter->a1) / transition_samples;
242 iir_filter->a2_incr = (a2_temp - iir_filter->a2) / transition_samples;
243 iir_filter->b02_incr = (b02_temp - iir_filter->b02) / transition_samples;
244 iir_filter->b1_incr = (b1_temp - iir_filter->b1) / transition_samples;
245 if (fabs(iir_filter->b02) > 0.0001) {
246 fluid_real_t quota = b02_temp / iir_filter->b02;
247 iir_filter->compensate_incr = quota < 0.5 || quota > 2;
248 }
249 /* Have to add the increments filter_coeff_incr_count times. */
250 iir_filter->filter_coeff_incr_count = transition_samples;
251 }
252 fluid_check_fpe ("voice_write filter calculation");
253 }
254
255
fluid_iir_filter_calc(fluid_iir_filter_t * iir_filter,fluid_real_t output_rate,fluid_real_t fres_mod)256 void fluid_iir_filter_calc(fluid_iir_filter_t* iir_filter,
257 fluid_real_t output_rate,
258 fluid_real_t fres_mod)
259 {
260 fluid_real_t fres;
261
262 /* calculate the frequency of the resonant filter in Hz */
263 fres = fluid_ct2hz(iir_filter->fres + fres_mod);
264
265 /* FIXME - Still potential for a click during turn on, can we interpolate
266 between 20khz cutoff and 0 Q? */
267
268 /* I removed the optimization of turning the filter off when the
269 * resonance frequence is above the maximum frequency. Instead, the
270 * filter frequency is set to a maximum of 0.45 times the sampling
271 * rate. For a 44100 kHz sampling rate, this amounts to 19845
272 * Hz. The reason is that there were problems with anti-aliasing when the
273 * synthesizer was run at lower sampling rates. Thanks to Stephan
274 * Tassart for pointing me to this bug. By turning the filter on and
275 * clipping the maximum filter frequency at 0.45*srate, the filter
276 * is used as an anti-aliasing filter. */
277
278 if (fres > 0.45f * output_rate)
279 fres = 0.45f * output_rate;
280 else if (fres < 5)
281 fres = 5;
282
283 /* if filter enabled and there is a significant frequency change.. */
284 if ((abs (fres - iir_filter->last_fres) > 0.01))
285 {
286 /* The filter coefficients have to be recalculated (filter
287 * parameters have changed). Recalculation for various reasons is
288 * forced by setting last_fres to -1. The flag filter_startup
289 * indicates, that the DSP loop runs for the first time, in this
290 * case, the filter is set directly, instead of smoothly fading
291 * between old and new settings. */
292 iir_filter->last_fres = fres;
293 fluid_iir_filter_calculate_coefficients(iir_filter, FLUID_BUFSIZE,
294 output_rate);
295 }
296
297
298 fluid_check_fpe ("voice_write DSP coefficients");
299
300 }
301
302