Little speed improvement
[calf.git] / src / audio_fx.cpp
blob36ada2819766a037de8bd3f8e187f12533653d93
1 /* Calf DSP Library
2 * Reusable audio effect classes - implementation.
4 * Copyright (C) 2001-2010 Krzysztof Foltman, Markus Schmidt, Thor Harald Johansen and others
6 * This program is free software; you can redistribute it and/or
7 * modify it under the terms of the GNU Lesser General Public
8 * License as published by the Free Software Foundation; either
9 * version 2 of the License, or (at your option) any later version.
11 * This program is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 * Lesser General Public License for more details.
16 * You should have received a copy of the GNU Lesser General
17 * Public License along with this program; if not, write to the
18 * Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
19 * Boston, MA 02110-1301 USA
22 #include <calf/audio_fx.h>
23 #include <calf/giface.h>
24 #include <limits.h>
25 #include <stdlib.h>
26 #include <time.h>
27 #include <math.h>
29 using namespace calf_plugins;
30 using namespace dsp;
32 simple_phaser::simple_phaser(int _max_stages, float *x1vals, float *y1vals)
34 max_stages = _max_stages;
35 x1 = x1vals;
36 y1 = y1vals;
38 set_base_frq(1000);
39 set_mod_depth(1000);
40 set_fb(0);
41 state = 0;
42 cnt = 0;
43 stages = 0;
44 set_stages(_max_stages);
47 void simple_phaser::set_stages(int _stages)
49 if (_stages > stages)
51 assert(_stages <= max_stages);
52 if (_stages > max_stages)
53 _stages = max_stages;
54 for (int i = stages; i < _stages; i++)
56 x1[i] = x1[stages-1];
57 y1[i] = y1[stages-1];
60 stages = _stages;
63 void simple_phaser::reset()
65 cnt = 0;
66 state = 0;
67 phase.set(0);
68 for (int i = 0; i < max_stages; i++)
69 x1[i] = y1[i] = 0;
70 control_step();
73 void simple_phaser::control_step()
75 cnt = 0;
76 int v = phase.get() + 0x40000000;
77 int sign = v >> 31;
78 v ^= sign;
79 // triangle wave, range from 0 to INT_MAX
80 double vf = (double)((v >> 16) * (1.0 / 16384.0) - 1);
82 float freq = base_frq * pow(2.0, vf * mod_depth / 1200.0);
83 freq = dsp::clip<float>(freq, 10.0, 0.49 * sample_rate);
84 stage1.set_ap_w(freq * (M_PI / 2.0) * odsr);
85 phase += dphase * 32;
86 for (int i = 0; i < stages; i++)
88 dsp::sanitize(x1[i]);
89 dsp::sanitize(y1[i]);
91 dsp::sanitize(state);
94 void simple_phaser::process(float *buf_out, float *buf_in, int nsamples)
96 for (int i=0; i<nsamples; i++) {
97 cnt++;
98 if (cnt == 32)
99 control_step();
100 float in = *buf_in++;
101 float fd = in + state * fb;
102 for (int j = 0; j < stages; j++)
103 fd = stage1.process_ap(fd, x1[j], y1[j]);
104 state = fd;
106 float sdry = in * gs_dry.get();
107 float swet = fd * gs_wet.get();
108 *buf_out++ = sdry + swet;
112 float simple_phaser::freq_gain(float freq, float sr) const
114 typedef std::complex<double> cfloat;
115 freq *= 2.0 * M_PI / sr;
116 cfloat z = 1.0 / exp(cfloat(0.0, freq)); // z^-1
118 cfloat p = cfloat(1.0);
119 cfloat stg = stage1.h_z(z);
121 for (int i = 0; i < stages; i++)
122 p = p * stg;
124 p = p / (cfloat(1.0) - cfloat(fb) * p);
125 return std::abs(cfloat(gs_dry.get_last()) + cfloat(gs_wet.get_last()) * p);
128 ///////////////////////////////////////////////////////////////////////////////////
130 void biquad_filter_module::calculate_filter(float freq, float q, int mode, float gain)
132 if (mode <= mode_36db_lp) {
133 order = mode + 1;
134 left[0].set_lp_rbj(freq, pow(q, 1.0 / order), srate, gain);
135 } else if ( mode_12db_hp <= mode && mode <= mode_36db_hp ) {
136 order = mode - mode_12db_hp + 1;
137 left[0].set_hp_rbj(freq, pow(q, 1.0 / order), srate, gain);
138 } else if ( mode_6db_bp <= mode && mode <= mode_18db_bp ) {
139 order = mode - mode_6db_bp + 1;
140 left[0].set_bp_rbj(freq, pow(q, 1.0 / order), srate, gain);
141 } else { // mode_6db_br <= mode <= mode_18db_br
142 order = mode - mode_6db_br + 1;
143 left[0].set_br_rbj(freq, order * 0.1 * q, srate, gain);
146 right[0].copy_coeffs(left[0]);
147 for (int i = 1; i < order; i++) {
148 left[i].copy_coeffs(left[0]);
149 right[i].copy_coeffs(left[0]);
153 void biquad_filter_module::filter_activate()
155 for (int i=0; i < order; i++) {
156 left[i].reset();
157 right[i].reset();
161 void biquad_filter_module::sanitize()
163 for (int i=0; i < order; i++) {
164 left[i].sanitize();
165 right[i].sanitize();
169 int biquad_filter_module::process_channel(uint16_t channel_no, const float *in, float *out, uint32_t numsamples, int inmask) {
170 dsp::biquad_d1 *filter;
171 switch (channel_no) {
172 case 0:
173 filter = left;
174 break;
176 case 1:
177 filter = right;
178 break;
180 default:
181 assert(false);
182 return 0;
185 if (inmask) {
186 switch(order) {
187 case 1:
188 for (uint32_t i = 0; i < numsamples; i++)
189 out[i] = filter[0].process(in[i]);
190 break;
191 case 2:
192 for (uint32_t i = 0; i < numsamples; i++)
193 out[i] = filter[1].process(filter[0].process(in[i]));
194 break;
195 case 3:
196 for (uint32_t i = 0; i < numsamples; i++)
197 out[i] = filter[2].process(filter[1].process(filter[0].process(in[i])));
198 break;
200 } else {
201 if (filter[order - 1].empty())
202 return 0;
203 switch(order) {
204 case 1:
205 for (uint32_t i = 0; i < numsamples; i++)
206 out[i] = filter[0].process_zeroin();
207 break;
208 case 2:
209 if (filter[0].empty())
210 for (uint32_t i = 0; i < numsamples; i++)
211 out[i] = filter[1].process_zeroin();
212 else
213 for (uint32_t i = 0; i < numsamples; i++)
214 out[i] = filter[1].process(filter[0].process_zeroin());
215 break;
216 case 3:
217 if (filter[1].empty())
218 for (uint32_t i = 0; i < numsamples; i++)
219 out[i] = filter[2].process_zeroin();
220 else
221 for (uint32_t i = 0; i < numsamples; i++)
222 out[i] = filter[2].process(filter[1].process(filter[0].process_zeroin()));
223 break;
226 for (int i = 0; i < order; i++)
227 filter[i].sanitize();
228 return filter[order - 1].empty() ? 0 : inmask;
231 float biquad_filter_module::freq_gain(int subindex, float freq, float srate) const
233 float level = 1.0;
234 for (int j = 0; j < order; j++)
235 level *= left[j].freq_gain(freq, srate);
236 return level;
239 /////////////////////////////////////////////////////////////////////////////////////////////////////////
241 void reverb::update_times()
243 switch(type)
245 case 0:
246 tl[0] = 397 << 16, tr[0] = 383 << 16;
247 tl[1] = 457 << 16, tr[1] = 429 << 16;
248 tl[2] = 549 << 16, tr[2] = 631 << 16;
249 tl[3] = 649 << 16, tr[3] = 756 << 16;
250 tl[4] = 773 << 16, tr[4] = 803 << 16;
251 tl[5] = 877 << 16, tr[5] = 901 << 16;
252 break;
253 case 1:
254 tl[0] = 697 << 16, tr[0] = 783 << 16;
255 tl[1] = 957 << 16, tr[1] = 929 << 16;
256 tl[2] = 649 << 16, tr[2] = 531 << 16;
257 tl[3] = 1049 << 16, tr[3] = 1177 << 16;
258 tl[4] = 473 << 16, tr[4] = 501 << 16;
259 tl[5] = 587 << 16, tr[5] = 681 << 16;
260 break;
261 case 2:
262 default:
263 tl[0] = 697 << 16, tr[0] = 783 << 16;
264 tl[1] = 957 << 16, tr[1] = 929 << 16;
265 tl[2] = 649 << 16, tr[2] = 531 << 16;
266 tl[3] = 1249 << 16, tr[3] = 1377 << 16;
267 tl[4] = 1573 << 16, tr[4] = 1671 << 16;
268 tl[5] = 1877 << 16, tr[5] = 1781 << 16;
269 break;
270 case 3:
271 tl[0] = 1097 << 16, tr[0] = 1087 << 16;
272 tl[1] = 1057 << 16, tr[1] = 1031 << 16;
273 tl[2] = 1049 << 16, tr[2] = 1039 << 16;
274 tl[3] = 1083 << 16, tr[3] = 1055 << 16;
275 tl[4] = 1075 << 16, tr[4] = 1099 << 16;
276 tl[5] = 1003 << 16, tr[5] = 1073 << 16;
277 break;
278 case 4:
279 tl[0] = 197 << 16, tr[0] = 133 << 16;
280 tl[1] = 357 << 16, tr[1] = 229 << 16;
281 tl[2] = 549 << 16, tr[2] = 431 << 16;
282 tl[3] = 949 << 16, tr[3] = 1277 << 16;
283 tl[4] = 1173 << 16, tr[4] = 1671 << 16;
284 tl[5] = 1477 << 16, tr[5] = 1881 << 16;
285 break;
286 case 5:
287 tl[0] = 197 << 16, tr[0] = 133 << 16;
288 tl[1] = 257 << 16, tr[1] = 179 << 16;
289 tl[2] = 549 << 16, tr[2] = 431 << 16;
290 tl[3] = 619 << 16, tr[3] = 497 << 16;
291 tl[4] = 1173 << 16, tr[4] = 1371 << 16;
292 tl[5] = 1577 << 16, tr[5] = 1881 << 16;
293 break;
296 float fDec=1000 + 2400.f * diffusion;
297 for (int i = 0 ; i < 6; i++) {
298 ldec[i]=exp(-float(tl[i] >> 16) / fDec),
299 rdec[i]=exp(-float(tr[i] >> 16) / fDec);
303 void reverb::reset()
305 apL1.reset();apR1.reset();
306 apL2.reset();apR2.reset();
307 apL3.reset();apR3.reset();
308 apL4.reset();apR4.reset();
309 apL5.reset();apR5.reset();
310 apL6.reset();apR6.reset();
311 lp_left.reset();lp_right.reset();
312 old_left = 0; old_right = 0;
315 void reverb::process(float &left, float &right)
317 unsigned int ipart = phase.ipart();
319 // the interpolated LFO might be an overkill here
320 int lfo = phase.lerp_by_fract_int<int, 14, int>(sine.data[ipart], sine.data[ipart+1]) >> 2;
321 phase += dphase;
323 left += old_right;
324 left = apL1.process_allpass_comb_lerp16(left, tl[0] - 45*lfo, ldec[0]);
325 left = apL2.process_allpass_comb_lerp16(left, tl[1] + 47*lfo, ldec[1]);
326 float out_left = left;
327 left = apL3.process_allpass_comb_lerp16(left, tl[2] + 54*lfo, ldec[2]);
328 left = apL4.process_allpass_comb_lerp16(left, tl[3] - 69*lfo, ldec[3]);
329 left = apL5.process_allpass_comb_lerp16(left, tl[4] + 69*lfo, ldec[4]);
330 left = apL6.process_allpass_comb_lerp16(left, tl[5] - 46*lfo, ldec[5]);
331 old_left = lp_left.process(left * fb);
332 sanitize(old_left);
334 right += old_left;
335 right = apR1.process_allpass_comb_lerp16(right, tr[0] - 45*lfo, rdec[0]);
336 right = apR2.process_allpass_comb_lerp16(right, tr[1] + 47*lfo, rdec[1]);
337 float out_right = right;
338 right = apR3.process_allpass_comb_lerp16(right, tr[2] + 54*lfo, rdec[2]);
339 right = apR4.process_allpass_comb_lerp16(right, tr[3] - 69*lfo, rdec[3]);
340 right = apR5.process_allpass_comb_lerp16(right, tr[4] + 69*lfo, rdec[4]);
341 right = apR6.process_allpass_comb_lerp16(right, tr[5] - 46*lfo, rdec[5]);
342 old_right = lp_right.process(right * fb);
343 sanitize(old_right);
345 left = out_left, right = out_right;
348 /// Distortion Module by Tom Szilagyi
350 /// This module provides a blendable saturation stage
351 ///////////////////////////////////////////////////////////////////////////////////////////////
353 tap_distortion::tap_distortion()
355 is_active = false;
356 srate = 0;
357 meter = 0.f;
358 rdrive = rbdr = kpa = kpb = kna = knb = ap = an = imr = kc = srct = sq = pwrq = prev_med = prev_out = 0.f;
359 drive_old = blend_old = -1.f;
360 over = 1;
363 void tap_distortion::activate()
365 is_active = true;
366 set_params(0.f, 0.f);
368 void tap_distortion::deactivate()
370 is_active = false;
373 void tap_distortion::set_params(float blend, float drive)
375 // set distortion coeffs
376 if ((drive_old != drive) || (blend_old != blend)) {
377 rdrive = 12.0f / drive;
378 rbdr = rdrive / (10.5f - blend) * 780.0f / 33.0f;
379 kpa = D(2.0f * (rdrive*rdrive) - 1.0f) + 1.0f;
380 kpb = (2.0f - kpa) / 2.0f;
381 ap = ((rdrive*rdrive) - kpa + 1.0f) / 2.0f;
382 kc = kpa / D(2.0f * D(2.0f * (rdrive*rdrive) - 1.0f) - 2.0f * rdrive*rdrive);
384 srct = (0.1f * srate) / (0.1f * srate + 1.0f);
385 sq = kc*kc + 1.0f;
386 knb = -1.0f * rbdr / D(sq);
387 kna = 2.0f * kc * rbdr / D(sq);
388 an = rbdr*rbdr / sq;
389 imr = 2.0f * knb + D(2.0f * kna + 4.0f * an - 1.0f);
390 pwrq = 2.0f / (imr + 1.0f);
392 drive_old = drive;
393 blend_old = blend;
397 void tap_distortion::set_sample_rate(uint32_t sr)
399 srate = sr;
400 over = srate * 2 > 96000 ? 1 : 2;
401 resampler.set_params(srate, over, 2);
404 float tap_distortion::process(float in)
406 double *samples = resampler.upsample((double)in);
407 meter = 0.f;
408 for (int o = 0; o < over; o++) {
409 float proc = samples[o];
410 float med;
411 if (proc >= 0.0f) {
412 med = (D(ap + proc * (kpa - proc)) + kpb) * pwrq;
413 } else {
414 med = (D(an - proc * (kna + proc)) + knb) * pwrq * -1.0f;
416 proc = srct * (med - prev_med + prev_out);
417 prev_med = M(med);
418 prev_out = M(proc);
419 samples[o] = proc;
420 meter = std::max(meter, proc);
422 float out = (float)resampler.downsample(samples);
423 return out;
426 float tap_distortion::get_distortion_level()
428 return meter;
431 ////////////////////////////////////////////////////////////////////////////////
433 simple_lfo::simple_lfo()
435 is_active = false;
436 phase = 0.f;
439 void simple_lfo::activate()
441 is_active = true;
442 phase = 0.f;
445 void simple_lfo::deactivate()
447 is_active = false;
450 float simple_lfo::get_value()
452 return get_value_from_phase(phase, offset) * amount;
455 float simple_lfo::get_value_from_phase(float ph, float off) const
457 float val = 0.f;
458 float phs = ph + off;
459 if (phs >= 1.0)
460 phs = fmod(phs, 1.f);
461 switch (mode) {
462 default:
463 case 0:
464 // sine
465 val = sin((phs * 360.f) * M_PI / 180);
466 break;
467 case 1:
468 // triangle
469 if(phs > 0.75)
470 val = (phs - 0.75) * 4 - 1;
471 else if(phs > 0.5)
472 val = (phs - 0.5) * 4 * -1;
473 else if(phs > 0.25)
474 val = 1 - (phs - 0.25) * 4;
475 else
476 val = phs * 4;
477 break;
478 case 2:
479 // square
480 val = (phs < 0.5) ? -1 : +1;
481 break;
482 case 3:
483 // saw up
484 val = phs * 2.f - 1;
485 break;
486 case 4:
487 // saw down
488 val = 1 - phs * 2.f;
489 break;
491 return val;
494 void simple_lfo::advance(uint32_t count)
496 //this function walks from 0.f to 1.f and starts all over again
497 phase += count * freq * (1.0 / srate);
498 if (phase >= 1.0)
499 phase = fmod(phase, 1.f);
502 void simple_lfo::set_phase(float ph)
504 //set the phase from outsinde
505 phase = fabs(ph);
506 if (phase >= 1.0)
507 phase = fmod(phase, 1.f);
510 void simple_lfo::set_params(float f, int m, float o, uint32_t sr, float a)
512 // freq: a value in Hz
513 // mode: sine=0, triangle=1, square=2, saw_up=3, saw_down=4
514 // offset: value between 0.f and 1.f to offset the lfo in time
515 freq = f;
516 mode = m;
517 offset = o;
518 srate = sr;
519 amount = a;
521 void simple_lfo::set_freq(float f)
523 freq = f;
525 void simple_lfo::set_mode(int m)
527 mode = m;
529 void simple_lfo::set_offset(float o)
531 offset = o;
533 void simple_lfo::set_amount(float a)
535 amount = a;
537 bool simple_lfo::get_graph(float *data, int points, cairo_iface *context, int *mode) const
539 if (!is_active)
540 return false;
541 for (int i = 0; i < points; i++) {
542 float ph = (float)i / (float)points;
543 data[i] = get_value_from_phase(ph, offset) * amount;
545 return true;
548 bool simple_lfo::get_dot(float &x, float &y, int &size, cairo_iface *context) const
550 if (!is_active)
551 return false;
552 float phs = phase + offset;
553 if (phs >= 1.0)
554 phs = fmod(phs, 1.f);
555 x = phase;
556 y = get_value_from_phase(phase, offset) * amount;
557 return true;
561 /// Lookahead Limiter by Christian Holschuh and Markus Schmidt
563 lookahead_limiter::lookahead_limiter() {
564 is_active = false;
565 channels = 2;
566 id = 0;
567 buffer_size = 0;
568 overall_buffer_size = 0;
569 att = 1.f;
570 att_max = 1.0;
571 pos = 0;
572 delta = 0.f;
573 _delta = 0.f;
574 peak = 0.f;
575 over_s = 0;
576 over_c = 1.f;
577 attack = 0.005;
578 use_multi = false;
579 weight = 1.f;
580 _sanitize = false;
581 auto_release = false;
582 asc_active = false;
583 nextiter = 0;
584 nextlen = 0;
585 asc = 0.f;
586 asc_c = 0;
587 asc_pos = -1;
588 asc_changed = false;
589 asc_coeff = 1.f;
591 lookahead_limiter::~lookahead_limiter()
593 free(buffer);
594 free(nextpos);
595 free(nextdelta);
598 void lookahead_limiter::activate()
600 is_active = true;
601 pos = 0;
605 void lookahead_limiter::set_multi(bool set) { use_multi = set; }
607 void lookahead_limiter::deactivate()
609 is_active = false;
612 float lookahead_limiter::get_attenuation()
614 float a = att_max;
615 att_max = 1.0;
616 return a;
619 void lookahead_limiter::set_sample_rate(uint32_t sr)
621 srate = sr;
622 // rebuild buffer
623 overall_buffer_size = (int)(srate * (100.f / 1000.f) * channels) + channels; // buffer size attack rate multiplied by 2 channels
624 buffer = (float*) calloc(overall_buffer_size, sizeof(float));
625 pos = 0;
627 nextdelta = (float*) calloc(overall_buffer_size, sizeof(float));
628 nextpos = (int*) malloc(overall_buffer_size * sizeof(int));
629 memset(nextpos, -1, overall_buffer_size * sizeof(int));
631 reset();
634 void lookahead_limiter::set_params(float l, float a, float r, float w, bool ar, float arc, bool d)
636 limit = l;
637 attack = a / 1000.f;
638 release = r / 1000.f;
639 auto_release = ar;
640 asc_coeff = arc;
641 debug = d;
642 weight = w;
645 void lookahead_limiter::reset() {
646 int bs = (int)(srate * attack * channels);
647 buffer_size = bs - bs % channels; // buffer size attack rate
648 _sanitize = true;
649 pos = 0;
650 nextpos[0] = -1;
651 nextlen = 0;
652 nextiter = 0;
653 delta = 0.f;
654 att = 1.f;
655 reset_asc();
658 void lookahead_limiter::reset_asc() {
659 asc = 0.f;
660 asc_c = 0;
661 asc_pos = pos;
662 asc_changed = true;
665 float lookahead_limiter::get_rdelta(float peak, float _limit, float _att, bool _asc) {
667 // calc the att for average input to walk to if we use asc (att of average signal)
668 float _a_att = (limit * weight) / (asc_coeff * asc) * (float)asc_c;
670 // calc a release delta from this attenuation
671 float _rdelta = (1.0 - _att) / (srate * release);
672 if(_asc and auto_release and asc_c > 0 and _a_att > _att) {
673 // check if releasing to average level of peaks is steeper than
674 // releasing to 1.f
675 float _delta = std::max((_a_att - _att) / (srate * release), _rdelta / 10);
676 if(_delta < _rdelta) {
677 asc_active = true;
678 _asc_used = true;
679 _rdelta = _delta;
682 return _rdelta;
685 void lookahead_limiter::process(float &left, float &right, float * multi_buffer)
687 // PROTIP: harming paying customers enough to make them develop a competing
688 // product may be considered an example of a less than sound business practice.
690 // fill lookahead buffer
691 if(_sanitize) {
692 // if we're sanitizing (zeroing) the buffer on attack time change,
693 // don't write the samples to the buffer
694 buffer[pos] = 0.f;
695 buffer[pos + 1] = 0.f;
696 } else {
697 buffer[pos] = left;
698 buffer[pos + 1] = right;
701 // are we using multiband? get the multiband coefficient or use 1.f
702 float multi_coeff = (use_multi) ? multi_buffer[pos] : 1.f;
704 // calc the real limit including weight and multi coeff
705 float _limit = limit * multi_coeff * weight;
707 // input peak - impact higher in left or right channel?
708 peak = fabs(left) > fabs(right) ? fabs(left) : fabs(right);
710 // add an eventually appearing peak to the asc fake buffer if asc active
711 if(auto_release and peak > _limit) {
712 asc += peak;
713 asc_c ++;
716 if(peak > _limit or multi_coeff < 1.0) {
717 float _multi_coeff = 1.f;
718 float _peak;
720 // calc the attenuation needed to reduce incoming peak
721 float _att = std::min(_limit / peak, 1.f);
722 // calc release without any asc to keep all relevant peaks
723 float _rdelta = get_rdelta(peak, _limit, _att, false);
725 // calc the delta for walking to incoming peak attenuation
726 float _delta = (_limit / peak - att) / buffer_size * channels;
728 if(_delta < delta) {
729 // is the delta more important than the actual one?
730 // if so, we can forget about all stored deltas (because they can't
731 // be more important - we already checked that earlier) and use this
732 // delta now. and we have to create a release delta in nextpos buffer
733 nextpos[0] = pos;
734 nextpos[1] = -1;
735 nextdelta[0] = _rdelta;
736 nextlen = 1;
737 nextiter = 0;
738 delta = _delta;
739 } else {
740 // we have a peak on input its delta is less important than the
741 // actual delta. But what about the stored deltas we're following?
742 bool _found = false;
743 int i = 0;
744 for(i = nextiter; i < nextiter + nextlen; i++) {
745 // walk through our nextpos buffer
746 int j = i % buffer_size;
747 // calculate a delta for the next stored peak
748 // are we using multiband? then get the multi_coeff for the
749 // stored position
750 _multi_coeff = (use_multi) ? multi_buffer[nextpos[j]] : 1.f;
751 // is the left or the right channel on this position more
752 // important?
753 _peak = fabs(buffer[nextpos[j]]) > fabs(buffer[nextpos[j] + 1]) ? fabs(buffer[nextpos[j]]) : fabs(buffer[nextpos[j] + 1]);
754 // calc a delta to use to reach our incoming peak from the
755 // stored position
756 _delta = (_limit / peak - (limit * _multi_coeff * weight) / _peak) / (((buffer_size - nextpos[j] + pos) % buffer_size) / channels);
757 if(_delta < nextdelta[j]) {
758 // if the buffered delta is more important than the delta
759 // used to reach our peak from the stored position, store
760 // the new delta at that position and stop the loop
761 nextdelta[j] = _delta;
762 _found = true;
763 break;
766 if(_found) {
767 // there was something more important in the next-buffer.
768 // throw away any position and delta after the important
769 // position and add a new release delta
770 nextlen = i - nextiter + 1;
771 nextpos[(nextiter + nextlen) % buffer_size] = pos;
772 nextdelta[(nextiter + nextlen) % buffer_size] = _rdelta;
773 // set the next following position value to -1 (cleaning up the
774 // nextpos buffer)
775 nextpos[(nextiter + nextlen + 1) % buffer_size] = -1;
776 // and raise the length of our nextpos buffer for keeping the
777 // release value
778 nextlen ++;
783 // switch left and right pointers in buffer to output position
784 left = buffer[(pos + channels) % buffer_size];
785 right = buffer[(pos + channels + 1) % buffer_size];
787 // if a peak leaves the buffer, remove it from asc fake buffer
788 // but only if we're not sanitizing asc buffer
789 float _peak = fabs(left) > fabs(right) ? fabs(left) : fabs(right);
790 float _multi_coeff = (use_multi) ? multi_buffer[(pos + channels) % buffer_size] : 1.f;
791 if(pos == asc_pos and !asc_changed) {
792 asc_pos = -1;
794 if(auto_release and asc_pos == -1 and _peak > (limit * weight * _multi_coeff)) {
795 asc -= _peak;
796 asc_c --;
799 // change the attenuation level
800 att += delta;
802 // ...and calculate outpout from it
803 left *= att;
804 right *= att;
806 if((pos + channels) % buffer_size == nextpos[nextiter]) {
807 // if we reach a buffered position, change the actual delta and erase
808 // this (the first) element from nextpos and nextdelta buffer
809 if(auto_release) {
810 // set delta to asc influenced release delta
811 delta = get_rdelta(_peak, (limit * weight * _multi_coeff), att);
812 if(nextlen > 1) {
813 // if there are more positions to walk to, calc delta to next
814 // position in buffer and compare it to release delta (keep
815 // changes between peaks below asc steepness)
816 int _nextpos = nextpos[(nextiter + 1) % buffer_size];
817 float __peak = fabs(buffer[_nextpos]) > fabs(buffer[_nextpos + 1]) ? fabs(buffer[_nextpos]) : fabs(buffer[_nextpos + 1]);
818 float __multi_coeff = (use_multi) ? multi_buffer[_nextpos] : 1.f;
819 float __delta = ((limit * __multi_coeff * weight) / __peak - att) / (((buffer_size + _nextpos - ((pos + channels) % buffer_size)) % buffer_size) / channels);
820 if(__delta < delta) {
821 delta = __delta;
824 } else {
825 // if no asc set delta from nextdelta buffer and fix the attenuation
826 delta = nextdelta[nextiter];
827 att = (limit * weight * _multi_coeff) / _peak;
829 // remove first element from circular nextpos buffer
830 nextlen -= 1;
831 nextpos[nextiter] = -1;
832 nextiter = (nextiter + 1) % buffer_size;
835 if (att > 1.0f) {
836 // release time seems over, reset attenuation and delta
837 att = 1.0f;
838 delta = 0.0f;
839 nextiter = 0;
840 nextlen = 0;
841 nextpos[0] = -1;
844 // main limiting party is over, let's cleanup the puke
846 if(_sanitize) {
847 // we're sanitizing? then send 0.f as output
848 left = 0.f;
849 right = 0.f;
852 // security personnel pawing your values
853 if(att <= 0.f) {
854 // if this happens we're doomed!!
855 // may happen on manually lowering attack
856 att = 0.0000000000001;
857 delta = (1.0f - att) / (srate * release);
860 if(att != 1.f and 1 - att < 0.0000000000001) {
861 // denormalize att
862 att = 1.f;
865 if(delta != 0.f and fabs(delta) < 0.00000000000001) {
866 // denormalize delta
867 delta = 0.f;
870 // post treatment (denormal, limit)
871 denormal(&left);
872 denormal(&right);
874 // store max attenuation for meter output
875 att_max = (att < att_max) ? att : att_max;
877 // step forward in our sample ring buffer
878 pos = (pos + channels) % buffer_size;
880 // sanitizing is always done after a full cycle through the lookahead buffer
881 if(_sanitize and pos == 0) _sanitize = false;
883 asc_changed = false;
886 bool lookahead_limiter::get_asc() {
887 if(!asc_active) return false;
888 asc_active = false;
889 return true;
893 ////////////////////////////////////////////////////////////////////////////////
895 transients::transients() {
896 envelope = 0.f;
897 attack = 0.f;
898 release = 0.f;
899 attack_coef = 0.f;
900 release_coef = 0.f;
901 att_time = 0.f;
902 att_level = 0.f;
903 rel_time = 0.f;
904 rel_level = 0.f;
905 sust_thres = 1.f;
906 maxdelta = 0.f;
907 new_return = 1.f;
908 old_return = 1.f;
909 lookahead = 0;
910 lookpos = 0;
911 channels = 1;
912 sustain_ended = false;
914 transients::~transients()
916 free(lookbuf);
918 void transients::set_channels(int ch) {
919 channels = ch;
920 lookbuf = (float*) calloc(looksize * channels, sizeof(float));
921 lookpos = 0;
923 void transients::set_sample_rate(uint32_t sr) {
924 srate = sr;
925 attack_coef = exp(log(0.01) / (0.001 * srate));
926 release_coef = exp(log(0.01) / (0.2f * srate));
927 // due to new calculation in attack, we sometimes get harsh
928 // gain reduction/boost.
929 // to prevent "clicks" a maxdelta is set, which allows the signal
930 // to raise/fall ~6dB/ms.
931 maxdelta = pow(4, 1.f / (0.001 * srate));
932 calc_relfac();
934 void transients::set_params(float att_t, float att_l, float rel_t, float rel_l, float sust_th, int look) {
935 lookahead = look;
936 sust_thres = sust_th;
937 att_time = att_t;
938 rel_time = rel_t;
939 att_level = att_l > 0 ? 0.25f * pow(att_l * 8, 2)
940 : -0.25f * pow(att_l * 4, 2);
941 rel_level = rel_l > 0 ? 0.5f * pow(rel_l * 8, 2)
942 : -0.25f * pow(rel_l * 4, 2);
943 calc_relfac();
945 void transients::calc_relfac()
947 relfac = pow(0.5f, 1.f / (0.001 * rel_time * srate));
949 void transients::process(float *in, float s) {
950 s = fabs(s);
951 // fill lookahead buffer
952 for (int i = 0; i < channels; i++) {
953 lookbuf[lookpos + i] = in[i];
956 // envelope follower
957 // this is the real envelope follower curve. It raises as
958 // fast as the signal is raising and falls much slower
959 // depending on the sample rate and the ffactor
960 // (the falling factor)
961 if(s > envelope)
962 envelope = attack_coef * (envelope - s) + s;
963 else
964 envelope = release_coef * (envelope - s) + s;
966 // attack follower
967 // this is a curve which follows the envelope slowly.
968 // It never can rise above the envelope. It reaches 70.7%
969 // of the envelope in a certain amount of time set by the user
971 double attdelta = (envelope - attack)
972 * 0.707
973 / (srate * att_time * 0.001);
974 if (sustain_ended == true and envelope / attack - 1 > 0.2f)
975 sustain_ended = false;
976 attack += attdelta;
978 // never raise above envelope
979 attack = std::min(envelope, attack);
981 // release follower
982 // this is a curve which is always above the envelope. It
983 // starts to fall when the envelope falls beneath the
984 // sustain threshold
986 if ((envelope / release) - sust_thres < 0 and sustain_ended == false)
987 sustain_ended = true;
988 double reldelta = sustain_ended ? relfac : 1;
990 // release delta can never raise above 0
991 release *= reldelta;
993 // never fall below envelope
994 release = std::max(envelope, release);
996 // difference between attack and envelope
997 double attdiff = attack > 0 ? log(envelope / attack) : 0;
999 // difference between release and envelope
1000 double reldiff = envelope > 0 ? log(release / envelope) : 0;
1002 // amplification factor from attack and release curve
1003 double ampfactor = attdiff * att_level + reldiff * rel_level;
1004 old_return = new_return;
1005 new_return = 1 + (ampfactor < 0 ? std::max(-1 + 1e-15, exp(ampfactor) - 1) : ampfactor);
1006 if (new_return / old_return > maxdelta)
1007 new_return = old_return * maxdelta;
1008 else if (new_return / old_return < 1 / maxdelta)
1009 new_return = old_return / maxdelta;
1011 int pos = (lookpos + looksize * channels - lookahead * channels) % (looksize * channels);
1013 for (int i = 0; i < channels; i++)
1014 in[i] = lookbuf[pos + i] * new_return;
1016 // advance lookpos
1017 lookpos = (lookpos + channels) % (looksize * channels);
1021 //////////////////////////////////////////////////////////////////
1023 crossover::crossover() {
1024 bands = -1;
1025 mode = -1;
1026 redraw_graph = 1;
1028 void crossover::set_sample_rate(uint32_t sr) {
1029 srate = sr;
1031 void crossover::init(int c, int b, uint32_t sr) {
1032 channels = std::min(8, c);
1033 bands = std::min(8, b);
1034 srate = sr;
1035 for(int b = 0; b < bands; b ++) {
1036 // reset frequency settings
1037 freq[b] = 1.0;
1038 active[b] = true;
1039 level[b] = 1.0;
1040 for (int c = 0; c < channels; c ++) {
1041 // reset outputs
1042 out[c][b] = 0.f;
1046 float crossover::set_filter(int b, float f, bool force) {
1047 // keep between neighbour bands
1048 if (b)
1049 f = std::max((float)freq[b-1] * 1.1f, f);
1050 if (b < bands - 2)
1051 f = std::min((float)freq[b+1] * 0.9f, f);
1052 // restrict to 10-20k
1053 f = std::max(10.f, std::min(20000.f, f));
1054 // nothing changed? return
1055 if (freq[b] == f and !force)
1056 return freq[b];
1057 freq[b] = f;
1058 float q;
1059 switch (mode) {
1060 case 0:
1061 default:
1062 q = 0.5;
1063 break;
1064 case 1:
1065 q = 0.7071068123730965;
1066 break;
1067 case 2:
1068 q = 0.54;
1069 break;
1071 for (int c = 0; c < channels; c ++) {
1072 if (!c) {
1073 lp[c][b][0].set_lp_rbj(freq[b], q, (float)srate);
1074 hp[c][b][0].set_hp_rbj(freq[b], q, (float)srate);
1075 } else {
1076 lp[c][b][0].copy_coeffs(lp[c-1][b][0]);
1077 hp[c][b][0].copy_coeffs(hp[c-1][b][0]);
1079 if (mode > 1) {
1080 if (!c) {
1081 lp[c][b][1].set_lp_rbj(freq[b], 1.34, (float)srate);
1082 hp[c][b][1].set_hp_rbj(freq[b], 1.34, (float)srate);
1083 } else {
1084 lp[c][b][1].copy_coeffs(lp[c-1][b][1]);
1085 hp[c][b][1].copy_coeffs(hp[c-1][b][1]);
1087 lp[c][b][2].copy_coeffs(lp[c][b][0]);
1088 hp[c][b][2].copy_coeffs(hp[c][b][0]);
1089 lp[c][b][3].copy_coeffs(lp[c][b][1]);
1090 hp[c][b][3].copy_coeffs(hp[c][b][1]);
1091 } else {
1092 lp[c][b][1].copy_coeffs(lp[c][b][0]);
1093 hp[c][b][1].copy_coeffs(hp[c][b][0]);
1096 redraw_graph = std::min(2, redraw_graph + 1);
1097 return freq[b];
1099 void crossover::set_mode(int m) {
1100 if(mode == m)
1101 return;
1102 mode = m;
1103 for(int i = 0; i < bands - 1; i ++) {
1104 set_filter(i, freq[i], true);
1106 redraw_graph = std::min(2, redraw_graph + 1);
1108 void crossover::set_active(int b, bool a) {
1109 if (active[b] == a)
1110 return;
1111 active[b] = a;
1112 redraw_graph = std::min(2, redraw_graph + 1);
1114 void crossover::set_level(int b, float l) {
1115 if (level[b] == l)
1116 return;
1117 level[b] = l;
1118 redraw_graph = std::min(2, redraw_graph + 1);
1120 void crossover::process(float *data) {
1121 for (int c = 0; c < channels; c++) {
1122 for(int b = 0; b < bands; b ++) {
1123 out[c][b] = data[c];
1124 for (int f = 0; f < get_filter_count(); f++){
1125 if(b + 1 < bands) {
1126 out[c][b] = lp[c][b][f].process(out[c][b]);
1127 lp[c][b][f].sanitize();
1129 if(b - 1 >= 0) {
1130 out[c][b] = hp[c][b - 1][f].process(out[c][b]);
1131 hp[c][b - 1][f].sanitize();
1134 out[c][b] *= level[b];
1138 float crossover::get_value(int c, int b) {
1139 return out[c][b];
1141 bool crossover::get_graph(int subindex, int phase, float *data, int points, cairo_iface *context, int *mode) const
1143 if (subindex >= bands) {
1144 redraw_graph = std::max(0, redraw_graph - 1);
1145 return false;
1147 float ret;
1148 double freq;
1149 for (int i = 0; i < points; i++) {
1150 ret = 1.f;
1151 freq = 20.0 * pow (20000.0 / 20.0, i * 1.0 / points);
1152 for(int f = 0; f < get_filter_count(); f ++) {
1153 if(subindex < bands -1)
1154 ret *= lp[0][subindex][f].freq_gain(freq, (float)srate);
1155 if(subindex > 0)
1156 ret *= hp[0][subindex - 1][f].freq_gain(freq, (float)srate);
1158 ret *= level[subindex];
1159 context->set_source_rgba(0.15, 0.2, 0.0, !active[subindex] ? 0.3 : 0.8);
1160 data[i] = dB_grid(ret);
1162 return true;
1164 bool crossover::get_layers(int index, int generation, unsigned int &layers) const
1166 layers = 0 | (redraw_graph or !generation ? LG_CACHE_GRAPH : 0) | (!generation ? LG_CACHE_GRID : 0);
1167 return redraw_graph or !generation;
1170 int crossover::get_filter_count() const
1172 switch (mode) {
1173 case 0:
1174 default:
1175 return 1;
1176 case 1:
1177 return 2;
1178 case 2:
1179 return 4;
1183 //////////////////////////////////////////////////////////////////
1185 bitreduction::bitreduction()
1187 coeff = 1;
1188 morph = 0;
1189 mode = 0;
1190 dc = 0;
1191 sqr = 0;
1192 aa = 0;
1193 aa1 = 0;
1194 redraw_graph = true;
1195 bypass = true;
1197 void bitreduction::set_sample_rate(uint32_t sr)
1199 srate = sr;
1201 void bitreduction::set_params(float b, float mo, bool bp, uint32_t md, float d, float a)
1203 morph = 1 - mo;
1204 bypass = bp;
1205 dc = d;
1206 aa = a;
1207 mode = md;
1208 coeff = powf(2.0f, b) - 1;
1209 sqr = sqrt(coeff / 2);
1210 aa1 = (1.f - aa) / 2.f;
1211 redraw_graph = true;
1213 float bitreduction::add_dc(float s, float dc) const
1215 return s > 0 ? s *= dc : s /= dc;
1217 float bitreduction::remove_dc(float s, float dc) const
1219 return s > 0 ? s /= dc : s *= dc;
1221 float bitreduction::waveshape(float in) const
1223 double y;
1224 double k;
1226 // add dc
1227 in = add_dc(in, dc);
1229 // main rounding calculation depending on mode
1231 // the idea for anti-aliasing:
1232 // you need a function f which brings you to the scale, where you want to round
1233 // and the function f_b (with f(f_b)=id) which brings you back to your original scale.
1235 // then you can use the logic below in the following way:
1236 // y = f(in) and k = roundf(y)
1237 // if (y > k + aa1)
1238 // k = f_b(k) + ( f_b(k+1) - f_b(k) ) *0.5 * (sin(x - PI/2) + 1)
1239 // if (y < k + aa1)
1240 // k = f_b(k) - ( f_b(k+1) - f_b(k) ) *0.5 * (sin(x - PI/2) + 1)
1242 // whereas x = (fabs(f(in) - k) - aa1) * PI / aa
1243 // for both cases.
1245 switch (mode) {
1246 case 0:
1247 default:
1248 // linear
1249 y = (in) * coeff;
1250 k = roundf(y);
1251 if (k - aa1 <= y and y <= k + aa1) {
1252 k /= coeff;
1253 } else if (y > k + aa1) {
1254 k = k / coeff + ((k + 1) / coeff - k / coeff) * 0.5 * (sin(M_PI * (fabs(y - k) - aa1) / aa - M_PI_2) + 1);
1255 } else {
1256 k = k / coeff - (k / coeff - (k - 1) / coeff) * 0.5 * (sin(M_PI * (fabs(y - k) - aa1) / aa - M_PI_2) + 1);
1258 break;
1259 case 1:
1260 // logarithmic
1261 y = sqr * log(fabs(in)) + sqr * sqr;
1262 k = roundf(y);
1263 if(!in) {
1264 k = 0;
1265 } else if (k - aa1 <= y and y <= k + aa1) {
1266 k = in / fabs(in) * exp(k / sqr - sqr);
1267 } else if (y > k + aa1) {
1268 k = in / fabs(in) * (exp(k / sqr - sqr) + (exp((k + 1) / sqr - sqr) - exp(k / sqr - sqr)) * 0.5 * (sin((fabs(y - k) - aa1) / aa * M_PI - M_PI_2) + 1));
1269 } else {
1270 k = in / fabs(in) * (exp(k / sqr - sqr) - (exp(k / sqr - sqr) - exp((k - 1) / sqr - sqr)) * 0.5 * (sin((fabs(y - k) - aa1) / aa * M_PI - M_PI_2) + 1));
1272 break;
1275 // morph between dry and wet signal
1276 k += (in - k) * morph;
1278 // remove dc
1279 k = remove_dc(k, dc);
1281 return k;
1283 float bitreduction::process(float in)
1285 return waveshape(in);
1288 bool bitreduction::get_layers(int index, int generation, unsigned int &layers) const
1290 layers = redraw_graph or !generation ? LG_CACHE_GRAPH | LG_CACHE_GRID : 0;
1291 return redraw_graph or !generation;
1293 bool bitreduction::get_graph(int subindex, int phase, float *data, int points, cairo_iface *context, int *mode) const
1295 if (subindex >= 2) {
1296 redraw_graph = false;
1297 return false;
1299 for (int i = 0; i < points; i++) {
1300 data[i] = sin(((float)i / (float)points * 360.) * M_PI / 180.);
1301 if (subindex and !bypass)
1302 data[i] = waveshape(data[i]);
1303 else {
1304 context->set_line_width(1);
1305 context->set_source_rgba(0.15, 0.2, 0.0, 0.15);
1308 return true;
1310 bool bitreduction::get_gridline(int subindex, int phase, float &pos, bool &vertical, std::string &legend, cairo_iface *context) const
1312 if (phase or subindex)
1313 return false;
1314 pos = 0;
1315 vertical = false;
1316 return true;
1319 //////////////////////////////////////////////////////////////////
1321 resampleN::resampleN()
1323 factor = 2;
1324 srate = 0;
1325 filters = 2;
1327 resampleN::~resampleN()
1329 free(tmp);
1331 void resampleN::set_params(uint32_t sr, int fctr = 2, int fltrs = 2)
1333 srate = sr;
1334 factor = std::min(16, std::max(1, fctr));
1335 filters = std::min(4, std::max(1, fltrs));
1336 // set all filters
1337 filter[0][0].set_lp_rbj(std::max(25000., (double)srate / 2), 0.8, (float)srate * factor);
1338 for (int i = 1; i < filters; i++) {
1339 filter[0][i].copy_coeffs(filter[0][0]);
1340 filter[1][i].copy_coeffs(filter[0][0]);
1343 double *resampleN::upsample(double sample)
1345 tmp[0] = sample;
1346 if (factor > 1) {
1347 for (int f = 0; f < filters; f++)
1348 tmp[0] = filter[0][f].process(sample);
1349 for (int i = 1; i < factor; i++) {
1350 tmp[i] = 0;
1351 for (int f = 0; f < filters; f++)
1352 tmp[i] = filter[0][f].process(sample);
1355 return tmp;
1357 double resampleN::downsample(double *sample)
1359 if (factor > 1) {
1360 for(int i = 0; i < factor; i++) {
1361 for (int f = 0; f < filters; f++) {
1362 sample[i] = filter[1][f].process(sample[i]);
1366 return sample[0];
1369 //////////////////////////////////////////////////////////////////
1371 samplereduction::samplereduction()
1373 target = 0;
1374 real = 0;
1375 samples = 0;
1376 last = 0;
1377 amount = 0;
1378 round = 0;
1380 void samplereduction::set_params(float am)
1382 amount = am;
1383 round = roundf(amount);
1384 //samples = round;
1386 double samplereduction::process(double in)
1388 samples ++;
1389 if (samples >= round) {
1390 target += amount;
1391 real += round;
1392 if (target + amount >= real + 1) {
1393 last = in;
1394 target = 0;
1395 real = 0;
1397 samples = 0;
1399 return last;