Sine hard sync

[rigatoni_modular]

Sorry to further hijack the math, I'm happy to start another topic if I should do that, but I'm honestly going insane at a place where I feel like I'm _so_ close to being able to release this project but blep is messing me up. The good news is I have THREE working implementations of it. A 0th order MinBLEP from VCV Rack's library, a 0th and 1st order 2-sample PolyBLEP from Mutable Instruments, and a 0th-3rd order BLEP using residuals based on this thread.

The bad news is that while I can confirm that all of these BLEPs work to anti-alias signals i.e. I think I'm offsetting them correctly at a sub-sample level and they almost completely remove the aliasing from a square and a synced square, they're hit or miss on whether they make a synced sine better or worse. This issue further compounds when using multiple summed sines, where I'm often getting worse aliasing when turning on BLEP.

I was under the (maybe naïve) understanding that antialiasing a synced signal composed of several sinusoids would be as effective as antialiasing a single synced sinusoid. I've confirmed some of the obvious things like approximate magnitude of the discontinuity being correct. Are there any known potential pitfalls with this approach or should I accept defeat?

Below I've included images of the spectrum and oscilloscope views of a few scenarios.

3500Hz square synced to 2330Hz square, BLEP off:

2330_3500_square_off.png

3500Hz square synced to 2330Hz square, BLEP on:

2330_3500_square_on.png

3500Hz sine synced to 2330Hz square, BLEP off:

2330_3500_sin_off.png

3500Hz sine synced to 2330Hz square, BLEP on:

2330_3500_sin_on.png

3500Hz 8-sine approximation of square synced to 2330Hz square, BLEP off:

2330_3500_square_additive_off.png

3500Hz 8-sine approximation of square synced to 2330Hz square, BLEP on:

2330_3500_square_additive_on.png

[Z1202]

Not sure what your exact problems are - don't have time to get into details ATM, but maybe you're running into the following. As sine frequency gets closer to Nyquist, adding further higher-order BLEPs becomes progressively more critical. Theoretically you need to add infinitely many BLEPs for a synced sine, it's just that at lower frequencies the higher-order corrections become small and could be ignored to an extent. This is different from sawtooth and square where only 0-th order BLEP is needed, or triangle, where additionally only 1st order BLEP is needed, no matter what the frequency. Some further ideas are shared here https://www.native-instruments.com/file ... neSync.pdf

[Z1202]

@mystran did you really compare the actual spectra of integrated windowed sync vs. windowed residual, or were your complaints just based on the derivative discontinuity? As I wrote earlier, for higher-order bleps I could imagine standard windows maybe being not the best fit, but I'd guess for 0th order Hann window should work pretty well.

[mystran]

@mystran did you really compare the actual spectra of integrated windowed sync vs. windowed residual, or were your complaints just based on the derivative discontinuity?

For the desmos plot? No.

I did however try this a few years ago and all I was able to get was significant (several dB) ripple in response, might have been with higher orders, can't remember exactly. It was around the same time I was messing with Lagrange-interpolators and figured out the maximum order you can get out of a BLEP seems related to how flat it is around DC.

With windowed sinc integrated using trapezoidal (first-half) and using symmetries, I can have "perfect" results up to piecewise cubic which has so far been sufficient for my purposes.

