Allpass delay line interpolator

[hibrasil]

I'm looking for some source code for a delay line / resonant comb filter / karpus strong algo using allpass interpolation rather than linear interpolation. Better still would be something that uses "warped" all pass interpolation as mentioned in Jon Dattoro's EffectDesign : Part2 paper :

http://www.stanford.edu/~dattorro/EffectDesignPart2.pdf

i have seen the example in the stk (DelayA.cpp). Does anyone know of any others?
I don't have the dsp knowledge to translate the diagrams etc from the papers into code, was hoping there might be some existing open source projects to look at

Cheers,

Oli

[xoxos]

at risk of being ignorant of the finer points of this inquiry having not implemented this myself, you've oc checked julius orion smith?? tons of interpolation material, some on allpasses i'm sure.

the other place i check for dsp is helsinki university of technology.. google hut and your subject..?? :p

[hibrasil]

yeah, i checked those, but they don't provide source code, just diagrams and equations . I guess I should just try and learn how to read this stuff.

oli

[MadBrain]

WARNING: Some of the parameters don't work (have been turned off for testing and debugging purposes)

This is a 2-tap crossfading delay line with 1-pole allpass interpolation. It also automatically compensates for the delay caused by the internal filter. When the pipe frequency changes, it calculates the new length and allpass parameters and fades from one pipe to the other.

It's intended for use in Synthedit so the parameters are setup for those ranges.

Code: Select all

class Pipe
{
private:
	float mix;
	float mix_delta;
	float *buffer;
	int mask;
	int buffer_size; // must be a power of 2
	int write_pointer;
	int read_pointer[2];
	float allpass_parameter[2];
	float allpass_memory[2];
	float allpass_memory2[2];

	int pipe_type;
	float feedback;
	float filter_parameter;
	float filter_parameter2;
	float filter_memory;
	float filter_memory2;
	float filter_gain;
	float filter_delay;
	float filter_old;
	float frequency;
	float fade_time;
	float glide;
	float glide_freq;
	float register_key;
	float delay_time_estimate;
	
public:
	Pipe();
	~Pipe();
	void set_max_length(int length); // Also clears the delay line!
	void try_max_length(int length); // Doesn't clear the delay line if the length is the same
	inline void set_frequency(float freq);
	void set_pipe_type(short in);
	void set_feedback(float in);
	void set_cutoff(float in);
	void set_fade_time(float in);
	void set_glide(float in);
	void set_register_key(float in);
	float output_signal();
	void input_signal(float input);
	float calculate_decay();
	void park();
};

Pipe::Pipe()
{
	// Initialize in a stable state
	mix = 0;
	mix_delta = 0;
	buffer = new float[1];
	mask = 0x1;
	buffer_size = 1;
	write_pointer = 0;
	read_pointer[0] = 0;
	read_pointer[1] = 0;
	allpass_parameter[0] = 0;
	allpass_parameter[1] = 0;
	allpass_memory[0] = 0;
	allpass_memory[1] = 0;
	feedback = 0;
	filter_parameter = 0;
	filter_parameter2 = 0;
	filter_memory = 0;
	filter_memory2 = 0;
	filter_gain = 0;
	filter_delay = 0;
	filter_old = 256124613613461.0f;
	frequency = -1000000000000000.0f;
	fade_time = 0.01f;
	glide = 1000;
	register_key = 2;
	glide_freq = 10000000000000000.0f;
	delay_time_estimate = 1000;
}

Pipe::~Pipe()
{
	delete [] buffer;
}

void Pipe::try_max_length(int length)
{
	if(length <= 0)
		length = 1;
	if(length >= 0x10000000)
		length = 0x10000000;
	int size;
	for(size=1; size<=length; size<<=1); // Get smallest pow(2) size in which length fits
	
	if(size!=buffer_size) // If buffer size different, reallocate
	{
		set_max_length(length);
	}
}

void Pipe::set_max_length(int length)
{
	

	if(length <= 0)
		length = 1;
	if(length >= 0x10000000)
		length = 0x10000000;
	int size;
	for(size=1; size<=length; size<<=1); // Get smallest pow(2) size in which length fits
	
	if(size!=buffer_size) // If buffer size different, reallocate
	{
		delete [] buffer;
		buffer = new float [size];
		buffer_size = size;	
	}

	// Reset pipe state
	write_pointer = 0;
	read_pointer[0] = 0;
	read_pointer[1] = 0;
	mask = size-1;
	allpass_parameter[0] = 0;
	allpass_parameter[1] = 0;
	allpass_memory[0] = 0;
	allpass_memory[1] = 0;
	allpass_memory2[0] = 0;
	allpass_memory2[1] = 0;
	mix = 0;
	mix_delta = 0;
	for(int i=0;i<buffer_size;i++) buffer[i]=0; return; // Clear buffer
	frequency = -1000000000000000.0f;
	return;
}

