Anklang 0.3.0-460-gc4ef46ba
ASE — Anklang Sound Engine (C++)

« « « Anklang Documentation
Loading...
Searching...
No Matches
bleposc.hh
Go to the documentation of this file.
1 // This Source Code Form is licensed MPL-2.0: http://mozilla.org/MPL/2.0
2#ifndef __ASE_DEVICES_BLEPOSC_HH__
3#define __ASE_DEVICES_BLEPOSC_HH__
4
5#include <ase/randomhash.hh>
6#include <ase/datautils.hh>
7
8namespace Ase {
9namespace BlepUtils {
10
11using std::max;
12
14{
15 double rate_;
16 double leaky_a;
17
18 static const int WIDTH = 13;
19 static const int WSHIFT = 6; // delay to align saw and impulse part of the signal
20 static const int OVERSAMPLE = 64;
21
22 static const float blep_table[WIDTH * OVERSAMPLE + 1];
23
24public:
25 double frequency_base = 440;
26 double frequency_factor = 1;
27
28 double freq_mod_octaves = 0;
29
30 double shape_base = 0; // 0 = saw, range [-1:1]
31 double shape_mod = 1.0;
32
33 double pulse_width_base = 0.5;
34 double pulse_width_mod = 0.0;
35
36 double sync_base = 0;
37 double sync_mod = 0;
38
39 double sub_base = 0;
40 double sub_mod = 0;
41
42 double sub_width_base = 0.5;
43 double sub_width_mod = 0.0;
44
45 bool need_reset_voice_state;
46
47 enum class State {
48 A,
49 B,
50 C,
51 D
52 };
53
55 {
56 double freq_factor = 1;
57 double left_factor = 1;
58 double right_factor = 0;
59
60 double master_phase = 0;
61 double slave_phase = 0;
62
63 double last_value = 0; /* leaky integrator state */
64 double current_level = 0; /* current position of the wave form (saw + jumps) */
65
66 double last_dc = 0; /* dc of previous parameters */
67 double dc_delta = 0;
68 int dc_steps = 0;
69
70 int future_pos = 0;
71
72
73 State state = State::A;
74
75 // could be shared under certain conditions
77
78 void
79 init_future()
80 {
81 future.fill (0);
82 future_pos = 0;
83 }
84 float
85 pop_future()
86 {
87 const float f = future[future_pos++];
88 if (future_pos == WIDTH)
89 {
90 // this loop was slightly faster than std::copy_n + std::fill_n
91 for (int i = 0; i < WIDTH; i++)
92 {
93 future[i] = future[WIDTH + i];
94 future[WIDTH + i] = 0;
95 }
96 future_pos = 0;
97 }
98 return f;
99 }
100 };
101
102 std::vector<UnisonVoice> unison_voices;
103
104 OscImpl()
105 {
106 set_unison (1, 0, 0); // default
107 }
108 void
109 reset()
110 {
111 const bool randomize_phase = unison_voices.size() > 1;
112
113 for (auto& voice : unison_voices)
114 {
115 voice.init_future();
116 voice.dc_steps = 0;
117 voice.dc_delta = 0;
118
119 if (randomize_phase) // randomize start phase for true unison
120 {
121 reset_master (voice, random_frange (0, 1));
122 }
123 else
124 {
125 reset_master (voice, 0);
126 }
127 }
128 }
129 void
130 reset_master (UnisonVoice& voice, double master_phase)
131 {
132 voice.master_phase = master_phase;
133 need_reset_voice_state = true;
134 }
135 void
136 set_unison (size_t n_voices, float detune, float stereo)
137 {
138 const bool unison_voices_changed = unison_voices.size() != n_voices;
139
140 unison_voices.resize (n_voices);
141
142 bool left_channel = true; /* start spreading voices at the left channel */
143 for (size_t i = 0; i < unison_voices.size(); i++)
144 {
145 if (n_voices == 1)
146 unison_voices[i].freq_factor = 1;
147 else
148 {
149 const float detune_cent = -detune / 2.0 + i / float (n_voices - 1) * detune;
150 unison_voices[i].freq_factor = pow (2, detune_cent / 1200);
151 }
152 /* stereo spread factors */
153 double left_factor, right_factor;
154 bool odd_n_voices = unison_voices.size() & 1;
155 if (odd_n_voices && i == unison_voices.size() / 2) // odd number of voices: this voice is centered
156 {
157 left_factor = (1 - stereo) + stereo * 0.5;
158 right_factor = (1 - stereo) + stereo * 0.5;
159 }
160 else if (left_channel) // alternate beween left and right voices
161 {
162 left_factor = 0.5 + stereo / 2;
163 right_factor = 0.5 - stereo / 2;
164 left_channel = false;
165 }
166 else
167 {
168 left_factor = 0.5 - stereo / 2;
169 right_factor = 0.5 + stereo / 2;
170 left_channel = true;
171 }
172 /* ensure constant total energy of left + right channel combined
173 *
174 * also take into account the more unison voices we add up, the louder the result
175 * will be, so compensate for this:
176 *
177 * -> each time the number of voices is doubled, the signal level is increased by
178 * a factor of sqrt (2)
179 */
180 const double norm = sqrt (left_factor * left_factor + right_factor * right_factor) * sqrt (n_voices / 2.0);
181 unison_voices[i].left_factor = left_factor / norm;
182 unison_voices[i].right_factor = right_factor / norm;
183 }
184 if (unison_voices_changed)
185 reset();
186 }
187 void
188 set_rate (double rate)
189 {
190 rate_ = rate;
191
192 /* get leaky integrator constant for sample rate from ms (half-life time) */
193 const double leaky_ms = 10;
194 leaky_a = pow (2.0, -1000.0 / (rate_ * leaky_ms));
195 }
196 double
197 rate()
198 {
199 return rate_;
200 }
201 double
202 estimate_dc (double shape,
203 double pulse_width,
204 double sub,
205 double sub_width,
206 double sync_factor)
207 {
208 const double bound_a = sub_width * pulse_width;
209 const double bound_b = 2 * sub_width * pulse_width + 1 - sub_width - pulse_width;
210 const double bound_c = sub_width * pulse_width + (1 - sub_width);
211 const double bound_d = 1.0;
212
213 const double saw_slope = -4.0 * (shape + 1) * (1 - sub);
214
215 const double a1 = 1;
216 const double a2 = a1 + saw_slope * bound_a;
217
218 const double b1 = a2 + 2.0 * (shape * (1 - sub) - sub);
219 const double b2 = b1 + saw_slope * (bound_b - bound_a);
220
221 const double c1 = b2 + 2 * (1 - sub);
222 const double c2 = c1 + saw_slope * (bound_c - bound_b);
223
224 const double d1 = c2 + 2.0 * (shape * (1 - sub) + sub);
225 const double d2 = d1 + saw_slope * (bound_d - bound_c);
226
227 /* dc without sync */
228 const double dc_base = (a1 + a2) / 2 * bound_a
229 + (b1 + b2) / 2 * (bound_b - bound_a)
230 + (c1 + c2) / 2 * (bound_c - bound_b)
231 + (d1 + d2) / 2 * (bound_d - bound_c);
232
233 /* quick path: no sync, no sync related dc computation */
234 if (sync_factor < 1.01)
235 return dc_base;
236
237 /* dc offset introduced by sync */
238 const double sync_phase = sync_factor - int (sync_factor);
239
240 double a_avg = (a1 + a2) / 2;
241 double b_avg = (b1 + b2) / 2;
242 double c_avg = (c1 + c2) / 2;
243 double d_avg = (d1 + d2) / 2;
244
245 if (sync_phase < bound_a)
246 {
247 const double frac = (bound_a - sync_phase) / bound_a;
248 const double sync_a2 = a1 * frac + a2 * (1 - frac);
249
250 a_avg = (1 - frac) * (a1 + sync_a2) / 2;
251 b_avg = c_avg = d_avg = 0;
252 }
253 else if (sync_phase < bound_b)
254 {
255 const double frac = (bound_b - sync_phase) / (bound_b - bound_a);
256 const double sync_b2 = b1 * frac + b2 * (1 - frac);
257
258 b_avg = (1 - frac) * (b1 + sync_b2) / 2;
259 c_avg = d_avg = 0;
260 }
261 else if (sync_phase < bound_c)
262 {
263 const double frac = (bound_c - sync_phase) / (bound_c - bound_b);
264 const double sync_c2 = c1 * frac + c2 * (1 - frac);
265
266 c_avg = (1 - frac) * (c1 + sync_c2) / 2;
267 d_avg = 0;
268 }
269 else
270 {
271 const double frac = (bound_d - sync_phase) / (bound_d - bound_c);
272 const double sync_d2 = d1 * frac + d2 * (1 - frac);
273
274 d_avg = (1 - frac) * (d1 + sync_d2) / 2;
275 }
276
277 /* dc sync part of the signal */
278 const double dc_sync = a_avg * bound_a
279 + b_avg * (bound_b - bound_a)
280 + c_avg * (bound_c - bound_b)
281 + d_avg * (bound_d - bound_c);
282 const double dc = (dc_base * (int) sync_factor + dc_sync) / sync_factor;
283
284 return dc;
285 }
286
287 void
288 reset_voice_state (double shape,
289 double pulse_width,
290 double sub,
291 double sub_width,
292 double sync_factor)
293 {
294 const double bound_a = sub_width * pulse_width;
295 const double bound_b = 2 * sub_width * pulse_width + 1 - sub_width - pulse_width;
296 const double bound_c = sub_width * pulse_width + (1 - sub_width);
297 const double bound_d = 1.0;
298
299 const double saw_slope = -4.0 * (shape + 1) * (1 - sub);
300
301 const double a1 = 1;
302 const double a2 = a1 + saw_slope * bound_a;
303
304 const double b1 = a2 + 2.0 * (shape * (1 - sub) - sub);
305 const double b2 = b1 + saw_slope * (bound_b - bound_a);
306
307 const double c1 = b2 + 2 * (1 - sub);
308 const double c2 = c1 + saw_slope * (bound_c - bound_b);
309
310 const double d1 = c2 + 2.0 * (shape * (1 - sub) + sub);
311 const double d2 = d1 + saw_slope * (bound_d - bound_c);
312
313 /* dc without sync */
314 const double dc_base = (a1 + a2) / 2 * bound_a
315 + (b1 + b2) / 2 * (bound_b - bound_a)
316 + (c1 + c2) / 2 * (bound_c - bound_b)
317 + (d1 + d2) / 2 * (bound_d - bound_c);
318
319 /* dc offset introduced by sync */
320 const double sync_phase = sync_factor - int (sync_factor);
321
322 double a_avg = (a1 + a2) / 2;
323 double b_avg = (b1 + b2) / 2;
324 double c_avg = (c1 + c2) / 2;
325 double d_avg = (d1 + d2) / 2;
326
327 if (sync_phase < bound_a)
328 {
329 const double frac = (bound_a - sync_phase) / bound_a;
330 const double sync_a2 = a1 * frac + a2 * (1 - frac);
331
332 a_avg = (1 - frac) * (a1 + sync_a2) / 2;
333 b_avg = c_avg = d_avg = 0;
334 }
335 else if (sync_phase < bound_b)
336 {
337 const double frac = (bound_b - sync_phase) / (bound_b - bound_a);
338 const double sync_b2 = b1 * frac + b2 * (1 - frac);
339
340 b_avg = (1 - frac) * (b1 + sync_b2) / 2;
341 c_avg = d_avg = 0;
342 }
343 else if (sync_phase < bound_c)
344 {
345 const double frac = (bound_c - sync_phase) / (bound_c - bound_b);
346 const double sync_c2 = c1 * frac + c2 * (1 - frac);
347
348 c_avg = (1 - frac) * (c1 + sync_c2) / 2;
349 d_avg = 0;
350 }
351 else
352 {
353 const double frac = (bound_d - sync_phase) / (bound_d - bound_c);
354 const double sync_d2 = d1 * frac + d2 * (1 - frac);
355
356 d_avg = (1 - frac) * (d1 + sync_d2) / 2;
357 }
358
359 /* dc sync part of the signal */
360 double dc_sync = a_avg * bound_a
361 + b_avg * (bound_b - bound_a)
362 + c_avg * (bound_c - bound_b)
363 + d_avg * (bound_d - bound_c);
364
365 const double dc = (dc_base * (int) sync_factor + dc_sync) / sync_factor;
366
367 for (auto& voice : unison_voices)
368 {
369 double dest_phase = voice.master_phase;
370
371 double last_value; /* leaky integrator state */
372
373 dest_phase *= sync_factor;
374 dest_phase -= (int) dest_phase;
375
376 voice.slave_phase = dest_phase;
377
378 /* compute voice state and initial value without dc */
379 if (dest_phase < bound_a)
380 {
381 double frac = (bound_a - dest_phase) / bound_a;
382 last_value = a1 * frac + a2 * (1 - frac);
383
384 voice.state = State::A;
385 }
386 else if (dest_phase < bound_b)
387 {
388 double frac = (bound_b - dest_phase) / (bound_b - bound_a);
389 last_value = b1 * frac + b2 * (1 - frac);
390
391 voice.state = State::B;
392 }
393 else if (dest_phase < bound_c)
394 {
395 double frac = (bound_c - dest_phase) / (bound_c - bound_b);
396 last_value = c1 * frac + c2 * (1 - frac);
397
398 voice.state = State::C;
399 }
400 else
401 {
402 double frac = (bound_d - dest_phase) / (bound_d - bound_c);
403 last_value = d1 * frac + d2 * (1 - frac);
404
405 voice.state = State::D;
406 }
407 voice.last_value = last_value - dc;
408 voice.last_dc = dc;
409 voice.current_level = last_value - 1;
410 }
411 }
412 void
413 insert_blep (UnisonVoice& voice, double frac, double weight)
414 {
415 int pos = frac * OVERSAMPLE;
416 const float inter_frac = frac * OVERSAMPLE - pos;
417 const float weight_left = (1 - inter_frac) * weight;
418 const float weight_right = inter_frac * weight;
419
420 pos = std::max (pos, 0);
421 pos = std::min (pos, OVERSAMPLE - 1);
422
423 for (int i = 0; i < WIDTH; i++)
424 {
425 voice.future[i + voice.future_pos] += blep_table[pos] * weight_left + blep_table[pos + 1] * weight_right;
426
427 pos += OVERSAMPLE;
428 }
429 }
430 void
431 insert_future_delta (UnisonVoice& voice, double weight)
432 {
433 voice.future[voice.future_pos + WSHIFT] += weight;
434 }
435
436 double
437 clamp (double d, double min, double max)
438 {
439 return ASE_CLAMP (d, min, max);
440 }
441
442 /* check if slave oscillator has reached target_phase
443 *
444 * when master oscillator sync occurs, only return true if this point in time is
445 * before master oscillator sync
446 */
447 bool
448 check_slave_before_master (UnisonVoice& voice, double target_phase, double sync_factor)
449 {
450 if (voice.slave_phase > target_phase)
451 {
452 if (voice.master_phase > 1)
453 {
454 const double slave_frac = (voice.slave_phase - target_phase) / sync_factor;
455 const double master_frac = voice.master_phase - 1;
456
457 return master_frac < slave_frac;
458 }
459 else
460 return true;
461 }
462 return false;
463 }
464 void
465 process_sample_stereo (float *left_out, float *right_out, unsigned int n_values,
466 const float *freq_in = nullptr,
467 const float *freq_mod_in = nullptr,
468 const float *shape_mod_in = nullptr,
469 const float *sub_mod_in = nullptr,
470 const float *sync_mod_in = nullptr,
471 const float *pulse_mod_in = nullptr,
472 const float *sub_width_mod_in = nullptr)
473 {
474 floatfill (left_out, 0.0, n_values);
475 floatfill (right_out, 0.0, n_values);
476
477 double master_freq = frequency_factor * frequency_base;
478 double pulse_width = clamp (pulse_width_base, 0.01, 0.99);
479 double sub = clamp (sub_base, 0.0, 1.0);
480 double sub_width = clamp (sub_width_base, 0.01, 0.99);
481 double shape = clamp (shape_base, -1.0, 1.0);
482 double sync_factor = fast_exp2 (clamp (sync_base, 0.0, 60.0) / 12);
483
484 /* dc substampling according to control frequency (cpu/quality trade off) */
485 const int dc_steps = max (irintf (rate_ / 4000), 1);
486
487 for (auto& voice : unison_voices)
488 {
489 const double master_freq2inc = 0.5 / rate_ * voice.freq_factor;
490
491 for (unsigned int n = 0; n < n_values; n++)
492 {
493 if (freq_in)
494 master_freq = frequency_factor * fast_voltage2hz (freq_in[n]);
495
496 double master_inc = master_freq * master_freq2inc;
497 if (freq_mod_in)
498 master_inc *= fast_exp2 (freq_mod_in[n] * freq_mod_octaves);
499
500 if (shape_mod_in)
501 shape = clamp (shape_base + shape_mod * shape_mod_in[n], -1.0, 1.0);
502
503 if (sub_mod_in)
504 sub = clamp (sub_base + sub_mod * sub_mod_in[n], 0.0, 1.0);
505
506 if (sync_mod_in)
507 sync_factor = fast_exp2 (clamp (sync_base + sync_mod * sync_mod_in[n], 0.0, 60.0) / 12);
508
509 if (pulse_mod_in)
510 pulse_width = clamp (pulse_width_base + pulse_width_mod * pulse_mod_in[n], 0.01, 0.99);
511
512 if (sub_width_mod_in)
513 sub_width = clamp (sub_width_base + sub_width_mod * sub_width_mod_in[n], 0.01, 0.99);
514
515 /* reset needs parameters, so we need to do it here */
516 if (need_reset_voice_state)
517 {
518 reset_voice_state (shape, pulse_width, sub, sub_width, sync_factor);
519 need_reset_voice_state = false;
520 }
521
522 const double slave_inc = master_inc * sync_factor;
523 const double saw_delta = -4.0 * slave_inc * (shape + 1) * (1 - sub);
524
525 voice.master_phase += master_inc;
526 voice.slave_phase += slave_inc;
527
528 bool state_changed;
529 do
530 {
531 state_changed = false;
532
533 if (voice.state == State::A)
534 {
535 const double bound_a = sub_width * pulse_width;
536
537 if (check_slave_before_master (voice, bound_a, sync_factor))
538 {
539 const double slave_frac = (voice.slave_phase - bound_a) / slave_inc;
540
541 const double jump_a = 2.0 * (shape * (1 - sub) - sub);
542 const double saw = -4.0 * (shape + 1) * (1 - sub) * bound_a;
543 const double blep_height = jump_a + saw - (voice.current_level + (1 - slave_frac) * saw_delta);
544
545 insert_blep (voice, slave_frac, blep_height);
546 voice.current_level += blep_height;
547 voice.state = State::B;
548 state_changed = true;
549 }
550 }
551 if (voice.state == State::B)
552 {
553 const double bound_b = 2 * sub_width * pulse_width + 1 - sub_width - pulse_width;
554
555 if (check_slave_before_master (voice, bound_b, sync_factor))
556 {
557 const double slave_frac = (voice.slave_phase - bound_b) / slave_inc;
558
559 const double jump_ab = 2.0 * ((shape + 1) * (1 - sub) - sub);
560 const double saw = -4.0 * (shape + 1) * (1 - sub) * bound_b;
561 const double blep_height = jump_ab + saw - (voice.current_level + (1 - slave_frac) * saw_delta);
562
563 insert_blep (voice, slave_frac, blep_height);
564 voice.current_level += blep_height;
565 voice.state = State::C;
566 state_changed = true;
567 }
568 }
569 if (voice.state == State::C)
570 {
571 const double bound_c = sub_width * pulse_width + (1 - sub_width);
572
573 if (check_slave_before_master (voice, bound_c, sync_factor))
574 {
575 const double slave_frac = (voice.slave_phase - bound_c) / slave_inc;
576
577 const double jump_abc = 2.0 * (2 * shape + 1) * (1 - sub);
578 const double saw = -4.0 * (shape + 1) * (1 - sub) * bound_c;
579 const double blep_height = jump_abc + saw - (voice.current_level + (1 - slave_frac) * saw_delta);
580
581 insert_blep (voice, slave_frac, blep_height);
582 voice.current_level += blep_height;
583 voice.state = State::D;
584 state_changed = true;
585 }
586 }
587 if (voice.state == State::D)
588 {
589 if (check_slave_before_master (voice, 1, sync_factor))
590 {
591 voice.slave_phase -= 1;
592
593 const double slave_frac = voice.slave_phase / slave_inc;
594
595 voice.current_level += (1 - slave_frac) * saw_delta;
596
597 insert_blep (voice, slave_frac, -voice.current_level);
598
599 voice.current_level = saw_delta * slave_frac - saw_delta;
600 voice.state = State::A;
601 state_changed = true;
602 }
603 }
604 if (!state_changed && voice.master_phase > 1)
605 {
606 voice.master_phase -= 1;
607
608 const double master_frac = voice.master_phase / master_inc;
609
610 const double new_slave_phase = voice.master_phase * sync_factor;
611
612 voice.current_level += (1 - master_frac) * saw_delta;
613
614 insert_blep (voice, master_frac, -voice.current_level);
615
616 voice.current_level = saw_delta * master_frac - saw_delta;
617 voice.slave_phase = new_slave_phase;
618
619 voice.state = State::A;
620 state_changed = true;
621 }
622 }
623 while (state_changed); // rerun all state checks if state was modified
624
625 if (voice.dc_steps > 0)
626 {
627 voice.dc_steps--;
628 }
629 else
630 {
631 const double dc = estimate_dc (shape, pulse_width, sub, sub_width, sync_factor);
632
633 voice.dc_steps = dc_steps - 1;
634 voice.dc_delta = (voice.last_dc - dc) / dc_steps;
635 voice.last_dc = dc;
636 }
637
638 voice.current_level += saw_delta;
639 insert_future_delta (voice, saw_delta + voice.dc_delta); // align with the impulses
640
641 /* leaky integration */
642 double value = leaky_a * voice.last_value + voice.pop_future();
643 voice.last_value = value;
644
645 left_out[n] += value * voice.left_factor;
646 right_out[n] += value * voice.right_factor;
647 }
648 }
649 }
650};
651
652class Osc /* simple interface to OscImpl */
653{
654 double
655 to_sync (double sync_factor)
656 {
657 return 12 * log (sync_factor) / log (2);
658 }
659public:
660 OscImpl osc_impl;
661
662 double rate;
663 double freq;
664 double pulse_width = 0.5;
665 double shape = 0; // 0 = saw, range [-1:1]
666 double sub = 0;
667 double sub_width = 0.5;
668
669 double master_freq;
670
671 void
672 set_unison (size_t n_voices, float detune, float stereo)
673 {
674 osc_impl.set_unison (n_voices, detune, stereo);
675 }
676 double
677 test_seek_to (double phase)
678 {
679 ASE_ASSERT_RETURN (osc_impl.unison_voices.size() > 0, 0);
680 osc_impl.reset();
681 osc_impl.reset_master (osc_impl.unison_voices[0], phase); // jump to phase
682
683 process_sample(); // propagate parameters & perform reset
684 return osc_impl.unison_voices[0].last_value;
685 }
686 double
687 process_sample()
688 {
689 float left, right;
690
691 process_sample_stereo (&left, &right);
692
693 return (left + right) / 2;
694 }
695 void
696 process_sample_stereo (float *left_out, float *right_out)
697 {
698 osc_impl.set_rate (rate);
699 osc_impl.sync_base = to_sync (freq / master_freq);
700 osc_impl.pulse_width_base = pulse_width;
701 osc_impl.shape_base = shape;
702 osc_impl.sub_base = sub;
703 osc_impl.sub_width_base = sub_width;
704 osc_impl.frequency_base = master_freq;
705
706 return osc_impl.process_sample_stereo (left_out, right_out, 1);
707 }
708};
709
710} // BlepUtils
711} // Ase
712
713#endif // __ASE_DEVICES_BLEPOSC_HH__
#define ASE_CLAMP(v, mi, ma)
Yield v clamped to [ mi .. ma ].
Definition cxxaux.hh:50
#define ASE_ASSERT_RETURN(expr,...)
Return from the current function if expr evaluates to false and issue an assertion warning.
Definition cxxaux.hh:82
typedef int
log
typedef float
T max(T... args)
T min(T... args)
The Anklang C++ API namespace.
Definition api.hh:9
void floatfill(float *dst, float f, size_t n)
Fill n values of dst with f.
Definition datautils.hh:29
float fast_voltage2hz(float x)
Float precision variant of voltage2hz using fast_exp2().
Definition signalmath.hh:80
float fast_exp2(float x)
Definition mathutils.hh:82
int irintf(float f)
Definition mathutils.hh:14
double random_frange(double begin, double end)
pow
sqrt