minBLEPS once and for all

[d_m_chambers_]

Hi guys,

I've tried to implement anti aliasing on a simple Saw, using minBLEP. I wonder if somebody could correct my code because the current output aliases like crazy and has a very strange spectrum, as shown below



So, in my oscillator's constructor I have the following, using Daniel's 'GenerateMinBLEP' with 24 zero crossings and 32 oversampling, giving minBLEP data 1537 in length with a circular buffer 48 in length :

Code: Select all

        circularBuffer = new float[48];
        for (int i = 0; i < 48; i++)
        {
            circularBuffer[i] = 0.0f;
        }
        minBlep = GenerateMinBLEP(24, 32);

Plotted, this looks right:



nextSample looks like this:

Code: Select all

#define LERP(A,B,F) ((B-A)*F+A)
double MinBlepOsc::nextSample() {
    
    mPhase += mPhaseIncrement;
    index = (index+1)%48;
    while(mPhase >= 1.0)
    {
        mPhase -= 1.0;
        double exactCrossTime = 1.0-((mPhaseIncrement-mPhase)/mPhaseIncrement);
        for(int i = 0; i < 47; i++)
        {
            double tempIndex = (exactCrossTime*32.0)+(i*32.0);
            double tempFraction = tempIndex-floor(tempIndex);
            circularBuffer[(index+i)%48] += (1.0-LERP(tempFraction, minBlep[(int)floor(tempIndex)], minBlep[(int)ceil(tempIndex)]));
        }
    }
    
    circularBuffer[index] += (mPhase-mPhaseIncrement);
    double output = 1.0*circularBuffer[index];
    circularBuffer[index] = 0.0;
    return output;

    
}

Clearly I'm doing something very wrong! It would be great if somebody could correct my code/ point me in the right direction.

Cheers, Dave

[stratum]

If you plot the result of that code I guess it would be easy to figure out what's wrong (given that the minBLEP is correct, that's what I would expect)

[d_m_chambers_]

Hi stratum. Thanks for the tip. I'll give that a try shortly..

[d_m_chambers_]

If you plot the result of that code I guess it would be easy to figure out what's wrong (given that the minBLEP is correct, that's what I would expect)

Well, when I looked at the waveform I saw that it's 'offset' in the y direction :



...changing that however didn't improve the weird spectrum and hideous aliasing but I did notice something rather strange, zooming in, you'll see that the minBLEP seems to be applied to the BOTTOM of the wave (y is zero) not the top (y is 1) where I think it should be applied:



This seems wrong and the minBLEP just looks wrong. I'm hoping somebody reading can correct the few lines of code written above. It must be that 48, 32 and 24 just don't work

[stratum]

Apparently there is some working code here http://www.musicdsp.org/showone.php?id=112

As far as I can see all that this does is to add the "blep" waveform at zero crossings, and its "amplitude" (sign actually, as it is either 1 or -1) changes according to the direction of signal during the zero crossing (upward or downward). The blep is mixed to the "naive" waveform but the mix control value (the amount of each waveform, "naive" vs "blep") is a linearly decaying fraction:

from the above link (MinBLEPS.zip / osc.c):

Code: Select all

// add impulse into buffer
void osc_AddBLEP(osc_t *lpO, double offset, double amp)
{
int i;
double *lpOut=lpO->buffer+lpO->iBuffer;
double *lpIn=gMinBLEP.lpTable+(int)(KTABLE*offset);
double frac=fmod(KTABLE*offset,1.0);
int cBLEP=(gMinBLEP.c/KTABLE)-1;
double *lpBufferEnd=lpO->buffer+lpO->cBuffer;
double f;

	// add
	for (i=0; i<lpO->nInit; i++, lpIn+=KTABLE, lpOut++)
	{
		if (lpOut>=lpBufferEnd) lpOut=lpO->buffer;
		f=LERP(lpIn[0],lpIn[1],frac);
		*lpOut+=amp*(1-f);
	}

	// copy
	for (; i<cBLEP; i++, lpIn+=KTABLE, lpOut++)
	{
		if (lpOut>=lpBufferEnd) lpOut=lpO->buffer;
		f=LERP(lpIn[0],lpIn[1],frac);
		*lpOut=amp*(1-f);
	}

	lpO->nInit=cBLEP;
}

The idea looks simple, it may be tricky to get right to cover all possible cases.

And this is how zero crossings are handled for square and saw (from the same file, osc.c):

