From 52233a9d63ff74e91060451063a8bd78d1bd94c8 Mon Sep 17 00:00:00 2001
From: Stefano D'Angelo <stefano.dangelo@orastron.com>
Date: Mon, 14 Aug 2023 17:34:22 +0200
Subject: [PATCH] improved bw_drive algorithm + added notes to
 bw_osc_{pulse,tri}

---
 TODO                   |  3 +--
 include/bw_drive.h     | 35 ++++++++++++++++++++---------------
 include/bw_osc_pulse.h |  4 ++--
 include/bw_osc_tri.h   |  4 ++--
 4 files changed, 25 insertions(+), 21 deletions(-)

diff --git a/TODO b/TODO
index a14cdb9..9c14fed 100644
--- a/TODO
+++ b/TODO
@@ -7,12 +7,11 @@ code:
 * osc post filter (and one pole init, slew rate, etc.) val from input? set state instead?
 * audio rate optional pulse width/slope inputs?
 * empty functions etc. to keep consistency and forward compatibility?
-* should clip slope in triangle?
 * common smoothing policy (as control rate as possible?) - smoothing control?
 * svf bandpass out polarity too confusing (inverted in mm2)?
 * get_y_z1, common strategy?
 * bw_satur gain compensation to divide by actual gain (derivative) rather than gain parameter?
-* add initial state (x0) to reset state of lp1, ap1, mm1, hs1, ls1, others? all?
+* add initial state (x0) to reset state of all modules?
 * bw_comb: should also modulate feedback?
 * bw_comb: integer target delay values?
 * bw_svf/lp1 automatically limited (max) prewarp frequency to avoid instability?
