From a09eca2f83cb664637855b444315d3394868096c Mon Sep 17 00:00:00 2001
From: Krzysztof Foltman <wdev@foltman.com>
Date: Sun, 2 Jun 2013 10:40:57 +0100
Subject: [PATCH] Sampler: new filter types, fix 6dB/oct lowpass/highpass (had
 broken pre-warping).

The whole math around 6dB/oct filters is probably off by some factor.
---
 biquad-float.h  | 45 ++++++++++++++++++++++++++++++---------------
 dspmath.h       |  1 +
 sampler.c       |  2 ++
 sampler_layer.h | 12 ++++++++++--
 sampler_voice.c |  8 ++++++--
 5 files changed, 49 insertions(+), 19 deletions(-)

diff --git a/biquad-float.h b/biquad-float.h
index 4ffaa14..9959f8f 100644
--- a/biquad-float.h
+++ b/biquad-float.h
@@ -149,6 +149,7 @@ static inline void cbox_biquadf_set_peakeq_rbj(struct cbox_biquadf_coeffs *coeff
     coeffs->b2 = ib0 * (1 - alpha/A);
 }
 
+// This is my math, and it's rather suspect
 static inline void cbox_biquadf_set_1plp(struct cbox_biquadf_coeffs *coeffs, float freq, float sr)
 {
     float w = hz2w(freq, sr);
@@ -167,7 +168,7 @@ static inline void cbox_biquadf_set_1plp(struct cbox_biquadf_coeffs *coeffs, flo
 static inline void cbox_biquadf_set_1php(struct cbox_biquadf_coeffs *coeffs, float freq, float sr)
 {
     float w = hz2w(freq, sr);
-    float x = tan (w * 0.5f); // XXXKF needs to be checked for accuracy
+    float x = tan (w * 0.5f); 
     float q = 1 / (1 + x);
     float a01 = x*q;
     float b1 = a01 - q;
@@ -179,32 +180,46 @@ static inline void cbox_biquadf_set_1php(struct cbox_biquadf_coeffs *coeffs, flo
     coeffs->b2 = 0;
 }
 
-static inline void cbox_biquadf_set_1plp_lookup(struct cbox_biquadf_coeffs *coeffs, struct cbox_sincos *sincos)
+static inline void cbox_biquadf_set_1p(struct cbox_biquadf_coeffs *coeffs, float a0, float a1, float b1, int two_copies)
 {
-    float x = sincos->sine / sincos->cosine;
+    if (two_copies)
+    {
+        // (a0 + a1z) * (a0 + a1z) = a0^2 + 2*a0*a1*z + a1^2*z^2
+        // (1 - b1z) * (1 - b1z) = 1 - 2b1*z + b1^2*z^2
+        coeffs->a0 = a0*a0;
+        coeffs->a1 = 2*a0*a1;
+        coeffs->b1 = 2 * b1;
+        coeffs->a2 = a1*a1;
+        coeffs->b2 = b1*b1;
+    }
+    else
+    {
+        coeffs->a0 = a0;
+        coeffs->a1 = a1;
+        coeffs->b1 = b1;
+        coeffs->a2 = 0;
+        coeffs->b2 = 0;
+    }
+}
+
+static inline void cbox_biquadf_set_1plp_lookup(struct cbox_biquadf_coeffs *coeffs, struct cbox_sincos *sincos, int two_copies)
+{
+    float x = sincos->prewarp;
     float q = 1 / (1 + x);
     float a01 = x*q;
     float b1 = a01 - q;
     
-    coeffs->a0 = a01;
-    coeffs->a1 = a01;
-    coeffs->b1 = b1;
-    coeffs->a2 = 0;
-    coeffs->b2 = 0;
+    cbox_biquadf_set_1p(coeffs, a01, a01, b1, two_copies);    
 }
 
-static inline void cbox_biquadf_set_1php_lookup(struct cbox_biquadf_coeffs *coeffs, struct cbox_sincos *sincos)
+static inline void cbox_biquadf_set_1php_lookup(struct cbox_biquadf_coeffs *coeffs, struct cbox_sincos *sincos, int two_copies)
 {
-    float x = sincos->sine / sincos->cosine;
+    float x = sincos->prewarp;
     float q = 1 / (1 + x);
     float a01 = x*q;
     float b1 = a01 - q;
     
-    coeffs->a0 = q;
-    coeffs->a1 = -q;
-    coeffs->b1 = b1;
-    coeffs->a2 = 0;
-    coeffs->b2 = 0;
+    cbox_biquadf_set_1p(coeffs, a01, -a01, b1, two_copies);    
 }
 
 #if USE_NEON
diff --git a/dspmath.h b/dspmath.h
index fb9c122..123b952 100644
--- a/dspmath.h
+++ b/dspmath.h
@@ -36,6 +36,7 @@ struct cbox_sincos
 {
     float sine;
     float cosine;
+    float prewarp;
 };
 
 static inline float hz2w(float hz, float sr)
diff --git a/sampler.c b/sampler.c
index 559cba3..be95056 100644
--- a/sampler.c
+++ b/sampler.c
@@ -593,6 +593,8 @@ MODULE_CREATE_FUNCTION(sampler)
         float omega=(float)(2*M_PI*freq/srate);
         m->sincos[i].sine = sinf(omega);
         m->sincos[i].cosine = cosf(omega);
+        m->sincos[i].prewarp = 2.0 * tan(hz2w(freq, srate) * 0.5f);
+
     }
     for (i = 0; ; i++)
     {
diff --git a/sampler_layer.h b/sampler_layer.h
index f37679c..5adf60d 100644
--- a/sampler_layer.h
+++ b/sampler_layer.h
@@ -91,6 +91,10 @@ enum sampler_filter_type
     sft_bp12,
     sft_lp6,
     sft_hp6,
+    sft_lp12nr,
+    sft_hp12nr,
+    sft_lp24nr,
+    sft_hp24nr,
 };
 
 #define ENUM_VALUES_sampler_filter_type(MACRO) \
@@ -101,7 +105,11 @@ enum sampler_filter_type
     MACRO("hpf_4p", sft_hp24) \
     MACRO("bpf_4p", sft_bp12) \
     MACRO("lpf_1p", sft_lp6)  \
-    MACRO("hpf_1p", sft_hp6)
+    MACRO("hpf_1p", sft_hp6)  \
+    MACRO("lpf_2p_nores", sft_lp12nr)  \
+    MACRO("hpf_2p_nores", sft_hp12nr)  \
+    MACRO("lpf_4p_nores", sft_lp24nr)  \
+    MACRO("hpf_4p_nores", sft_hp24nr)  \
 
 #define ENUM_LIST(MACRO) \
     MACRO(sampler_loop_mode) \
@@ -392,7 +400,7 @@ static inline gboolean sampler_layer_data_is_4pole(struct sampler_layer_data *v)
 {
     if (v->cutoff == -1)
         return FALSE;
-    return v->fil_type == sft_lp24 || v->fil_type == sft_hp24 || v->fil_type == sft_bp12;
+    return v->fil_type == sft_lp24 || v->fil_type == sft_lp24nr || v->fil_type == sft_hp24 || v->fil_type == sft_hp24nr || v->fil_type == sft_bp12;
 }
 
 #endif
diff --git a/sampler_voice.c b/sampler_voice.c
index f4c2844..6228788 100644
--- a/sampler_voice.c
+++ b/sampler_voice.c
@@ -576,10 +576,14 @@ void sampler_voice_process(struct sampler_voice *v, struct sampler_module *m, cb
             cbox_biquadf_set_bp_rbj_lookup(&v->filter_coeffs, &m->sincos[(int)logcutoff], resonance);
             break;
         case sft_lp6:
-            cbox_biquadf_set_1plp_lookup(&v->filter_coeffs, &m->sincos[(int)logcutoff]);
+        case sft_lp12nr:
+        case sft_lp24nr:
+            cbox_biquadf_set_1plp_lookup(&v->filter_coeffs, &m->sincos[(int)logcutoff], l->fil_type != sft_lp6);
             break;
         case sft_hp6:
-            cbox_biquadf_set_1php_lookup(&v->filter_coeffs, &m->sincos[(int)logcutoff]);
+        case sft_hp12nr:
+        case sft_hp24nr:
+            cbox_biquadf_set_1php_lookup(&v->filter_coeffs, &m->sincos[(int)logcutoff], l->fil_type != sft_hp6);
             break;
         default:
             assert(0);
-- 
2.11.4.GIT