Implement overdrive in resonant filter

This patch implements nonlinear distortion in the resonant ladder
filter. It's based on the differential equations in the Huovilainen '04
DAFx paper, but using matrix exponential to compute the evolution of the
state variables. The implementation is scalar (and calls into sqrtf),
but designed to be implemented in very efficient SIMD.

android/res/layout/piano2.xml

@@ -19,7 +19,9 @@
         android:id="@+id/resonanceLabel"
         android:gravity="center_horizontal" />
     <TextView
-        android:text=""
+        android:text="@+string/overdrive"
+        android:width="30sp"
+        android:id="@+id/overdriveLabel"
         android:gravity="center_horizontal" />
     <TextView
         android:text=""
@@ -44,9 +46,15 @@
         app:min="0"
         app:max="1"
         android:layout_margin="2dp" />
+    <com.levien.synthesizer.android.widgets.knob.KnobView
+        android:id="@+id/overdriveKnob"
+        app:value="0.0"
+        app:min="0"
+        app:max="1"
+        android:layout_margin="2dp" />
     <LinearLayout
         android:orientation="vertical"
-        android:layout_span="4" >
+        android:layout_span="3" >
     <Spinner
         android:id="@+id/presetSpinner"
         android:layout_width="wrap_content"

android/res/values/strings.xml

@@ -48,6 +48,7 @@
   <!-- Low Pass Filter -->
   <string name="cutoff">Cutoff</string>
   <string name="resonance">Resonance</string>
+  <string name="overdrive">Overdrive</string>
 
   <!-- Amp -->
   <string name="attack">Attack</string>
@@ -79,4 +80,7 @@
   </string-array>
 
   <string name="capture">Capture</string>
+
+  <string name="settings">Settings</string>
+
 </resources>

android/src/com/levien/synthesizer/android/ui/PianoActivity2.java

@@ -68,6 +68,7 @@ public class PianoActivity2 extends Activity {
     piano_ = (PianoView)findViewById(R.id.piano);
     cutoffKnob_ = (KnobView)findViewById(R.id.cutoffKnob);
     resonanceKnob_ = (KnobView)findViewById(R.id.resonanceKnob);
+    overdriveKnob_ = (KnobView)findViewById(R.id.overdriveKnob);
     presetSpinner_ = (Spinner)findViewById(R.id.presetSpinner);
 
     AudioParams params = new AudioParams(44100, 64);
@@ -116,6 +117,12 @@ public class PianoActivity2 extends Activity {
         androidGlue_.onController(0, 2, value);
       }
     });