ps. I would also like to add that I have not spent too much time on trying to hard-sync sines, because I feel like that problem is somewhat academic (ie. having a solution that actually gives the full result would be interesting, but just brute-forcing a bunch of derivatives until it's "good enough" is not terribly interesting to me). You can hard-sync them perfectly if you approximate them as piece-wise polynomial. You can also theoretically oversample by 2x (to allow for bandwidth expansion) and window each sine-segment with a rectangular window anti-alised by simple BLEP-steps... though this gets a bit expensive obviously if a lot of those windows (including the BLEP part) overlap.

[Z1202]

having a solution that actually gives the full result would be interesting

The paper I linked has a "full solution", although I forgot all the details, since I haven't been using it ever since. What I remember though is that it involves one of exponential integrals, which (lacking a realtime function to evaluate it) needs to be tabulated, and the required table length grows with the sine frequency (not unlike the number of BLEP orders grows), as the function needs to be evaluated at larger and larger argument values, IIRC.

[rigatoni_modular]

Not sure what your exact problems are - don't have time to get into details ATM, but maybe you're running into the following. As sine frequency gets closer to Nyquist, adding further higher-order BLEPs becomes progressively more critical. Theoretically you need to add infinitely many BLEPs for a synced sine, it's just that at lower frequencies the higher-order corrections become small and could be ignored to an extent. This is different from sawtooth and square where only 0-th order BLEP is needed, or triangle, where additionally only 1st order BLEP is needed, no matter what the frequency. Some further ideas are shared here https://www.native-instruments.com/file ... neSync.pdf

I'm following the native instruments paper for my implementation. I guess I would've just expected the first blep to at least reduce the aliasing by ~6db/oct as advertised instead of making the aliasing worse...

My sine is currently a 9th order Chebyshev approximation stolen from: https://web.archive.org/web/20200628195 ... oximation/

Maybe if I used a piecewise cubic sine approximation I'd have better results?

[juha_p]

Does that matter? Oh.. it does.. https://www.desmos.com/calculator/ospvpzbijq

Well.. it's closer now, but still not quite a windowed sinc... and in fact we see that if we do subtract a windowed sinc, that there appears to be a discontinuity in first derivative at zero.

Legal to manipulate w(x) equation (by using bigger divider in cos() you'll get h(x) = w(x)k(x) better match)?

[mystran]

Does that matter? Oh.. it does.. https://www.desmos.com/calculator/ospvpzbijq

Well.. it's closer now, but still not quite a windowed sinc... and in fact we see that if we do subtract a windowed sinc, that there appears to be a discontinuity in first derivative at zero.

Legal to manipulate w(x) equation (by using bigger divider in cos() you'll get h(x) = w(x)k(x) better match)?

Yes. The error goes down when the window gets longer. I intentionally left the window short to emphasize.

[mystran]

I'm following the native instruments paper for my implementation. I guess I would've just expected the first blep to at least reduce the aliasing by ~6db/oct as advertised instead of making the aliasing worse...

BLEPs should never make things worse if implemented correctly. Every additional derivative you deal with should result in at least some improvement. They are quite finicky though 'cos you're essentially trying to cancel out aliasing by adding an inverted copy and if something is slightly off, the result is quite often worse than what you started with.

So.. don't panic, the method isn't fundamentally broken, most likely you just have some subtle error somewhere... and frankly I highly doubt anyone gets this stuff right on the first try.

[rigatoni_modular]

So.. don't panic, the method isn't fundamentally broken, most likely you just have some subtle error somewhere... and frankly I highly doubt anyone gets this stuff right on the first try.

Thanks for the encouragement I'm going to log all the blepped samples I create so I can graph them and see if things look smooth. The one thing that's confusing me is that I'm using the exact same sub-sample offset for my square and for my sine so unless my scaling is wrong or something I'm at a loss. Hopefully the logging will help.

[rigatoni_modular]

So.. don't panic, the method isn't fundamentally broken, most likely you just have some subtle error somewhere... and frankly I highly doubt anyone gets this stuff right on the first try.

Following up here, I'm finally in a place where I can properly debug this and would love a sanity check. I've measured the raw synced sine values as well as my anti-aliased values and the residual I'm adding. For each sync I've also measured:

- The discontinuity magnitude we're correcting for (the value the sine would be at this sample without sync subtracted from the value it's at after the sync)
- The sub-sample offset into the residual table (for this I'm using the estimated # of samples (<1) ago that the sync happened which seems correct).

I'm still seeing some discontinuities in the output and I was curious if this was expected or not. I would naïvely expect the difference between the residual value 1 sample before the discontinuity and at the discontinuity to be the size of the discontinuity to make it disappear completely. That's not what I'm seeing, which makes some sense given that the interpolated residual values will never be the extrama on either side of the ideal residual.

I'm pretty sure that my sub-sample offset is correct so is it possible that I have a magnitude or interpolation issue or does this look about right?

AntiAliasingView1.png

[mystran]

It's quite hard to say from this picture what might be wrong.. but you can make things easier for yourself, by slowing everything down: make the frequency lower (so longer period between syncs), make the BLEP cutoff lower (by a few octaves, so it becomes easier to see the actual shape) and then make it longer (to compensate for the lower frequency).

With a lower BLEP cutoff, you also gain another advantage: now you can see what happens to the spectrum above the cutoff without having to infer this from aliasing.

[rigatoni_modular]

It's quite hard to say from this picture what might be wrong.. but you can make things easier for yourself, by slowing everything down: make the frequency lower (so longer period between syncs), make the BLEP cutoff lower (by a few octaves, so it becomes easier to see the actual shape) and then make it longer (to compensate for the lower frequency).

With a lower BLEP cutoff, you also gain another advantage: now you can see what happens to the spectrum above the cutoff without having to infer this from aliasing.

Great idea! I've been doing these tests at higher frequencies because it's hard to reproduce the weird stuff at lower frequencies so I'll try your suggestion next. Me being a DSP noob, please forgive this question: is lowering the "ringing frequency" of the sinc I'm integrating to get the 0th order BLEP analogous to lowering the BLEP cutoff?

[2DaT]

Following up here, I'm finally in a place where I can properly debug this and would love a sanity check. I've measured the raw synced sine values as well as my anti-aliased values and the residual I'm adding. For each sync I've also measured:

- The discontinuity magnitude we're correcting for (the value the sine would be at this sample without sync subtracted from the value it's at after the sync)
- The sub-sample offset into the residual table (for this I'm using the estimated # of samples (<1) ago that the sync happened which seems correct).

First of all, it's recommended to make a basic sawtooth oscillator first to make sure everything works as expected. Then add hard sync.

There are 2 conventions one can use for computing fractions: you can use (0,1], where 0 means 1 sample behind and 1 means it's on current sample. But I find it to be very inconvenient, because it makes blep table lookup harder (need to (1-frac), because in blep table moving right means that we delay the transition so it's the other way around).

I prefer (1,0] where 0 means it's on current sample and 1 means it's 1 sample behind. It's also simpler to calculate "inverted" fraction when doing a phasor reset:

Code: Select all

phase += delta;
if ( phase >= 1.0f )
{
	phase -= 1.0f;
	fraction = phase / delta;
}

Then for hard-sync things get more complicated, because we can get either of two: a phase reset then a hard sync reset, or hard sync reset alone, so extra care must be taken to determine phase at hard sync precisely. Also not to insert a phase reset transition if it happened after hard sync.

Code: Select all

if ( phase >= 1.0f )
{
	phase -= 1.0f;
	fraction = phase / delta;
	if ( sync_frac >= fraction )
	{
		phase += 1.0f;
	}
}
float phase_after_sync = phase;
float phase_at_sync;
if ( sync_frac >= 0 )
{
	phase_after_sync = sync_frac * delta;
	phase_at_sync = phase - phase_after_sync;
}

That's more or less everything you need to compute hard-synced saw, sine and triangle. Pulse waves with pwm on the other hand are much more complicated, because pwm signal must be interpolated, and there can be multiple transitions per sample, as well as transitions after hard sync which must be handled properly. Also, under extreme pwm modulation, pulse can back swing (i.e go from triggered back to untriggered, but without phase going through 1), which needs to be handled. Although "backswing" isn't extremely noticeable until audio rate pwm, so it isn't a very crucial thing to do. But it can happen.

A tip: phase-pw, and same difference one sample earlier is enough to calculate everything without using extra state.

[rigatoni_modular]

Code: Select all

if ( phase >= 1.0f )
{
	phase -= 1.0f;
	fraction = phase / delta;
	if ( sync_frac >= fraction )
	{
		phase += 1.0f;
	}
}
float phase_after_sync = phase;
float phase_at_sync;
if ( sync_frac >= 0 )
{
	phase_after_sync = sync_frac * delta;
	phase_at_sync = phase - phase_after_sync;
}

I’m using the 2nd sync offset scheme you mentioned fortunately. A little confused about the 2nd code sample though - shouldn’t phase_at_sync just be the old phase plus delta? The subtraction there is puzzling.