Code: Select all

// play waveform
double osc_Play(osc_t *lpO)
{
double v;
double fs=lpO->f/44100;

	// create waveform
	lpO->p=lpO->p+fs;

	// add BLEP at end of waveform
	if (lpO->p>=1)
	{
		lpO->p=lpO->p-1.0;
		lpO->v=0.0f;
		osc_AddBLEP(lpO, lpO->p/fs,1.0f);
	}

	// add BLEP in middle of wavefor for squarewave
	if (!lpO->v && lpO->p>lpO->fPWM && lpO->type==OT_SQUARE)
	{
		lpO->v=1.0f;
		osc_AddBLEP(lpO, (lpO->p-lpO->fPWM)/fs,-1.0f);
	}

	// sample value
	if (lpO->type==OT_SAW)
		v=lpO->p;
	else
		v=lpO->v;

	// add BLEP buffer contents
	if (lpO->nInit)
	{
		v+=lpO->buffer[lpO->iBuffer];
		lpO->nInit--;
		if (++lpO->iBuffer>=lpO->cBuffer) lpO->iBuffer=0;
	}

	return v;
}

The blep table is loaded from a matlab matrix file in function minBLEP_Init, and you can replace this with your own blep generator:

Code: Select all

BOOL minBLEP_Init()
{
	// load table
	FILE *fp=fopen("minblep.mat","rb");
	unsigned int iSize;

	if (!fp) return FALSE;

	fseek(fp,0x134,SEEK_SET);

	fread(&iSize,sizeof(int),1,fp);
	gMinBLEP.c=iSize/sizeof(double);

	gMinBLEP.lpTable=(double*)malloc(iSize);
	if (!gMinBLEP.lpTable) return FALSE;

	fread(gMinBLEP.lpTable,iSize,1,fp);

	fclose(fp);

	return TRUE;
}

edit 1: From the replies in http://www.musicdsp.org/showone.php?id=112 it seems like this code is not entirely correct, especially for the square wave, but nevertheless it looks better than a text-only description.

edit 2: KTABLE is oversampling ratio of the blep table (i.e. blep sampling rate / oscillator sampling rate) and gMinBLEP.c is the number of samples in the blep table. osc_AddBLEP simply reads values corresponding to the current oscillator position from the table without doing any operation other than simply decaying the mix ratio (the reason is that the wave in the table already properly band-limited.)

edit 3: Hard sync does not seem to be implemented, probably would require some modification to osc_Play (i.e. handle zero crossings due to an extra oscillator) Apparently this can complicate matters enough to require modification of osc_AddBLEP too (I'm not sure about that, since there is an extra mixing step in osc_Play, osc_AddBLEP could possibly stay as it is, perhaps that extra mixing step would just be replicated for the second oscillator) https://www.keithmcmillen.com/blog/simp ... ator-sync/

[d_m_chambers_]

Thanks stratum, I'll give that a try soon and post if I had any luck.

[kryptonaut]

Dave (d_m_chambers) - in your code, there is a discrepancy between the definition and usage of LERP. You define it with the fraction as the third parameter, but call it with the fraction first.
Try:

Code: Select all

#define LERP(F,A,B) ((B-A)*F+A)

Apart from that the BLEP stuff looks ok to me.

[d_m_chambers_]

kryptonaut, I can't thank you enough. You're certainly eagle-eyed. That was the issue. I now get a clean, alias free signal.

Thanks so much!!

[kryptonaut]

Great, glad I could help. Hope the rest of your project goes well - MinBLEPs are nice when they work, but the details can be tricky to get right!

[stratum]

but the details can be tricky to get right!

This picture shows why
https://www.keithmcmillen.com/wp-conten ... 6/sync.png

[stratum]

Actually looking at that picture, seems like all that needs to be done is to insert another blep with a different amplitude at that hard sync point since just because we have a step with different amplitude does not mean that its width should also be different. I'd just mix them all.

[kryptonaut]

Yes for hard sync you need to add in another blep, scaled according to the size of the step it's replacing. Care needs to be taken when the slave reset happens in the same single-sample timestep as the master reset.

[stratum]

Care needs to be taken when the slave reset happens in the same single-sample timestep as the master reset.

Thanks for the tip.

[d_m_chambers_]

Thanks guys but I worked out my problem so there is no need for my original post in this space any more. Happy coding!

[d_m_chambers_]

ignore - I accidentally replied to my own post