Non-working second order system code

Saxo · July 7, 2022, 7:53pm

inspired by this video: Giving Personality to Procedural Animations using Math - YouTube in which the poster uses a second order system to add personality to animations, I wanted to use second order systems in order to add some gesture to control signals (i.e. taking in a signal used to control pitch, and causing it to overshoot/anticipate/resonate/lag in real time in a parametrized way). I wrote this function based on the video, but plotting it just results in DC at 0 any advice would be greatly appreciated!

(
{
	var x, y, dxdt, dydt, d2ydt2, k1, k2, k3, f, zeta, r;
	f = 1;
	zeta = 1;
	r = 1;
	t = 1/s.sampleRate;

	k1 = (zeta) / (pi*f);
	k2 = (1) / ((2*pi*f).squared);
	k3 = (r*zeta) / (2*pi*f);

	x = LFPulse.ar(1, 0.5);//step to measure step response
	dxdt = Slope.ar(x);

	y = LocalIn.ar;
	dydt = Slope.ar(y);
	d2ydt2 = Slope.ar(dydt);

	y = y + (t*dydt);
	dydt = dydt + (t * (x + (k3 * dxdt) - y - (k1 * dydt)) / k2);

	LocalOut.ar(y);
	y;
}.plot(1)
)

droptableuser · July 7, 2022, 7:58pm

i don’t have any specific help but it made me think of this quark

semiquaver · July 7, 2022, 10:11pm

hmmmm

I think you have the order of your calculations wrong.

so y = Local in
dydt = Slope.ar(y)
y=y + dydt/s.sampleRate
dydt you now redefine according to complicated formula but
now you send y out before using your new dydt for anything.

another issue that might keep you from getting results you want is dxdt is 0 (except one sample/ second where it is inf and one where it is -inf!)

finally be aware that signal sent out to LocalOut returns to LocalIn 1 block later (not one sample!)

Saxo · July 8, 2022, 12:30am

Attempted to fix what you mentioned, still getting 0 DC

(
~fbBusY =Bus.audio(s, 1);
~fbBusDY =Bus.audio(s, 1);
{
	var x, y, dxdt, dydt, k1, k2, k3, f, zeta, r, block, fbBusY = ~fbBusY, fbBusDY = ~fbBusDY;
	f = 0.59;
	zeta = 0.15;
	r = 2;
	t = 1/s.sampleRate;

	k1 = (zeta) / (pi*f);
	k2 = (1) / ((2*pi*f).squared);
	k3 = (r*zeta) / (2*pi*f);

	x = LFPulse.ar(1, 0.5).lag(0.01);//step to measure step response
	dxdt = Slope.ar(x);

	y = In.ar(fbBusY);
	dydt = In.ar(fbBusDY);
	y = y + (t*dydt);
	dydt = dydt + (t * (x + (k3 * dxdt) - y - (k1 * dydt)) / k2);
	Out.ar(fbBusY, Delay1.ar(y));
	Out.ar(fbBusDY, Delay1.ar(dydt));
	y;
}.plot(5);
)

semiquaver · July 8, 2022, 12:51am

use InFeedback.ar instead of In.ar !!

Sam_Pluta · July 8, 2022, 3:15pm

This at least doesn’t sit at 0:

(
{
	var x, y, dxdt, dydt, d2ydt2, k1, k2, k3, f, zeta, r;
	f = 1;
	zeta = 1;
	r = 1;
	t = 1/s.sampleRate;

	k1 = (zeta) / (pi*f);
	k2 = (1) / ((2*pi*f).squared);
	k3 = (r*zeta) / (2*pi*f);

	x = LFPulse.ar(1, 0.5);//step to measure step response
	dxdt = Slope.ar(x);

	#y, dydt = LocalIn.ar(2);
	d2ydt2 = Slope.ar(dydt);

	y = y + (t*dydt);
	dydt = dydt + (t * (x + (k3 * dxdt) - y - (k1 * dydt)) / k2);

	LocalOut.ar([y, dydt]);
	y;
}.plot(5)
)

Sam

jordan · July 8, 2022, 3:45pm

I don’t think this is possible to implement in sclang for the same reason you can’t implement a one pole filter - you’d need a block size of 1. It should be possible at control rate though.

When you calculate dy you’re calculating the entire block’s slope, and use that to calculate the next value, but you want to then use that value to alter the derivative next sample, which isn’t possible unless it happens to lie on a block boundary.

