Wild ZamTube DSP!

Signed-off-by: Damien Zammit <damien@zamaudio.com>

11 files changed, 298 insertions, 499 deletions

diff --git a/plugins/ZamTube/Makefile b/plugins/ZamTube/Makefile
index 71ae929..d641ac5 100644
--- a/plugins/ZamTube/Makefile
+++ b/plugins/ZamTube/Makefile
@@ -14,7 +14,7 @@ NAME = ZamTube
 
 OBJS_DSP = \
 	ZamTubePlugin.cpp.o \
-	wdf.cpp.o
+	triode.cpp.o
 
 OBJS_UI  = \
 	ZamTubeArtwork.cpp.o \
diff --git a/plugins/ZamTube/ZamTubePlugin.cpp b/plugins/ZamTube/ZamTubePlugin.cpp
index e4aa101..da5bf1d 100644
--- a/plugins/ZamTube/ZamTubePlugin.cpp
+++ b/plugins/ZamTube/ZamTubePlugin.cpp
@@ -19,17 +19,14 @@
  */
 
 #include "ZamTubePlugin.hpp"
+#include <stdio.h>
 
 START_NAMESPACE_DISTRHO
 
 // -----------------------------------------------------------------------
 
 ZamTubePlugin::ZamTubePlugin()
-    : Plugin(paramCount, 1, 0), // 1 program, 0 states
-		Vi(0.0,10000.0), Ci(0.0000001,48000.0), Ck(0.00001,48000.0), Co(0.00000001,48000.0), Ro(1000000.0), Rg(20000.0), 
-		Ri(1000000.0), Rk(1000.0), E(250.0,100000.0), 
-		S0(&Ci,&Vi), I0(&S0), P0(&I0,&Ri), S1(&Rg,&P0), 
-		I1(&S1), I3(&Ck,&Rk), S2(&Co,&Ro), I4(&S2), P2(&I4,&E) 
+    : Plugin(paramCount, 1, 0) // 1 program, 0 states
 {
     // set default values
     loadProgram(0);
@@ -48,7 +45,7 @@ void ZamTubePlugin::initParameter(uint32_t index, Parameter& parameter)
         parameter.symbol     = "tubedrive";
         parameter.unit       = " ";
         parameter.ranges.def = 0.0f;
-        parameter.ranges.min = -30.0f;
+        parameter.ranges.min = 0.0f;
         parameter.ranges.max = 30.0f;
         break;
     case paramBass:
@@ -56,18 +53,18 @@ void ZamTubePlugin::initParameter(uint32_t index, Parameter& parameter)
         parameter.name       = "Bass";
         parameter.symbol     = "bass";
         parameter.unit       = " ";
-        parameter.ranges.def = 0.5f;
+        parameter.ranges.def = 0.0f;
         parameter.ranges.min = 0.0f;
-        parameter.ranges.max = 1.0f;
+        parameter.ranges.max = 0.5f;
         break;
     case paramMiddle:
         parameter.hints      = kParameterIsAutomable;
         parameter.name       = "Mids";
         parameter.symbol     = "mids";
         parameter.unit       = " ";
-        parameter.ranges.def = 0.5f;
+        parameter.ranges.def = 0.0f;
         parameter.ranges.min = 0.0f;
-        parameter.ranges.max = 1.0f;
+        parameter.ranges.max = 0.5f;
         break;
     case paramTreble:
         parameter.hints      = kParameterIsAutomable;
@@ -76,7 +73,7 @@ void ZamTubePlugin::initParameter(uint32_t index, Parameter& parameter)
         parameter.unit       = " ";
         parameter.ranges.def = 0.0f;
         parameter.ranges.min = 0.0f;
-        parameter.ranges.max = 1.0f;
+        parameter.ranges.max = 0.5f;
         break;
     case paramToneStack:
         parameter.hints      = kParameterIsAutomable | kParameterIsInteger;
@@ -184,8 +181,8 @@ void ZamTubePlugin::loadProgram(uint32_t index)
 
     /* Default parameter values */
     tubedrive = 0.0f;
-    bass = 0.5f;
-    middle = 0.5f;
+    bass = 0.0f;
+    middle = 0.0f;
     treble = 0.0f;
     tonestack = 0.0f;
     mastergain = 0.0f;
@@ -206,51 +203,25 @@ void ZamTubePlugin::activate()
 	T Fs = getSampleRate();
 	
 	// Passive components
-	T ci = 0.0000001;	//100nF
-	T ck = 0.00001;		//10uF
-	T co = 0.00000001;	//10nF
-	T ro = 1000000.0;	//1Mohm
-	T rp = 100000.0;	//100kohm
-	T rg = 20000.0;		//20kohm
-	T ri = 1000000.0;	//1Mohm
-	T rk = 1000.0;		//1kohm
+	ci = 0.0000001	;	//100nF
+	ck = 0.00001;		//10uF
+	co = 0.00000001;	//10nF
+	ro = 1000000.0;		//1Mohm
+	rp = 100000.0;		//100kohm
+	rg = 20000.0;		//20kohm
+	ri = 1000000.0;		//1Mohm
+	rk = 1000.0;		//1kohm
 	e = 250.0;		//250V
 
-	Vi = V(0.0,10000.0);	//1kohm internal resistance
-	Ci = C(ci, Fs);
-	Ck = C(ck, Fs);
-	Co = C(co, Fs);
-	Ro = R(ro);
-	Rg = R(rg);
-	Ri = R(ri);
-	Rk = R(rk);
-	E = V(e, rp);
-
-	//Official
-	//->Gate
-	S0 = ser(&Ci, &Vi);
-	I0 = inv(&S0);
-	P0 = par(&I0, &Ri);
-	S1 = ser(&Rg, &P0);
-	I1 = inv(&S1);
-
-	//->Cathode
-	I3 = par(&Ck, &Rk);
-
-	//->Plate
-	S2 = ser(&Co, &Ro);
-	I4 = inv(&S2);
-	P2 = par(&I4, &E);
-
 	// 12AX7 triode model RSD-1
-	v.g1 = 2.242e-3;
-	v.mu1 = 103.2;
-	v.gamma1 = 1.26;
-	v.c1 = 3.4;
-	v.gg1 = 6.177e-4;
-	v.e1 = 1.314;
-	v.cg1 = 9.901;
-	v.ig01 = 8.025e-8;
+	v.g = 2.242e-3;
+	v.mu = 103.2;
+	v.gamma = 1.26;
+	v.c = 3.4;
+	v.gg = 6.177e-4;
+	v.e = 1.314;
+	v.cg = 9.901;
+	v.ig0 = 8.025e-8;
 
 	// 12AX7 triode model EHX-1
 	v.g2 = 1.371e-3;
@@ -262,8 +233,11 @@ void ZamTubePlugin::activate()
 	v.cg2 = 11.99;
 	v.ig02 = 3.917e-8;
 
+	ckt.updateRValues(ci, ck, co, e, 0.0, rp, rg, ri, rk, ro, 10000.0, Fs, v);
+	ckt.warmup_tubes();
+
         fSamplingFreq = Fs;
-        fConst0 = float(min(192000, max(1, fSamplingFreq)));
+        fConst0 = float(std::fminf(192000.0, std::fmaxf(1.0, fSamplingFreq)));
         fConst1 = (2 * fConst0);
         fConst2 = faustpower<2>(fConst1);
         fConst3 = (3 * fConst1);
@@ -299,18 +273,6 @@ void ZamTubePlugin::activate()
 
 void ZamTubePlugin::run(const float** inputs, float** outputs, uint32_t frames)
 {
-	T tubetone = 0.f;
-	
-	//12AX7 triode tube mod
-	v.g = v.g1 + (v.g2-v.g1)*tubetone;
-	v.mu = v.mu1 + (v.mu2-v.mu1)*tubetone;
-	v.gamma = v.gamma1 + (v.gamma2-v.gamma1)*tubetone;
-	v.c = v.c1 + (v.c2-v.c1)*tubetone;
-	v.gg = v.gg1 + (v.gg2-v.gg1)*tubetone;
-	v.e = v.e1 + (v.e2-v.e1)*tubetone;
-	v.cg = v.cg1 + (v.cg2-v.cg1)*tubetone;
-	v.ig0 = v.ig01 + (v.ig02-v.ig01)*tubetone;
-
 	float 	fSlow0 = middle;
 	float 	fSlow1 = (1.3784375e-06f * fSlow0);
 	float 	fSlow2 = expf((3.4f * (bass - 1)));
@@ -978,63 +940,25 @@ void ZamTubePlugin::run(const float** inputs, float** outputs, uint32_t frames)
 	float 	fSlow664 = ((fSlow31 == 0) / fSlow651);
 
 	float tubeout = 0.f;
+	ckt.t.insane = insane;
 	
 	for (uint32_t i = 0; i < frames; ++i) {
 
 		//Step 1: read input sample as voltage for the source
-		float in = sanitize_denormal(inputs[0][i]);
+		double in = sanitize_denormal(inputs[0][i]);
 
 		// protect against overflowing circuit
 		in = fabs(in) < DANGER ? in : 0.f; 
 
-		Vi.e = in*from_dB(tubedrive);
-
-		//Step 2: propagate waves up to the triode and push values into triode element
-		I1.waveUp();
-		I3.waveUp();
-		P2.waveUp();
-		v.G.WD = sanitize_denormal(I1.WU);
-		v.K.WD = sanitize_denormal(I3.WU); 
-		v.P.WD = sanitize_denormal(P2.WU);
-		v.vg = v.G.WD;
-		v.vk = v.K.WD;
-		v.vp = v.P.WD;
-		v.G.PortRes = sanitize_denormal(I1.PortRes);
-		v.K.PortRes = sanitize_denormal(I3.PortRes);
-		v.P.PortRes = sanitize_denormal(P2.PortRes);
-
-		//Step 3: compute wave reflections inside the triode
-		T vg0, vg1, vp0, vp1;
-
-		vg0 = -10.0;
-		vg1 = 10.0;
-		v.vg = sanitize_denormal(v.zeroffg(vg0,vg1,TOLERANCE));
-		//v.vg = v.secantfg(&vg0,&vg1);
-
-		vp0 = e;
-		vp1 = 0.0;
-		if (insane > 0.5f) {
-			v.vp = sanitize_denormal(v.zeroffp_insane(vp0,vp1,TOLERANCE));
+		double ViE = in*from_dB(tubedrive);
+		tubeout = 10. * ckt.advanc(ViE) * from_dB(30. - tubedrive);
+		if (!ckt.on) {
+			tubeout = 0.0;
 		} else {
-			v.vp = sanitize_denormal(v.zeroffp(vp0,vp1,TOLERANCE));
+			tubeout = sanitize_denormal(tubeout);
 		}
-		//v.vp = v.secantfp(&vp0,&vp1);
-
-		v.vk = sanitize_denormal(v.ffk());
-
-		v.G.WU = sanitize_denormal(2.0*v.vg-v.G.WD);
-		v.K.WU = sanitize_denormal(2.0*v.vk-v.K.WD);
-		v.P.WU = sanitize_denormal(2.0*v.vp-v.P.WD);
-		
 		outputs[0][i] = in;
-		
-		tubeout = -Ro.Voltage()/e; //invert signal and rescale
-		tubeout = sanitize_denormal(tubeout);
 
-		P2.setWD(v.P.WU); 
-		I1.setWD(v.G.WU);
-		I3.setWD(v.K.WU);
-		
 		//Tone Stack sim
 		fRec0[0] = sanitize_denormal((float)tubeout - (fSlow16 * (((fSlow14 * fRec0[2]) + (fSlow13 * fRec0[1])) + (fSlow11 * fRec0[3]))));
 		fRec1[0] = sanitize_denormal((float)tubeout - (fSlow47 * (((fSlow45 * fRec1[2]) + (fSlow44 * fRec1[1])) + (fSlow42 * fRec1[3]))));
diff --git a/plugins/ZamTube/ZamTubePlugin.hpp b/plugins/ZamTube/ZamTubePlugin.hpp
index 4f3c9c7..343cf2a 100644
--- a/plugins/ZamTube/ZamTubePlugin.hpp
+++ b/plugins/ZamTube/ZamTubePlugin.hpp
@@ -17,7 +17,8 @@
 #define ZAMTUBEPLUGIN_HPP_INCLUDED
 
 #include "DistrhoPlugin.hpp"
-#include "wdf.h"
+#include "triode.h"
+#include "wdfcircuits.h"
 
 START_NAMESPACE_DISTRHO
 
@@ -39,33 +40,16 @@ class ZamTubePlugin : public Plugin
 {
 public:
 	Triode v;
-	V Vi;
-	C Ci;
-	C Ck;
-	C Co;
-	R Ro;
-	R Rg;
-	R Ri;
-	R Rk;
-	V E;
-
+	TubeStageCircuit ckt;
 	T e;
-
-	//Official
-	//->Gate
-	ser S0;
-	inv I0;
-	par P0;
-	ser S1;
-	inv I1;
-
-	//->Cathode
-	par I3;
-
-	//->Plate
-	ser S2;
-	inv I4;
-	par P2;
+	T ci;
+	T ck;
+	T co;
+	T ro;
+	T rp;
+	T rg;
+	T ri;
+	T rk;
 
         float   fConst0;
         float   fConst1;
diff --git a/plugins/ZamTube/ZamTubeUI.cpp b/plugins/ZamTube/ZamTubeUI.cpp
index efc7fe1..5caf9d2 100644
--- a/plugins/ZamTube/ZamTubeUI.cpp
+++ b/plugins/ZamTube/ZamTubeUI.cpp
@@ -58,25 +58,25 @@ ZamTubeUI::ZamTubeUI()
 
     fKnobBass = new ZamKnob(this, knobImage);
     fKnobBass->setAbsolutePos(63, 140.5);
-    fKnobBass->setRange(0.f, 1.0f);
+    fKnobBass->setRange(0.f, 0.5f);
     fKnobBass->setLabel(true);
     fKnobBass->setScrollStep(0.1f);
-    fKnobBass->setDefault(0.5f);
+    fKnobBass->setDefault(0.0f);
     fKnobBass->setRotationAngle(240);
     fKnobBass->setCallback(this);
 
     fKnobMids = new ZamKnob(this, knobImage);
     fKnobMids->setAbsolutePos(63, 87);
-    fKnobMids->setRange(0.f, 1.0f);
+    fKnobMids->setRange(0.f, 0.5f);
     fKnobMids->setLabel(true);
     fKnobMids->setScrollStep(0.1f);
-    fKnobMids->setDefault(0.5f);
+    fKnobMids->setDefault(0.0f);
     fKnobMids->setRotationAngle(240);
     fKnobMids->setCallback(this);
 
     fKnobTreb = new ZamKnob(this, knobImage);
     fKnobTreb->setAbsolutePos(63, 33);
-    fKnobTreb->setRange(0.f, 1.0f);
+    fKnobTreb->setRange(0.f, 0.5f);
     fKnobTreb->setLabel(true);
     fKnobTreb->setScrollStep(0.1f);
     fKnobTreb->setDefault(0.0f);
@@ -151,8 +151,8 @@ void ZamTubeUI::programLoaded(uint32_t index)
 
     // Default values
     fKnobTube->setValue(0.0f);
-    fKnobBass->setValue(0.5f);
-    fKnobMids->setValue(0.5f);
+    fKnobBass->setValue(0.0f);
+    fKnobMids->setValue(0.0f);
     fKnobTreb->setValue(0.0f);
     fKnobGain->setValue(0.0f);
     fToggleInsane->setDown(false);
diff --git a/plugins/ZamTube/glue.h b/plugins/ZamTube/glue.h
new file mode 100644
index 0000000..2ee702b
--- /dev/null
+++ b/plugins/ZamTube/glue.h
@@ -0,0 +1,15 @@
+
+#include <vector>
+#include <cmath>
+#include <stdint.h>
+#include <stdlib.h>
+
+#define LOG_ERROR(x)
+#define LOG_INFO(x)
+#define LOG_SAMPLE2(x)
+#define LOG_INNER_LOOP(x)
+#define SCOPE(x, ...)
+
+#define Assert(x); if (!(x)) { printf("ASSERT FAILED!\n"); }
+typedef double Real;
+typedef uint32_t uint;
diff --git a/plugins/ZamTube/inv.h b/plugins/ZamTube/inv.h
deleted file mode 100644
index bb63f41..0000000
--- a/plugins/ZamTube/inv.h
+++ /dev/null
@@ -1,21 +0,0 @@
-//inv.h
-#ifndef __INV_H__
-#define __INV_H__
-
-class ser;
-class par;
-
-#include "wdf.h"
-#include "ser.h"
-#include "par.h"
-
-class inv : public Adaptor { 
-public:
-	template <class Port>inv(Port *left);
-	inv(ser *l);
-	T waveUp();
-	void setWD(T waveparent);
-};
-
-#else
-#endif
diff --git a/plugins/ZamTube/par.h b/plugins/ZamTube/par.h
deleted file mode 100644
index 8373971..0000000
--- a/plugins/ZamTube/par.h
+++ /dev/null
@@ -1,23 +0,0 @@
-//par.h
-#ifndef __PAR_H__
-#define __PAR_H__
-
-class inv;
-class ser;
-
-#include "wdf.h"
-#include "inv.h"
-#include "ser.h"
-
-class par : public Adaptor {
-public:
-	template <class Port1, class Port2>par(Port1 *left, Port2 *right);
-	par(inv *l, R* r);
-	par(inv *l, V* r);
-	par(C *l, R* r);
-	T waveUp();
-	void setWD(T waveparent);
-};
-
-#else
-#endif
diff --git a/plugins/ZamTube/ser.h b/plugins/ZamTube/ser.h
deleted file mode 100644
index 319c3eb..0000000
--- a/plugins/ZamTube/ser.h
+++ /dev/null
@@ -1,22 +0,0 @@
-//ser.h
-#ifndef __SER_H__
-#define __SER_H__
-
-class inv;
-class par;
-
-#include "inv.h"
-#include "par.h"
-
-class ser : public Adaptor {
-public:
-	template <class Port1, class Port2>ser(Port1 *left, Port2 *right);
-	ser(R* l, par* r);
-	ser(C* l, R* r);
-	ser(C* l, V* r);
-	T waveUp();
-	void setWD(T waveparent);	
-};
-
-#else
-#endif
diff --git a/plugins/ZamTube/wdf.cpp b/plugins/ZamTube/triode.cpp
index 0aa0359..d310aca 100644
--- a/plugins/ZamTube/wdf.cpp
+++ b/plugins/ZamTube/triode.cpp
@@ -1,202 +1,62 @@
-//wdf.cpp
-
 #include <stdio.h>
 #include <inttypes.h>
 #include <cmath>
-#include "wdf.h"
+#include "triode.h"
 using std::abs;
 
-WDF::WDF() {}
-
-void WDF::setWD(T val) {
-	WD = val;
-	state = val;
-	DUMP(printf("DOWN\tWDF\t%c\tWD=%f\tWU=%f\tV=%f\n",type,WD,WU,(WD+WU)/2.0));	
-}
-
-void OnePort::setWD(T val) {
-	WD = val;
-	state = val;
-	DUMP(printf("DOWN\tOneport\t%c\tWD=%f\tWU=%f\tV=%f\n",type,WD,WU,(WD+WU)/2.0));	
-}
-
-T WDF::Voltage() {
-	T Volts = (WU + WD) / 2.0;
-	return Volts;
-}
-
-T WDF::Current() {
-	T Amps = (WU - WD) / (2.0*PortRes);
-	return Amps;
-}
-
-template <class Port1, class Port2>ser::ser(Port1 *l, Port2 *r) : Adaptor(THREEPORT) {
-	left = l;
-	right = r;
-	PortRes = l->PortRes + r->PortRes;
-	type = 'S';
-}
-
-ser::ser(R* l, par* r) : Adaptor(THREEPORT) {
-	left = l;
-	right = r;
-	PortRes = l->PortRes + r->PortRes;
-	type = 'S';
-}
-
-ser::ser(C* l, R* r) : Adaptor(THREEPORT) {
-	left = l;
-	right = r;
-	PortRes = l->PortRes + r->PortRes;
-	type = 'S';
-}
-
-ser::ser(C* l, V* r) : Adaptor(THREEPORT) {
-	left = l;
-	right = r;
-	PortRes = l->PortRes + r->PortRes;
-	type = 'S';
-}
-
-template <class Port>inv::inv(Port *l) : Adaptor(PASSTHROUGH) {
-	left = l;
-	PortRes = l->PortRes;
-	type = 'I';
-}
-
-inv::inv(ser *l) : Adaptor(PASSTHROUGH) {
-	left = l;
-	PortRes = l->PortRes;
-	type = 'I';
-}
-
-T ser::waveUp() {
-	//Adaptor::WU = -left->waveUp() - right->waveUp();
-	WDF::WU = -left->waveUp() - right->waveUp();
-	DUMP(printf("UP\tser\tWU=%f\tWD=%f\tV=%f\n",WU,WD,(WD+WU)/2.0));
-	return WU;
-}
-
-template <class Port1, class Port2>par::par(Port1 *l, Port2 *r) : Adaptor(THREEPORT) {
-	left = l;
-	right = r;
-	PortRes = 1.0 / (1.0 / l->PortRes + 1.0 / r->PortRes);
-	type = 'P';
-}
-
-par::par(inv* l, R* r) : Adaptor(THREEPORT) {
-	left = l;
-	right = r;
-	PortRes = 1.0 / (1.0 / l->PortRes + 1.0 / r->PortRes);
-	type = 'P';
-}
-
-par::par(inv* l, V* r) : Adaptor(THREEPORT) {
-	left = l;
-	right = r;
-	PortRes = 1.0 / (1.0 / l->PortRes + 1.0 / r->PortRes);
-	type = 'P';
+T Triode::getC(void)
+{
+	return Kb;
 }
 
-par::par(C* l, R* r) : Adaptor(THREEPORT) {
-	left = l;
-	right = r;
-	PortRes = 1.0 / (1.0 / l->PortRes + 1.0 / r->PortRes);
-	type = 'P';
+T Triode::getP(void)
+{
+	return Pb;
 }
 
-T par::waveUp() {
-	T G23 = 1.0 / left->PortRes + 1.0 / right->PortRes;
-	WDF::WU = (1.0 / left->PortRes)/G23*left->waveUp() + (1.0 / right->PortRes)/G23*right->waveUp();
-	DUMP(printf("UP\tpar\tWU=%f\tWD=%f\tV=%f\n",WU,WD,(WD+WU)/2.0));
-	return WU;
+T Triode::getG(void)
+{
+	return Gb;
 }
 
-Adaptor::Adaptor(int flag) {
-	WU = 0.0;
-	WD = 0.0;
-	switch (flag) {
-		case ONEPORT:
-			left = NULL;
-			right = NULL;
-			break;
-		case PASSTHROUGH:
-			right = NULL;
-			break;
-		default:
-		case THREEPORT:
-			break;
+void Triode::compute(T Pbb, T Gbb, T Kbb)
+{
+//	T Kb_o = Kb;
+//	T Gb_o = Gb;
+//	T Pb_o = Pb;
+
+//	Kb = (2.0*vk-Kbb);
+//	Gb = (2.0*vg-Gbb);
+//	Pb = (2.0*vp-Pbb);
+
+	Kb = Kbb;
+	Gb = Gbb;
+	Pb = Pbb;
+
+	//Step 3: compute wave reflections inside the triode
+	T vg0, vg1, vp0, vp1;
+
+	vg0 = -10.0;
+	vg1 = 10.0;
+	vg = sanitize_denormal(zeroffg(vg0, vg1, TOLERANCE));
+	//v.vg = v.secantfg(&vg0,&vg1);
+
+	vp0 = e;
+	vp1 = 0.0;
+	if (insane) {
+		vp = sanitize_denormal(zeroffp_insane(vp0,vp1,TOLERANCE));
+	} else {
+		vp = sanitize_denormal(zeroffp(vp0,vp1,TOLERANCE));
 	}
-}
-
-void ser::setWD(T waveparent) {
-	Adaptor::setWD(waveparent);
-	DUMP(printf("SER WP=%f\n",waveparent));
-	left->setWD(left->WU-(2.0*left->PortRes/(PortRes+left->PortRes+right->PortRes))*(waveparent+left->WU+right->WU));
-	right->setWD(right->WU-(2.0*right->PortRes/(PortRes+left->PortRes+right->PortRes))*(waveparent+left->WU+right->WU));
-}
+	//v.vp = v.secantfp(&vp0,&vp1);
 
-void par::setWD(T waveparent) {
-	Adaptor::setWD(waveparent);
-	DUMP(printf("PAR WP=%f\n",waveparent));
-	T p = 2.0*(waveparent/PortRes + left->WU/left->PortRes + right->WU/right->PortRes)/(1.0/PortRes + 1.0/left->PortRes + 1.0/right->PortRes);
+	vk = sanitize_denormal(ffk());
 
-	left->setWD((p - left->WU));
-	right->setWD((p - right->WU));
-}
-
-T inv::waveUp() {
-	///////////WD = -left->WD;
-	WU = -left->waveUp(); 	//-
-	DUMP(printf("UP\tinv\tWU=%f\tWD=%f\tV=%f\n",WU,WD,(WD+WU)/2.0));
-	return WU;
-}
-
-void inv::setWD(T waveparent) {
-	WDF::setWD(waveparent);
-	DUMP(printf("INV WP=%f\n",waveparent));
-	//left->WD = -waveparent;		//-
-	///////////left->WU = -WU;
-	left->setWD(-waveparent);	//-
-	
-}
-
-R::R(T res) : Adaptor(ONEPORT) {
-	PortRes = res;
-	type = 'R';
-}
+	Kb = (2.0*vk-Kb);
+	Gb = (2.0*vg-Gb);
+	Pb = (2.0*vp-Pb);
 
-T R::waveUp() {
-	WU = 0.0;
-	DUMP(printf("UP\tR\tWU=%f\tWD=%f\tV=%f\n",WU, WD,(WD+WU)/2.0));
-	return WU;
-}
-
-C::C(T c, T fs) : Adaptor(ONEPORT) {
-	PortRes = 1.0/(2.0*c*fs);
-	state = 0.0;
-	type = 'C';
-}
-
-T C::waveUp() {
-	WU = state;
-	DUMP(printf("UP\tC\tWU=%f\tWD=%f\tV=%f\n",WU,WD,(WD+WU)/2.0));
-	return WU;
-}
-
-V::V(T ee, T r) : Adaptor(ONEPORT) {
-	e = ee;
-	PortRes = r;
-	WD = 0.0;  //always?
-	type = 'V';
-}
-
-T V::waveUp() {
-	T watts = 100.0;
-	WU = 2.0*e - WD;
-	if (Voltage()*Current() > watts) WU *= 0.995;//0.9955;
-	DUMP(printf("UP\tV\tWU=%f\tWD=%f\tV=%f\n",WU, WD,(WD+WU)/2.0));
-	return WU;
+	//printf("vg = %f  vp = %f  vk = %f  Gb = %f  Kb = %f  Pb = %f\n", vg, vp, vk, Gb, Kb, Pb);
 }
 
 inline T _exp(const T x)
@@ -216,55 +76,50 @@ inline T _log(const T x)
     return 2.0*(a+a*a*a/3.0+a*a*a*a*a/5.0+a*a*a*a*a*a*a/7.0+a*a*a*a*a*a*a*a*a/9.0); 
 }
 
-inline T _pow(const T a, const T b)
-{
-    return pow(a,b);
-}
-
 T Triode::ffg(T VG) {
-        return (G.WD-G.PortRes*(gg*_pow(_log(1.0+_exp(cg*VG))/cg,e)+ig0)-VG);
+        return (Gb-Gr*(gg*pow(log(1.0+exp(cg*VG))/cg,e)+ig0)-VG);
 }
 
 T Triode::fgdash(T VG) {
         T a1 = exp(cg*VG);
         T b1 = -e*pow(log(a1+1.0)/cg,e-1.0);
-        T c1 = a1/(a1+1.0)*gg*G.PortRes;
+        T c1 = a1/(a1+1.0)*gg*Gr;
         return (b1*c1);
 }
 
-T Triode::ffp(T VP) { 
-    static bool prepared = false;
-    static double scale;
-    static double coeff[4];
-    if(!prepared) {
-        //go go series expansion
-        const double L2 = log(2.0);
-
-        const double scale = pow(L2,gamma-2)/(8.0*pow(c,gamma));
-        coeff[0] = 8.0*L2*L2*scale;
-        coeff[1] = gamma*c*L2*4*scale;
-        coeff[2] = (c*c*gamma*gamma+L2*c*c*gamma-c*c*gamma)*scale;
-        coeff[3] = 0.0;
-        prepared = true;
-    }
-
-    double A = VP/mu+vg;
-    return (P.WD+P.PortRes*((g*(coeff[0]+coeff[1]*A+coeff[2]*A*A))+(G.WD-vg)/G.PortRes)-VP);
+T Triode::ffp(T VP) {
+	static bool prepared = false;
+	static double scale;
+	static double coeff[4];
+	if(!prepared) {
+		//go go series expansion
+		const double L2 = log(2.0);
+
+		const double scale = pow(L2,gamma-2)/(8.0*pow(c,gamma));
+		coeff[0] = 8.0*L2*L2*scale;
+		coeff[1] = gamma*c*L2*4*scale;
+		coeff[2] = (c*c*gamma*gamma+L2*c*c*gamma-c*c*gamma)*scale;
+		coeff[3] = 0.0;
+		prepared = true;
+	}
+
+	double A = VP/mu+vg;
+	return (Pb+Pr*((g*(coeff[0]+coeff[1]*A+coeff[2]*A*A))+(Gb-vg)/Gr)-VP);
 }
 
 T Triode::ffp_insane(T VP) {
-	return (P.WD+P.PortRes*((g*pow(log(1.0+exp(c*(VP/mu+vg)))/c,gamma))+(G.WD-vg)/G.PortRes)-VP);
+	return (Pb+Pr*((g*pow(log(1.0+exp(c*(VP/mu+vg)))/c,gamma))+(Gb-vg)/Gr)-VP);
 }
 
 T Triode::fpdash(T VP) {
         T a1 = exp(c*(vg+VP/mu));
         T b1 = a1/(mu*(a1+1.0));
-        T c1 = g*gamma*P.PortRes*pow(log(a1+1.0)/c,gamma-1.0);
+        T c1 = g*gamma*Pr*pow(log(a1+1.0)/c,gamma-1.0);
         return (c1*b1);
 }
 
 T Triode::ffk() {
-        return (K.WD - K.PortRes*(g*pow(log(1.0+exp(c*(vp/mu+vg)))/c,gamma)));
+        return (Kb - Kr*(g*pow(log(1.0+exp(c*(vp/mu+vg)))/c,gamma)));
 }
 /*
 T Triode::secantfg(T *i1, T *i2) {
@@ -362,6 +217,11 @@ T Triode::r8_abs ( T x )
 
 Triode::Triode()
 {
+	vg = 0.0;
+	vk = 0.0;
+	vp = 0.0;
+	insane = false;
+
 	T r = 1.0;
 
 	while ( 1.0 < ( T ) ( 1.0 + r )	)
diff --git a/plugins/ZamTube/wdf.h b/plugins/ZamTube/triode.h
index a315f97..fd2da51 100644
--- a/plugins/ZamTube/wdf.h
+++ b/plugins/ZamTube/triode.h
@@ -21,13 +21,6 @@ Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
 
 #include <math.h>
 
-#define ONEPORT 0
-#define PASSTHROUGH 1
-#define THREEPORT 2
-
-#define max(x,y) (( (x) > (y) ) ? x : y )
-#define min(x,y) (( (x) < (y) ) ? x : y )
-#define sign(x) ( (x) >= 0.0 ? 1.0 : -1.0 )
 #define BIG 1e12
 #define SMALL 1e-14
 #define EPSILON 1e-9
@@ -36,44 +29,26 @@ Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
 #define SWAP_PN(x,y) {T tmp=y; y=-x; x=tmp;}
 #define SWAP_NP(x,y) {T tmp=y; y=x; x=-tmp;}
 #define SWAP_NN(x,y) {T tmp=y; y=-x; x=-tmp;}
-
-#if 0
-#include <stdio.h>
-#define DUMP(x) x
-#else
-#define DUMP(x)
-#endif
+#define TOLERANCE 1e-6
 
 typedef double T;
 
-class WDF {
-public:
-	T WD;
-	T WU;
-	T PortRes;
-	WDF();
-	T Voltage();
-	T Current();
-	T state;
-	char type;
-
-	virtual T waveUp() { return 0.0; }
-	virtual void setWD(T waveparent);
-};
-
-class OnePort : public WDF {
-public:
-	void setWD(T waveparent);
-};
+static inline float
+sanitize_denormal(float v) {
+        if(!std::isnormal(v) || !std::isfinite(v))
+                return 0.f;
+        return v;
+}
 
 class Triode {
 public:
-	WDF G, K, P;
+	T Kb, Gb, Pb;
+	T Kr, Gr, Pr;
 
-        T vg, vk, vp;
-        T g, mu, gamma, c, gg, e, cg, ig0;
-        T g1, mu1, gamma1, c1, gg1, e1, cg1, ig01;
-        T g2, mu2, gamma2, c2, gg2, e2, cg2, ig02;
+	T vg, vk, vp;
+	T g, mu, gamma, c, gg, e, cg, ig0;
+	T g1, mu1, gamma1, c1, gg1, e1, cg1, ig01;
+	T g2, mu2, gamma2, c2, gg2, e2, cg2, ig02;
 
 	T ffg(T VG);
 	T fgdash(T VG);
@@ -85,8 +60,13 @@ public:
 	T newtonfg(T *i1);
 	T secantfp(T *i1, T *i2);
 	T newtonfp(T *i1);
+	bool insane;
 	
 	Triode();
+	void compute(T Kbb, T Gbb, T Pbb);
+	T getC(void);
+	T getG(void);
+	T getP(void);
 
 	//Brent's method
 	T r8_abs ( T x );
@@ -96,41 +76,6 @@ public:
 	T zeroffp ( T a, T b, T t );
 	T zeroffp_insane ( T a, T b, T t );
 	T zeroffg ( T a, T b, T t );
-
 };
 
-class Adaptor : public OnePort {
-public:
-	WDF *left;
-	WDF *right;
-	Adaptor(int flag);
-};
-
-class R : public Adaptor {
-public:
-	R(T res);
-	T waveUp();
-};
-
-class C : public Adaptor {
-public:
-	C(T c, T fs);
-	T waveUp();
-};
-
-class V : public Adaptor {
-public:
-	T e;
-	V(T ee, T r);
-	T waveUp();
-};
-
-class inv;
-class ser;
-class par;
-#include "inv.h"
-#include "ser.h"
-#include "par.h"
-
-#else
 #endif
diff --git a/plugins/ZamTube/wdfcircuits.h b/plugins/ZamTube/wdfcircuits.h
new file mode 100644
index 0000000..6be9f2a
--- /dev/null
+++ b/plugins/ZamTube/wdfcircuits.h
@@ -0,0 +1,137 @@
+#ifndef WDFCIRCUITS_H
+#define WDFCIRCUITS_H
+#include "glue.h"
+#include "triode.h"
+
+class TubeStageCircuit {
+	/*Tube Preamp*/
+public:
+	Triode t;
+	bool on;
+
+	TubeStageCircuit() {
+		Cia = 0.0;
+		Cka = 0.0;
+		Coa = 0.0;
+		on = false;
+	}
+	TubeStageCircuit(Real C_Ci, Real C_Ck, Real C_Co, Real E_E250, Real E_Vi, Real R_E250, Real R_Rg, Real R_Ri, Real R_Rk, Real R_Ro, Real R_Vi, Real sampleRate, Triode& tube) {
+		Cia = 0.0;
+		Cka = 0.0;
+		Coa = 0.0;
+		on = false;
+		updateRValues(C_Ci, C_Ck, C_Co, E_E250, E_Vi, R_E250, R_Rg, R_Ri, R_Rk, R_Ro, R_Vi, sampleRate, tube);
+	}
+
+	void updateRValues(Real C_Ci, Real C_Ck, Real C_Co, Real E_E250, Real E_Vi, Real R_E250, Real R_Rg, Real R_Ri, Real R_Rk, Real R_Ro, Real R_Vi, Real sampleRate, Triode& tube) {
+		t = tube;
+		Real ViR = R_Vi;
+		ViE = E_Vi;
+		Real CiR = 1.0 / (2.0*C_Ci*sampleRate);
+		Real CkR = 1.0 / (2.0*C_Ck*sampleRate);
+		Real CoR = 1.0 / (2.0*C_Co*sampleRate);
+		Real RoR = R_Ro;
+		Real RgR = R_Rg;
+		Real RiR = R_Ri;
+		Real RkR = R_Rk;
+		Real E250R = R_E250;
+		E250E = E_E250;
+		Real S0_3R = (CiR + ViR);
+		S0_3Gamma1 = CiR/(CiR + ViR);
+		Assert(S0_3Gamma1 >= 0.0 && S0_3Gamma1 <= 1.0);
+		Real P0_1R = S0_3R;
+		Real P0_2R = RiR;
+		Real P0_3R = 1.0 /(1.0 / P0_1R + 1.0 / P0_2R);
+		P0_3Gamma1 = 1.0 / P0_1R/(1.0 / P0_1R + 1.0 / P0_2R);
+		Assert(P0_3Gamma1 >= 0.0 && P0_3Gamma1 <= 1.0);
+		S1_3Gamma1 = RgR/(RgR + P0_3R);
+		Assert(S1_3Gamma1 >= 0.0 && S1_3Gamma1 <= 1.0);
+		Real I3_1R = CkR;
+		Real I3_2R = RkR;
+		I3_3Gamma1 = 1.0 / I3_1R/(1.0 / I3_1R + 1.0 / I3_2R);
+		Assert(I3_3Gamma1 >= 0.0 && I3_3Gamma1 <= 1.0);
+		Real S2_3R = (CoR + RoR);
+		S2_3Gamma1 = CoR/(CoR + RoR);
+		Assert(S2_3Gamma1 >= 0.0 && S2_3Gamma1 <= 1.0);
+		Real P2_1R = S2_3R;
+		Real P2_2R = E250R;
+		P2_3Gamma1 = 1.0 / P2_1R/(1.0 / P2_1R + 1.0 / P2_2R);
+		Assert(P2_3Gamma1 >= 0.0 && P2_3Gamma1 <= 1.0);
+		t.Kr = sanitize_denormal(I3_3Gamma1);
+		t.Pr = sanitize_denormal(S2_3Gamma1);
+		t.Gr = sanitize_denormal(S1_3Gamma1);
+		//printf("Kr = %f  Pr = %f  Gr = %f\n", t.Kr, t.Pr, t.Gr);
+	}
+
+	void warmup_tubes(void) {
+		int i;
+		on = false;
+		for (i = 0; i < 8000; i++) {
+			advanc(0.0);
+		}
+		on = true;
+	}
+
+	Real advanc(Real VE){
+		ViE = VE;
+		Real Ckb = Cka;
+		Real I3_3b3 = -I3_3Gamma1*(-Ckb);
+		Real Cib = Cia;
+		Real S0_3b3 = -( Cib) + -( ViE);
+		Real P0_3b3 = -P0_3Gamma1*(-S0_3b3);
+		Real S1_3b3 = -( 0.0) + -( P0_3b3);
+		Real Cob = Coa;
+		Real S2_3b3 = -( Cob) + -( 0.0);
+		Real P2_3b3 = E250E - P2_3Gamma1*(E250E - S2_3b3);
+		//Tube:    K       G      P
+		t.compute(I3_3b3,S1_3b3,P2_3b3);
+		Real b1 = t.getC();
+		Real b2 = t.getG();
+		Real b3 = t.getP();
+		//Set As
+		Real I3_3b1 = b1  - Ckb - I3_3Gamma1*(-Ckb);
+		Cka = I3_3b1;
+		Real S1_3b2 = -( -( 0.0) + -( b2) - S1_3Gamma1*(-( 0.0) + -( P0_3b3) + -( b2)));
+		Real P0_3b1 = S1_3b2  - S0_3b3 - P0_3Gamma1*(-S0_3b3);
+		Real S0_3b1 = -( -( Cib) - S0_3Gamma1*(-( Cib) + -( ViE) + -( P0_3b1)));
+		Cia = S0_3b1;
+		Real P2_3b1 = b3 + E250E - S2_3b3 - P2_3Gamma1*(E250E - S2_3b3);
+		Real S2_3b1 = -( -( Cob) - S2_3Gamma1*(-( Cob) + -( 0.0) + -( P2_3b1)));
+		Coa = S2_3b1;
+		Real S2_3b2 = -( -( Cob) + -( P2_3b1) - S2_3Gamma1*(-( Cob) + -( 0.0) + -( P2_3b1)));
+		Real Roa = S2_3b2;
+		return -(Roa);
+	}
+
+	std::vector<Real> getState(){
+		std::vector<Real> state(3, 0.0);
+		state[0] = Cia;
+		state[1] = Cka;
+		state[2] = Coa;
+		return state;
+	}
+	void setState(std::vector<Real> state) {
+		Assert(state.size() == 3);
+		Cia = state[0];
+		Cka = state[1];
+		Coa = state[2];
+	}
+
+private:
+	//State variables
+	Real Cia;
+	Real Cka;
+	Real Coa;
+
+	//R values
+	Real I3_3Gamma1;
+	Real E250E;
+	Real ViE;
+	Real P0_3Gamma1;
+	Real S0_3Gamma1;
+	Real S2_3Gamma1;
+	Real S1_3Gamma1;
+	Real P2_3Gamma1;
+};
+
+#endif