inline void Pipe::set_frequency(float freq)
{
	int dl;
	float freq_in = freq;
	int j;
	
	if(freq==frequency)
		return;

	if(mix<=0)
	{
		if(glide_freq >= 10000000000000000)
			j=1;
		else
			j=0;
		
		if (mix<-1)
			mix = 1;
		else
			mix += 1;
		mix_delta = pow(2,12.0-fade_time*12)/(float)(getSampleRate());
		read_pointer[0] = read_pointer[1];
		allpass_parameter[0] = allpass_parameter[1];
		allpass_memory[0] = allpass_memory[1];
		allpass_memory2[0] = allpass_memory2[1];
		
		freq *= 10;
		freq -= 5;
		
		freq -= register_key;
		switch(pipe_type)
		{
			case 1: // brass
				if(freq>=2)
				{
					freq-=(int)(freq);
					freq-=2;
					if(freq >= 0.5833)
						freq -= 0.584962500721156181453738943947817f;
					else if(freq >= 0.32)
						freq -= 0.321928094887362347870319429489390f;
					else if(freq >= 0.1666)
						freq -= 0.169925001442312362907477887895633f;
					freq-=1;
				}
				else if(freq>=1)
				{
					freq-=1;
					if(freq >= 0.5833)
						freq -= 0.584962500721156181453738943947817f;
					else if(freq >= 0.32)
						freq -= 0.321928094887362347870319429489390f;
					freq-=1;
				}
				else if(freq>=0)
				{
					if(freq >= 0.5833)
						freq -= 0.584962500721156181453738943947817f;
					freq-=1;
				}
			break;
			case 2: // clarinet
				if(freq >= 3.90689059560851852932405837343721) // 45th harmonic!
					freq -= 5.49185309632967471077779731738502f;
				else if(freq >= 3.1666) // 27th harmonic
					freq -= 4.75488750216346854436121683184345f;
				else if(freq >= 2.66296501272242931233978047895902f) // 19th
					freq -= 4.24792751344358549379351942290683f;
				else if(freq >= 1.5833) // 9th
					freq -= 3.16992500144231236290747788789563f;
				else if(freq >= 0.736965594166206166416580485541574f) // 5th
					freq -= 2.32192809488736234787031942948939f;
				else if(freq >= 0) // 3rd
					freq -= 1.58496250072115618145373894394782f;
			break;
		}
		freq += register_key;
		
		freq = (float)(getSampleRate())*pow(0.5,freq)/440.0;
		
		delay_time_estimate = freq;
		
		freq -= filter_delay;
//		if(filter_parameter!=0 && filter_parameter!=1)
//			freq -= 1/(1/(1-filter_parameter)-1);;
		
/*		if(glide_freq >= 10000000000000000)
			glide_freq = freq;
			
		if(freq > glide_freq)
		{
			glide_freq += glide / (pow(2,12.0-fade_time*10)/getSampleRate());
			if (glide_freq > freq)
			{
				frequency = freq_in;
				glide_freq = freq;
			}
			else
			{
				frequency = -1000000000000000.0f;
				freq = glide_freq;
			}
		}
		else
		{
			glide_freq -= glide / (pow(2,12.0-fade_time*10)/getSampleRate());
			if (glide_freq < freq)
			{
				frequency = freq_in;
				glide_freq = freq;
			}
			else
			{
				frequency = -1000000000000000.0f;
				freq = glide_freq;
			}
		}
*/
		dl = (int)(freq-0.5);
		// CLAMP freq to min & max delay length!
		if (dl<1)
			dl = 1;
		if (dl>buffer_size)
			dl = buffer_size;
		read_pointer[1] = (write_pointer+dl) & mask;
		
		freq -= dl;
		if(freq<=0.1)
			freq = 0.1f;
		if(freq>=21.5)
			freq = 21.5f;

		if(freq == 1)
			allpass_parameter[1] = 0;
		else
			allpass_parameter[1] = 1 / (2/(freq - 1) + 1);

		if(j==1)
		{
			mix = 0;
			mix_delta = 0;
			read_pointer[0] = read_pointer[1];
			allpass_parameter[0] = allpass_parameter[1];
			allpass_memory[0] = allpass_memory[1];
			allpass_memory2[0] = allpass_memory2[1];
		}
		
	}
}

void Pipe::set_pipe_type(short in)
{
	pipe_type = in;
	if(pipe_type == 2 && feedback>0)
		feedback *= -1;
	else if(pipe_type != 2 && feedback<0)
		feedback *= -1;
}

void Pipe::set_feedback(float in)
{
	if(in>=0.9999)
		in = 0.9999f;
	else if(in<=-0.9999)
		in = -0.9999f;
	feedback = in;

	if(pipe_type == 2 && feedback>0)
		feedback *= -1;
	else if(pipe_type != 2 && feedback<0)
		feedback *= -1;
}