diff --git a/include/bw_drive.h b/include/bw_drive.h
index 8fbfe89..85737ca 100644
--- a/include/bw_drive.h
+++ b/include/bw_drive.h
@@ -22,7 +22,8 @@
  *  module_type {{{ dsp }}}
  *  version {{{ 1.0.0 }}}
  *  requires {{{
- *    bw_common bw_gain bw_math bw_mm2 bw_one_pole bw_peak bw_satur bw_svf
+ *    bw_common bw_gain bw_lp1 bw_math bw_mm2 bw_one_pole bw_peak bw_satur
+ *    bw_svf
  *  }}}
  *  description {{{
  *    Overdrive effect.
@@ -33,6 +34,8 @@
  *    <ul>
  *      <li>Version <strong>1.0.0</strong>:
  *        <ul>
+ *          <li>Improved algorithm to be a bit more faithful to the
+ *              original.</li>
  *          <li><code>bw_drive_process()</code> and
  *              <code>bw_drive_process_multi()</code> now use
  *              <code>size_t</code> to count samples and channels.</li>
@@ -179,6 +182,7 @@ static inline void bw_drive_set_volume(bw_drive_coeffs *BW_RESTRICT coeffs, floa
 #include <bw_svf.h>
 #include <bw_peak.h>
 #include <bw_satur.h>
+#include <bw_lp1.h>
 #include <bw_gain.h>
 
 #ifdef __cplusplus
@@ -190,7 +194,7 @@ struct bw_drive_coeffs {
 	bw_svf_coeffs	hp2_coeffs;
 	bw_peak_coeffs	peak_coeffs;
 	bw_satur_coeffs	satur_coeffs;
-	bw_svf_coeffs	lp2_coeffs;
+	bw_lp1_coeffs	lp1_coeffs;
 	bw_gain_coeffs	gain_coeffs;
 };
 
@@ -199,27 +203,28 @@ struct bw_drive_state {
 	bw_svf_state	hp2_state;
 	bw_peak_state	peak_state;
 	bw_satur_state	satur_state;
-	bw_svf_state	lp2_state;
+	bw_lp1_state	lp1_state;
 };
 
 static inline void bw_drive_init(bw_drive_coeffs *BW_RESTRICT coeffs) {
 	bw_svf_init(&coeffs->hp2_coeffs);
 	bw_peak_init(&coeffs->peak_coeffs);
 	bw_satur_init(&coeffs->satur_coeffs);
-	bw_svf_init(&coeffs->lp2_coeffs);
+	bw_lp1_init(&coeffs->lp1_coeffs);
 	bw_gain_init(&coeffs->gain_coeffs);
 	bw_svf_set_cutoff(&coeffs->hp2_coeffs, 16.f);
+	bw_peak_set_peak_gain_dB(&coeffs->peak_coeffs, 20.f);
 	bw_peak_set_cutoff(&coeffs->peak_coeffs, 2e3f);
 	bw_peak_set_bandwidth(&coeffs->peak_coeffs, 10.f);
 	bw_satur_set_gain(&coeffs->satur_coeffs, 1.5f);
-	bw_svf_set_cutoff(&coeffs->lp2_coeffs, 600.f + (7.5e3f - 600.f) * 0.125f);
+	bw_lp1_set_cutoff(&coeffs->lp1_coeffs, 400.f + (5e3f - 400.f) * 0.125f);
 }
 
 static inline void bw_drive_set_sample_rate(bw_drive_coeffs *BW_RESTRICT coeffs, float sample_rate) {
 	bw_svf_set_sample_rate(&coeffs->hp2_coeffs, sample_rate);
 	bw_peak_set_sample_rate(&coeffs->peak_coeffs, sample_rate);
 	bw_satur_set_sample_rate(&coeffs->satur_coeffs, sample_rate);
-	bw_svf_set_sample_rate(&coeffs->lp2_coeffs, sample_rate);
+	bw_lp1_set_sample_rate(&coeffs->lp1_coeffs, sample_rate);
 	bw_gain_set_sample_rate(&coeffs->gain_coeffs, sample_rate);
 	bw_svf_reset_coeffs(&coeffs->hp2_coeffs);
 	bw_satur_reset_coeffs(&coeffs->satur_coeffs);
@@ -227,7 +232,7 @@ static inline void bw_drive_set_sample_rate(bw_drive_coeffs *BW_RESTRICT coeffs,
 
 static inline void bw_drive_reset_coeffs(bw_drive_coeffs *BW_RESTRICT coeffs) {
 	bw_peak_reset_coeffs(&coeffs->peak_coeffs);
-	bw_svf_reset_coeffs(&coeffs->lp2_coeffs);
+	bw_lp1_reset_coeffs(&coeffs->lp1_coeffs);
 	bw_gain_reset_coeffs(&coeffs->gain_coeffs);
 }
 
@@ -235,18 +240,18 @@ static inline void bw_drive_reset_state(const bw_drive_coeffs *BW_RESTRICT coeff
 	bw_svf_reset_state(&coeffs->hp2_coeffs, &state->hp2_state, 0.f);
 	bw_peak_reset_state(&coeffs->peak_coeffs, &state->peak_state, 0.f);
 	bw_satur_reset_state(&coeffs->satur_coeffs, &state->satur_state);
-	bw_svf_reset_state(&coeffs->lp2_coeffs, &state->lp2_state, 0.f);
+	bw_lp1_reset_state(&coeffs->lp1_coeffs, &state->lp1_state, 0.f);
 }
 
 static inline void bw_drive_update_coeffs_ctrl(bw_drive_coeffs *BW_RESTRICT coeffs) {
 	bw_peak_update_coeffs_ctrl(&coeffs->peak_coeffs);
-	bw_svf_update_coeffs_ctrl(&coeffs->lp2_coeffs);
+	bw_lp1_update_coeffs_ctrl(&coeffs->lp1_coeffs);
 	bw_gain_update_coeffs_ctrl(&coeffs->gain_coeffs);
 }
 
 static inline void bw_drive_update_coeffs_audio(bw_drive_coeffs *BW_RESTRICT coeffs) {
 	bw_peak_update_coeffs_audio(&coeffs->peak_coeffs);
-	bw_svf_update_coeffs_audio(&coeffs->lp2_coeffs);
+	bw_lp1_update_coeffs_audio(&coeffs->lp1_coeffs);
 	bw_gain_update_coeffs_audio(&coeffs->gain_coeffs);
 }
 
@@ -254,9 +259,9 @@ static inline float bw_drive_process1(const bw_drive_coeffs *BW_RESTRICT coeffs,
 	float v_lp, v_hp, v_bp;
 	bw_svf_process1(&coeffs->hp2_coeffs, &state->hp2_state, x, &v_lp, &v_bp, &v_hp);
 	float y = bw_peak_process1(&coeffs->peak_coeffs, &state->peak_state, v_hp);
-	y = v_hp + 0.5f * bw_satur_process1_comp(&coeffs->satur_coeffs, &state->satur_state, y);
-	bw_svf_process1(&coeffs->lp2_coeffs, &state->lp2_state, y, &v_lp, &v_bp, &v_hp);
-	return bw_gain_process1(&coeffs->gain_coeffs, v_lp);
+	y = v_hp + bw_satur_process1_comp(&coeffs->satur_coeffs, &state->satur_state, y - v_hp);
+	y = bw_lp1_process1(&coeffs->lp1_coeffs, &state->lp1_state, y);
+	return bw_gain_process1(&coeffs->gain_coeffs, y);
 }
 
 static inline void bw_drive_process(bw_drive_coeffs *BW_RESTRICT coeffs, bw_drive_state *BW_RESTRICT state, const float *x, float *y, size_t n_samples) {
@@ -277,11 +282,11 @@ static inline void bw_drive_process_multi(bw_drive_coeffs *BW_RESTRICT coeffs, b
 }
 
 static inline void bw_drive_set_drive(bw_drive_coeffs *BW_RESTRICT coeffs, float value) {
-	bw_peak_set_peak_gain_dB(&coeffs->peak_coeffs, 40.f * value);
+	bw_peak_set_peak_gain_dB(&coeffs->peak_coeffs, 20.f + 20.f * value);
 }
 
 static inline void bw_drive_set_tone(bw_drive_coeffs *BW_RESTRICT coeffs, float value) {
-	bw_svf_set_cutoff(&coeffs->lp2_coeffs, 600.f + (7.5e3f - 600.f) * value * value * value);
+	bw_lp1_set_cutoff(&coeffs->lp1_coeffs, 400.f + (5e3f - 400.f) * value * value * value);
 }
 
 static inline void bw_drive_set_volume(bw_drive_coeffs *BW_RESTRICT coeffs, float value) {
diff --git a/include/bw_osc_pulse.h b/include/bw_osc_pulse.h
index 256326e..a84b193 100644
--- a/include/bw_osc_pulse.h
+++ b/include/bw_osc_pulse.h
@@ -234,7 +234,7 @@ static inline void bw_osc_pulse_update_coeffs_audio(bw_osc_pulse_coeffs *BW_REST
 }
 
 static inline float bw_osc_pulse_process1(const bw_osc_pulse_coeffs *BW_RESTRICT coeffs, float x) {
-	const float pw = bw_one_pole_get_y_z1(&coeffs->smooth_state);
+	const float pw = bw_one_pole_get_y_z1(&coeffs->smooth_state); // output should always be in bounds given good input
 	return bw_signf(pw - x);
 }
 
@@ -246,7 +246,7 @@ static inline float bw_osc_pulse_blep_diff(float x) {
 }
 
 static inline float bw_osc_pulse_process1_antialias(const bw_osc_pulse_coeffs *BW_RESTRICT coeffs, float x, float x_phase_inc) {
-	const float pw = bw_one_pole_get_y_z1(&coeffs->smooth_state);
+	const float pw = bw_one_pole_get_y_z1(&coeffs->smooth_state); // output should always be in bounds given good input
 	const float pw_m_phase = pw - x;
 	float v = bw_copysignf(1.f, pw_m_phase); // pw = phase case should be properly compensated by the AA
 	if (x_phase_inc != 0.f) {
diff --git a/include/bw_osc_tri.h b/include/bw_osc_tri.h
index 2020cfa..e442e1c 100644
--- a/include/bw_osc_tri.h
+++ b/include/bw_osc_tri.h
@@ -236,7 +236,7 @@ static inline void bw_osc_tri_update_coeffs_audio(bw_osc_tri_coeffs *BW_RESTRICT
 }
 
 static inline float bw_osc_tri_process1(const bw_osc_tri_coeffs *BW_RESTRICT coeffs, float x) {
-	const float slope = bw_one_pole_get_y_z1(&coeffs->smooth_state);
+	const float slope = bw_one_pole_get_y_z1(&coeffs->smooth_state); // output should always be in bounds given good input
 	const float phase_d = x + x;
 	return x < slope ? (phase_d - slope) * bw_rcpf(slope) : (1.f + slope - phase_d) * bw_rcpf(1.f - slope);
 }
@@ -249,7 +249,7 @@ static inline float bw_osc_tri_blamp_diff(float x) {
 }
 
 static inline float bw_osc_tri_process1_antialias(const bw_osc_tri_coeffs *BW_RESTRICT coeffs, float x, float x_phase_inc) {
-	const float slope = bw_one_pole_get_y_z1(&coeffs->smooth_state);
+	const float slope = bw_one_pole_get_y_z1(&coeffs->smooth_state); // output should always be in bounds given good input
 	const float s_1_p_pw = 1.f + slope;
 	const float s_1_m_pw = 1.f - slope;
 	const float phase_d = x + x;