// draw the shape of a probability density function
// generate a time-series and play as audification
//
// C. Frauenberger, April 2006
//
// see http://swiki.hfbk-hamburg.de:8888/MusicTechnology/709#3

(
q = q ? ();	
p = p ?? { ProxySpace.new }; p.push;	
)

(
	// size for mapping, histogramm and pdf 
q.tabSize = 1024;
	// alloc Buffer for mapping table
q.mapBuf = Buffer.alloc(s, 2 * q.tabSize, 1);
	// alloc Buffer for histogram
q.histoBuf = Buffer.alloc(s, q.tabSize, 1);
	// alloc Buffer for spectrum
q.spectrumBuf = Buffer.alloc(s, q.tabSize, 1);
	// functions for the pdf
q.dfunc = { arg i; i = i / q.tabSize; (sin(i * 80) * 0.3 + sin(i * 23)).max(0) };
)

(
	// the pdf and mapping table
q.pdf = Array.fill(q.tabSize, q.at(\dfunc)).normalize;
q.mapTab = q.pdf.asRandomTable;

	// the views:
q.width = 400;
q.height = 150;
q.win = SCWindow("pdf / map table / histogram / spectrum", Rect(100, 600, q.width+8, 4*q.height + 20));
q.win.view.decorator =  FlowLayout(q.win.view.bounds);
q.pdfPlot = SCMultiSliderView(q.win, Rect(0, 0, q.width, q.height)); 
q.mapPlot = SCMultiSliderView(q.win, Rect(0, q.height + 5, q.width, q.height));
q.histoPlot = SCMultiSliderView(q.win, Rect(0, 2*q.height + 5, q.width, q.height));
q.spectrumPlot = SCMultiSliderView(q.win, Rect(0, 3*q.height + 5, q.width, q.height));
[q.pdfPlot, q.mapPlot, q.histoPlot, q.spectrumPlot].do{ |i|
	i.drawLines_(true);
	i.indexThumbSize_(1);
	i.valueThumbSize_(1);
	i.xOffset_(0);
};
q.pdfPlot.value_(q.pdf.resamp1(q.width));

q.plotUpdate = { 
	var sig, linMap;
	q.pdf = q.pdfPlot.value.resamp1(q.tabSize);
	q.mapTab = q.pdf.asRandomTable;
	q.mapPlot.value_(q.mapTab.resamp1(q.width));
	linMap = q.mapTab.linlin(0,1,-1,1);
	sig = Signal.fill(q.tabSize, { |i| linMap[i] });
	q.mapBuf.sendCollection(sig.asWavetable);
};

q.pdfPlot.mouseUpAction = { q.plotUpdate.value };
q.plotUpdate.value;
q.win.front;	
)

	// The audification
(
~player.ar(2);
~histo.ar(1);
~spectrum.ar(1);
)
(
~player = { 
	//var source = TRand.ar(0.0, 0.5,Dust.ar(rate));
	var source = WhiteNoise.ar;
	Shaper.ar(q.mapBuf.bufnum, source) ! 2;
	//source ! 2;
};

~histo = { BufWr.ar(~player.ar(1,0), 
			q.histoBuf.bufnum, 
			Phasor.ar(0, BufRateScale.kr(q.histoBuf.bufnum), 0, BufFrames.kr(q.histoBuf.bufnum)), 
			1);
		0.0; //quiet
};

~spectrum = { FFT(q.spectrumBuf.bufnum, ~player.ar(1,0)); 0.0 };

)
ProxyMixer(p);

	// Histogram and Spectrum