semiquaver · July 8, 2022, 5:49pm

I believe is an Fb1 class in miSCellaneous_lib which might offer single-sample feedback… @dkmayer

Saxo · July 8, 2022, 6:11pm

replacing In with InFeedback yields what appears to be a jittery exponential curve that just keeps increasing…

jordan · July 8, 2022, 7:28pm

So I also saw this video, thought it was very good and wanted something like this in supercollider.

So I made a plugin this evening! The supercollider side of things i particularly bare and you must wrap all args right now … PIDContoller.ar(t, DC.ar(4), DC.ar(0.3), DC.ar(1.4)).

There is a release here but its only for linux as thats where I am. Anyone wanna do the rest?

Huge thanks to this guide by madskjeldgaard nice and easy to follow!

dkmayer · July 8, 2022, 9:12pm

Interesting thread, thanks for the link! Here is how it can be done with Fb1_ODE from miSCellaneous_lib.


// reduce blockSize for faster SynthDef build, not necessary though

(
s.options.blockSize = 8;
s.reboot;
)


// define first order system with substitution y' = w from
// y + k1 y' + k2 y'' = x + k3 x'


(
Fb1_ODEdef(\sys, { |t, y, k1 = 0, k2 = 1, k3 = 1, x = 1, dxdt = 0|
    [
        { y[1] },
		{ (x / k2) + (k3 * dxdt / k2) - (y[0] / k2) - (k1 * y[1] / k2) }
    ]
// didn't think about the init values yet ...
}, 0, [0, 0], 1, 1);
)


(
SynthDef(\freqCtrlSmoother, { |sinFreq = 1000, f = 3, zeta = 0.5, r = 2, amp = 0.1|
	// variable transformation
	var k1 = zeta / (pi * f);
	var k2 = 1 / ((2 * pi * f) ** 2);
	var k3 = r * zeta / (2 * pi * f);

	// feed freq control into system
	var x = K2A.ar(sinFreq);
	var dxdt = Slope.ar(x);

	// important: we want no automatic DC leaker on the system solution
	var ode_solved = Fb1_ODE.ar(\sys, [k1, k2, k3, x, dxdt], leakDC: false);

	// get dampened control
	var freq = ode_solved[0];
	var audio = SinOsc.ar(freq * [1, 1.01], mul: amp);

	Out.ar(0, audio)
}, metadata: (
	specs: (
		sinFreq: [50, 2000, \exp, 0, 1000],
		f: [0, 100, 5, 0, 3],
		zeta: [0, 5, 3, 0, 0.15],
		r: [-5, 5, \lin, 0, 2],
		amp: [0, 0.5, \db, 0, 0.1]
	)
)).add
)

// change values of sinFreq
// check influence of params f, zeta, and r on sinFreq tracking

\freqCtrlSmoother.sVarGui.gui

dkmayer · July 9, 2022, 9:04pm

Two further remarks on this:

.) ar ODE solving might be an unnecessary luxury in this case. At least for low f values, Fb1_ODE.kr worked fine in my tests.

.) as such kind of control can be used on many occasions, a wrapper makes sense.

////////// save as .sc file in Extensions folder and recompile

SmoothControl : UGen {

	*initClass {
		Fb1_ODEdef(\smoothControl, { |t, y, k1 = 0, k2 = 1, k3 = 1, x = 1, dxdt = 0|
			[
				{ y[1] },
				{ (x / k2) + (k3 * dxdt / k2) - (y[0] / k2) - (k1 * y[1] / k2) }
			]
		}, 0, [0, 0], 1, 1);
	}

	*ar { |in, f = 3, zeta = 0.5, r = 2|
		var k1 = zeta / (pi * f);
		var k2 = 1 / ((2 * pi * f) ** 2);
		var k3 = r * zeta / (2 * pi * f);
		^Fb1_ODE.ar(\smoothControl, [k1, k2, k3, in, Slope.ar(in)], leakDC: false)[0]
	}

	*kr { |in, f = 3, zeta = 0.5, r = 2|
		var k1 = zeta / (pi * f);
		var k2 = 1 / ((2 * pi * f) ** 2);
		var k3 = r * zeta / (2 * pi * f);
		^Fb1_ODE.kr(\smoothControl, [k1, k2, k3, in, Slope.kr(in)], leakDC: false)[0]
	}

}

//////////

Then:


// significantly less ugens and CPU usage (especially with default blockSize 64)