void Pipe::set_cutoff(float in)
{

	if(in<0.0001)
		in=0.0001;
	if(in>1)
		in=1;
	filter_parameter = in;
	if(in!=1)
		filter_delay = 2/(1/(1-in)-1);
	else
		filter_delay = 0;
/*	
	frq = pow(2,in*11-5.5)/getSampleRate()*440;
	bandwidth = pow(2,in*24-12)/getSampleRate()*440;
	if (bandwidth>=0.99999f)
		bandwidth = 0.99999f;

	filter_parameter = cos(frq*2*FLOAT_PI) * (1-bandwidth);
	filter_parameter2 = sin(frq*2*FLOAT_PI) * (1-bandwidth);
	filter_gain = bandwidth;
	
	filter_delay = 1/(1/bandwidth-1);
*/
}

void Pipe::set_fade_time(float in)
{
	fade_time = in;
}

void Pipe::set_glide(float in)
{
	if(in<0)
		in=0;
	glide = in;
}

void Pipe::set_register_key(float in)
{
	register_key = in;
}


void Pipe::park()
{
	frequency = -1000000000000000.0f;
	glide_freq = 10000000000000000.0f;
	mix=0;
}
float Pipe::calculate_decay()
{
	float temp;
	temp = fabs(feedback);
	if (temp<0.01)
		temp = 0.01f;
	temp = log(1000000.0f) / -log(temp);
	
	return (temp * delay_time_estimate + 100);
}

inline float Pipe::output_signal()
{
	float pipe0, pipe1, temp;
	float tf1, tf2;
	
	mix -= mix_delta;
	if(mix<0)
		mix=0;
	pipe0 = buffer[read_pointer[0]];
	allpass_memory[0] = pipe0 + allpass_memory[0] * allpass_parameter[0];
	pipe0 = allpass_memory[0];
	temp = pipe0;
	pipe0 = allpass_memory2[0] - pipe0 * allpass_parameter[0];
	allpass_memory2[0] = temp;
	
	pipe1 = buffer[read_pointer[1]];
	allpass_memory[1] = pipe1 + allpass_memory[1] * allpass_parameter[1];
	pipe1 = allpass_memory[1];
	temp = pipe1;
	pipe1 = allpass_memory2[1] - pipe1 * allpass_parameter[1];
	allpass_memory2[1] = temp;
	
	temp = (mix)*pipe0 + (1-mix)*pipe1;
	temp *= feedback;
	filter_memory  += (temp-filter_memory) *filter_parameter;
	filter_memory2 += (filter_memory-filter_memory2)*filter_parameter;
	return ( filter_memory2 );
}

inline void Pipe::input_signal(float input)
{
	buffer[write_pointer] = input;
	read_pointer[0] = (read_pointer[0]-1)&mask;
	read_pointer[1] = (read_pointer[1]-1)&mask;
	write_pointer = (write_pointer-1)&mask;
}

Usage:
pipe.set_frequency(freq); // this will not use too much CPU
float temp = pipe.output_signal();
// Do stuff here
pipe.input_signal(pipe_input);

[hibrasil]

wow, thanks .

Is this synthedit module available somewhere?

oli

[MadBrain]

Well, here's the VST export from synthedit:
http://www.engramstudio.com/uploader/src/model2.zip

WARNING: very alpha test version prototype. Has some weaknesses I have to fix (low range sounds awful, needs better filters, parameter ranges are all over the place).

[xoxos]

are you planning on releasing a freeware derived from this? i don't feel like looking at the code atm, your 'lips' are better than mine i think a more versatile instrument than milestone can be produced from your script.

thanks for sharing regardless, i'm sure i'll learn from it.

i don't know who or where you are in your pm dev, my 'reed5c' (my forum quickies thread) is an *awful* instrument.. it adds a register hole to a simple waveguide model, and while it doesn't play in tune and my techniques are primitive, you can hear from it that the addition of the reg hole adds overblow harmonics and realism for reeds.

i would be happy to share what i have with anyone working on freeware brass and reed candidates.

[Fred e]

Take a look at the SM synth kaplus plus, that one has a Synthmaker code module in it implementing a allpass interpolated delay line.

go to www.synthmaker.com and dig around the forum

[JonasNorberg]

I'm looking for some source code for a delay line / resonant comb filter / karpus strong algo using allpass interpolation rather than linear interpolation. I don't have the dsp knowledge to translate the diagrams etc from the papers into code

Helix is open source check out the string model:

In this file:
http://audjoo-helix.svn.sourceforge.net ... iew=markup

line 1914

The function is called
void OscillatorInstance::processAccumStringModel( Synth* synth, OscillatorSettings* osc, float* controlbuffer, float* audiobuffer, int count )

And you might want to look for "// 4-point, 3rd-order Lagrange (x-form)"

This code is not ready to all copy and paste, but feel free to take parts or read it for reference.

Cheers.

PS. Also have a look at

"Polynomial Interpolators for High-Quality Resampling of Oversampled Audio"
http://www.audjoo.com/docs/Polynomial%2 ... lators.pdf

DS.


Jonas