(
var histMax = 20;
var specMax = 1000;
var numBins = q.width;
var histogram = Array.fill(numBins, 0);
var rn;
Tdef(\histogram, {
	loop{
		q.histoBuf.getToFloatArray(action: { |randomNumbers|
			rn = randomNumbers;
			histogram = randomNumbers.histo(q.width, -1, 1);
			{ q.histoPlot.value_(histogram/histMax) }.defer;
		});
		q.spectrumBuf.getToFloatArray(action: { |spec|
			var z, x, y;
			z = spec.clump(2).flop;
			z = [Signal.newFrom(z[0]), Signal.newFrom(z[1])];
			x = Complex(z[0], z[1]);
			{ q.spectrumPlot.value_(x.magnitude.pow(2).resamp1(q.width)/specMax) }.defer;
		});
		// calculate moments
		q.win.drawHook = {
			("Mean: " ++ rn.mean.round(0.01).asString).drawAtPoint(10@(2*q.height + 20));			("Variance: " ++ rn.variance.round(0.01).asString).drawAtPoint(90@(2*q.height + 20));
			("Skew: " ++ rn.skew.round(0.01).asString).drawAtPoint(190@(2*q.height + 20));
			("Kurtosis: " ++ rn.kurtosis.round(0.01).asString).drawAtPoint(270@(2*q.height + 20));
			("Mean: " ++ rn.mean.round(0.01).asString).drawAtPoint(10@20);			("Variance: " ++ rn.variance.round(0.01).asString).drawAtPoint(90@20);
			("Skew: " ++ rn.skew.round(0.01).asString).drawAtPoint(190@(2*q.height + 20));
			("Kurtosis: " ++ rn.kurtosis.round(0.01).asString).drawAtPoint(270@(2*q.height + 20));

		};
		{ q.win.refresh; }.defer;
		0.5.wait;
	}
});
)
Tdef(\histogram).play;
Tdef(\histogram).stop;


// different pdfs

(	// normal distribution
q.pdfPlot.value_(Array.fill(q.width, { |i|
	var x = 2*i/q.width - 1;
	var v = 0.25;
	1/((2pi*v).sqrt) * (-1 * (x-0).squared / (2 * v)).exp;
}).normalize);
q.plotUpdate
)

(	// exponential
q.pdfPlot.value_(Array.fill(q.width, { |i|
	var x = 2*i/q.width - 1;
	var l = 3;
	l * (-1*x.abs*l).exp;
}).normalize);
q.plotUpdate
)

(	// skewy exponential
q.pdfPlot.value_(Array.fill(q.width, { |i|
	var x = 2*i/q.width - 1.8;
	var l = 3;
	l * (-1*x.abs*l).exp;
}).normalize);
q.plotUpdate
)

(	// skewy cos
q.pdfPlot.value_(Array.fill(q.width, { |i|
	var x = pi*i/q.width;
	1- cos(x);
}).normalize);
q.plotUpdate
)

(	// levels
q.pdfPlot.value_(
	Array.fill(q.width/2, 0.15) ++ 
	Array.fill(q.width/4, 0.4) ++
	Array.fill(q.width/4, 0.03)
);
q.plotUpdate
)


(	// uniform 
q.pdfPlot.value_(Array.fill(q.width, 0.1));
q.plotUpdate
)


// IFFT approach 

	// t.. time, u[4].. raw moments
(
q.cfunc = { |ts, u| 
	ts.collect({ |t| 
		Complex.new(1 - (1/2*t.pow(2)*u[1].pow(2)) + (1/24*t.pow(4)*u[3].pow(4)), 
				(t*u[0])-(1/6*t.pow(3)*u[2].pow(3)));
	});
};

	// pdf by ifft of cfunc
q.calcPdf = { |u| 
	var t = Array.fill(q.tabSize, { |i| i / 44100});
	var b = Buffer.alloc(s, 2 * q.tabSize, 1);
	var c = q.at(\cfunc).value(t,u).collect({ |x| [x.real, x.imag] }).flatten;
	b.sendCollection(c);
	IFFT(b);
	b.getToFloatArray(action:{|a| 
		q.cpdf = a.clump(2).collect({|x| Complex.new(x[0], x[1])});
		{ q.cpdf.magnitude.plot; q.cpdf.phase.plot; }.defer;
	});
};
)
	// experiment
q.at(\calcPdf).value([0, 0.1, 0.5, 0]);
q.pdfPlot.value_(q.cpdf.magnitude.resamp1(q.width).normalize);
q.plotUpdate;