(
SynthDef(\freqCtrlSmoother_kr, { |sinFreq = 1000, f = 3, zeta = 0.5, r = 2, amp = 0.1|
	var freq = SmoothControl.kr(sinFreq, f, zeta, r);
	var audio = SinOsc.ar(freq * [1, 1.01], mul: amp);
	Out.ar(0, audio)
}, metadata: (
	specs: (
		sinFreq: [50, 2000, \exp, 0, 1000],
		f: [0, 100, 5, 0, 3],
		zeta: [0, 5, 3, 0, 0.15],
		r: [-5, 5, \lin, 0, 2],
		amp: [0, 0.5, \db, 0, 0.1]
	)
)).add
)

// change values of sinFreq
// check influence of params f, zeta, and r on sinFreq tracking

\freqCtrlSmoother_kr.sVarGui.gui

For the ar variant you’d replace the line with SmoothControl:

	var freq = SmoothControl.ar(K2A.ar(sinFreq), f, zeta, r);

Saxo · July 9, 2022, 11:27pm

Thank you so much for this! having lots of fun patching with it right now. One strange thing I am running into: I tried replacing the unsmoothed frequency arg with a pitch tracker, but for whatever reason, it only works with ZeroCrossing and not Pitch! With Pitch the smoothed output migrates to an equilibrium and stays there, no matter what the changing input does. This only happens with Pitch, and when the unsmoothed freq is replaced by say a SinOsc.kr, or by ZeroCrossing it functions fine. Wondering if you may have an idea of why this is.

(//Pitch Ugen as input not working
SynthDef(\freqCtrlSmoother, { |f = 3, zeta = 0.5, r = 2, amp = 0.1|
	// variable transformation
	var k1 = zeta / (pi * f);
	var k2 = 1 / ((2 * pi * f) ** 2);
	var k3 = r * zeta / (2 * pi * f);

	var x = K2A.ar(Pitch.kr(SoundIn.ar())[0]);
	//var x = K2A.ar(ZeroCrossing.kr(A2K.kr(SoundIn.ar())));//works fine
	var dxdt = Slope.ar(x);

	var ode_solved = Fb1_ODE.ar(\sys, [k1, k2, k3, x, dxdt], leakDC: false);
	var freq = ode_solved[0];
	var audio = SinOsc.ar(freq * [1, 1.01], mul: amp);

	freq.poll; //goes to 440 and stays there with Pitch.kr

	Out.ar(0, audio)
}, metadata: (
	specs: (
		f: [0, 100, 5, 0, 3],
		zeta: [0, 5, 3, 0, 0.15],
		r: [-5, 5, \lin, 0, 2],
		amp: [0, 0.5, \db, 0, 0.1]
	)
)).add
)

\freqCtrlSmoother.sVarGui.gui

dkmayer · July 10, 2022, 12:05am

This issue appears if Pitch doesn’t detect a frequency. A workaround would be selecting a number (arbitrary or somehow with sample and hold) in this case.


(
SynthDef(\freqCtrlSmoother_pitch, { |f = 3, zeta = 0.5, r = 2, amp = 0.1|
	// variable transformation
	var k1 = zeta / (pi * f);
	var k2 = 1 / ((2 * pi * f) ** 2);
	var k3 = r * zeta / (2 * pi * f);
	var freqData = Pitch.kr(SoundIn.ar());
	var x = K2A.ar(Select.kr(freqData[1], [freqData[0], DC.kr(400)]));
	var dxdt = Slope.ar(x);

	var ode_solved = Fb1_ODE.ar(\sys, [k1, k2, k3, x, dxdt], leakDC: false);
	var freq = ode_solved[0];
	var audio = SinOsc.ar(freq * [1, 1.01], mul: amp);

	freq.poll; 

	Out.ar(0, audio)
}, metadata: (
	specs: (
		f: [0, 100, 5, 0, 3],
		zeta: [0, 5, 3, 0, 0.15],
		r: [-5, 5, \lin, 0, 2],
		amp: [0, 0.5, \db, 0, 0.1]
	)
)).add
)

\freqCtrlSmoother_pitch.sVarGui.gui

A rough workaround indeed, but Pitch is also a rather raw pitchtracker (would probably work better with Tartini from the plugins).

dkmayer · July 10, 2022, 12:10am

… oops, it should be

var x = K2A.ar(Select.kr(freqData[1], [DC.kr(400), freqData[0]]));