chris@164: /*
chris@164:  *
chris@164:  * Simple 1-Dimensional Waveguide
chris@164:  *
chris@164:  */
chris@164: 
chris@164: #include "Waveguide.h"
chris@164: #include "../include/Utilities.h"
chris@164: #include <rtdk.h>
chris@164: #include <cmath>
chris@164: #include <stdio.h>
chris@164: #include <cstdlib>
chris@164: 
chris@164: #define DECAY 0.995//0.999
chris@164: #define DAMPING 0.01//0.05
chris@164: 
chris@164: // TODO: make damping and decay parametrisable
chris@164: 
chris@164: Waveguide::Waveguide()	{
chris@164: 
chris@164: 	// initialize variables
chris@164: 	a1_lp = 0;
chris@164: 	a2_lp = 0;
chris@164: 	b0_lp = 0;
chris@164: 	b1_lp = 0;
chris@164: 	b2_lp = 0;
chris@164: 	_dt = 1.0/44100.0;
chris@164: 	setFrequency(440);
chris@164: 	updateFilterCoeffs(8000);
chris@164: 	_filterReadPtr=0;
chris@164: 	for(int i=0;i<FILTER_BUFFER_SIZE;i++)	{
chris@164: 		_filterBuffer_x[i] = 0;
chris@164: 		_filterBuffer_y[i] = 0;
chris@164: 	}
chris@164: 	for(int i=0;i<WG_BUFFER_SIZE;i++)	{
chris@164: 		_buffer[i] = 0;
chris@164: 	}
chris@164: 	_lastX = 0;
chris@164: 	_lastY = 0;
chris@164: 	_readPtr = 0;
chris@164: 
chris@164: }
chris@164: 
chris@164: void Waveguide::setup()	{
chris@164: 
chris@164: }
chris@164: 
chris@164: float Waveguide::update(float in)	{
chris@164: 
chris@164: 	// 1. advance delay buffer read pointer
chris@164: 
chris@164: 	if(++_readPtr>=WG_BUFFER_SIZE)
chris@164: 		_readPtr=0;
chris@164: 
chris@164: 	// 2. write input into buffer
chris@164: 
chris@164: 	_buffer[_readPtr] = in;
chris@164: 
chris@164: 	// 3. read delayed sample from buffer
chris@164: 
chris@164: 	float out = _buffer[(_readPtr-_periodInSamples+WG_BUFFER_SIZE)%WG_BUFFER_SIZE];
chris@164: 
chris@164: 	// 4. apply damping (low-pass) filter to output
chris@164: 
chris@164: 	if(++_filterReadPtr>=FILTER_BUFFER_SIZE)
chris@164: 			_filterReadPtr=0;
chris@164: 
chris@164: 	out = b0_lp*out +
chris@164: 					b1_lp*_filterBuffer_x[(_filterReadPtr-1+FILTER_BUFFER_SIZE)%FILTER_BUFFER_SIZE] +
chris@164: 					b2_lp*_filterBuffer_x[(_filterReadPtr-2+FILTER_BUFFER_SIZE)%FILTER_BUFFER_SIZE] -
chris@164: 					a1_lp*_filterBuffer_y[(_filterReadPtr-1+FILTER_BUFFER_SIZE)%FILTER_BUFFER_SIZE] -
chris@164: 					a2_lp*_filterBuffer_y[(_filterReadPtr-2+FILTER_BUFFER_SIZE)%FILTER_BUFFER_SIZE];
chris@164: 
chris@164: 	// 5. Simple high-pass filter to block DC-offset
chris@164: 	// y[n] = x[n] - x[n-1] + a * y[n-1]
chris@164: 	float gain = 0.9999;
chris@164: 	float temp = out;
chris@164: 	out = out - _lastX + gain * _lastY;
chris@164: 	_lastY = out;
chris@164: 	_lastX = temp;
chris@164: 
chris@164: 	// 6. Apply intensity damping
chris@164: 	out *= DECAY;
chris@164: 
chris@164: 	_filterBuffer_x[_filterReadPtr] = in;
chris@164: 	_filterBuffer_y[_filterReadPtr] = out;
chris@164: 
chris@164: 	return out;
chris@164: 
chris@164: }
chris@164: 
chris@164: void Waveguide::setFrequency(float frequency)	{
chris@164: 
chris@164: 	// NB: currently no interpolation, so may not be ideal for dynamically changing waveguide frequency
chris@164: 	_periodInMilliseconds = 1000.0/frequency;
chris@164: 	_periodInSamples = (int)(_periodInMilliseconds * 44.1);
chris@164: 
chris@164: }
chris@164: 
chris@164: void Waveguide::updateFilterCoeffs(float frequency)	{
chris@164: 
chris@164: 	// FIXME: Butterworth filter doesn't work very well,
chris@164: 	//		  using simple FIR in the meantime
chris@164: 
chris@164: 	a1_lp = 0;
chris@164: 	a2_lp = 0;
chris@164: 	b0_lp = 1.0 - DAMPING;
chris@164: 	b1_lp = DAMPING;
chris@164: 	b2_lp = 0;
chris@164: 
chris@164: 	/*
chris@164: 	// 'w' for sake of resembling lower-case 'omega'
chris@164: 	float w = 2.0 * M_PI * frequency;
chris@164: 	float t = _dt;
chris@164: 	// The Q for a 2nd-order Butterworth is sqrt(2)/2
chris@164: 	float q = 0.707;//sqrt(2.0)/2.0;
chris@164: 
chris@164: 	// low-pass filter coefficients
chris@164: 	float a0_lp = w*w*t*t + 2*(w/q)*t + 4.0;
chris@164: 	float k = 1.0/a0_lp;
chris@164: 	a1_lp = (2.0*w*w*t*t - 8.0) * k;
chris@164: 	a2_lp = (4.0 - (w/q)*2.0*t + w*w*t*t) * k;
chris@164: 	b0_lp = (w*w*t*t) * k;
chris@164: 	b1_lp = (2.0*w*w*t*t) * k;
chris@164: 	b2_lp = (w*w*t*t) * k;
chris@164: 	*/
chris@164: 
chris@164: }