+    overdriveKnob_.setKnobListener(new KnobListener() {
+      public void onKnobChanged(double newValue) {
+        int value = (int)Math.round(newValue * 127);
+        androidGlue_.onController(0, 3, value);
+      }
+    });
     
     piano_.bindTo(androidGlue_);
     if (Build.VERSION.SDK_INT >= Build.VERSION_CODES.ICE_CREAM_SANDWICH) {
@@ -289,6 +296,7 @@ public class PianoActivity2 extends Activity {
   private PianoView piano_;
   private KnobView cutoffKnob_;
   private KnobView resonanceKnob_;
+  private KnobView overdriveKnob_;
   private Spinner presetSpinner_;
   private Handler statusHandler_;
   private Runnable statusRunnable_;

cpp/src/resofilter.cc

@@ -23,10 +23,12 @@
 
 #include <string.h>
 #include <stdio.h>
+#include <math.h>
 
 #include "synth.h"
 #include "freqlut.h"
 #include "exp2.h"
+#include "aligned_buf.h"
 #include "resofilter.h"
 
 double this_sample_rate;
@@ -173,21 +175,51 @@ void test_matrix() {
 }
 
 #if defined(USE_MATRIX)
+static float sigmoid(float x, float overdrive) {
+  float xs = overdrive * x * (1.0 / (1 << 24));
+  float isq = 1.0 / sqrtf(1 + xs * xs);
+  return x * isq;
+}
+
 void ResoFilter::process(const int32_t **inbufs, const int32_t *control_in,
                          const int32_t *control_last, int32_t **outbufs) {
-  float a[20];
-  make_state_transition(a, compute_alpha(control_in[0]), control_in[1]);
+  AlignedBuf<float, 20> a;
+  make_state_transition(a.get(), compute_alpha(control_in[0]), control_in[1]);
+  float overdrive = control_in[2] * (1.0 / (1 << 24));
   const int32_t *ibuf = inbufs[0];
   int32_t *obuf = outbufs[0];
+  if (overdrive == 0) {
     for (int i = 0; i < n; i++) {
       float signal = ibuf[i];
       float tmp[4];
-    matvec4(tmp, a + 4, x);
+      matvec4(tmp, a.get() + 4, x);
       for (int k = 0; k < 4; k++) {
-      x[k] = tmp[k] + signal * a[k];
+        x[k] = tmp[k] + signal * a.get()[k];
         obuf[i] = x[3];
       }
     }
+  } else {
+    float ogain = 1 + overdrive;
+    float k = control_in[1] * (1.0 / (1<<24));
+    for (int i = 0; i < 4; i++) {
+      a.get()[4 + 5 * i] -= 1.0;
+      a.get()[16 + i] += k * a.get()[i];
+    }
+    for (int i = 0; i < n; i++) {
+      float signal = ibuf[i];
+      float tmp[4];
+      float tx[4];
+      for (int j = 0; j < 4; j++) {
+        tx[j] = sigmoid(x[j], overdrive);
+      }
+      matvec4(tmp, a.get() + 4, tx);
+      float xin = sigmoid(signal - k * x[3], overdrive);
+      for (int j = 0; j < 4; j++) {
+        x[j] += tmp[j] + xin * a.get()[j];
+        obuf[i] = x[3] * ogain;
+      }
+    }
+  }
 }
 #else
 void ResoFilter::process(const int32_t **inbufs, const int32_t *control_in,

cpp/src/synth_unit.cc

@@ -49,6 +49,7 @@ char epiano[] = {
 void SynthUnit::Init(double sample_rate) {
   Freqlut::init(sample_rate);
   Exp2::init();
+  Tanh::init();
   Sin::init();
   Lfo::init(sample_rate);
   PitchEnv::init(sample_rate);
@@ -68,6 +69,7 @@ SynthUnit::SynthUnit(RingBuffer *ring_buffer) {
   current_note_ = 0;
   filter_control_[0] = 258847126;
   filter_control_[1] = 0;
+  filter_control_[2] = 0;
   controllers_.values_[kControllerPitch] = 0x2000;
   sustain_ = false;
   extra_buf_size_ = 0;
@@ -121,6 +123,7 @@ void SynthUnit::SetController(int controller, int value) {
 int SynthUnit::ProcessMidiMessage(const uint8_t *buf, int buf_size) {
   uint8_t cmd = buf[0];
   uint8_t cmd_type = cmd & 0xf0;
+  //LOGI("got %d midi: %02x %02x %02x", buf_size, buf[0], buf[1], buf[2]);
   if (cmd_type == 0x80 || (cmd_type == 0x90 && buf[2] == 0)) {
     if (buf_size >= 3) {
       // note off
@@ -163,6 +166,8 @@ int SynthUnit::ProcessMidiMessage(const uint8_t *buf, int buf_size) {
         filter_control_[0] = 142365917 + value * 917175;
       } else if (controller == 2) {
         filter_control_[1] = value * 528416;
+      } else if (controller == 3) {
+        filter_control_[2] = value * 528416;
       } else if (controller == 64) {
         sustain_ = value != 0;
         if (!sustain_) {

cpp/src/synth_unit.h

@@ -70,7 +70,7 @@ class SynthUnit {
   Controllers controllers_;
 
   ResoFilter filter_;
-  int32_t filter_control_[2];
+  int32_t filter_control_[3];
   bool sustain_;
 
   // Extra buffering for when GetSamples wants a buffer not a multiple of N