jpcima / bbd-delay-experimental

Testing grounds for newer BBD delay implementation
Boost Software License 1.0
11 stars 0 forks source link

Making the BBD implementation correct #1

Open jpcima opened 5 years ago

jpcima commented 5 years ago

@b3ll, I'll describe in some detail what I'm doing here, and I let you experiment what may work or not. Considering the paper is not really clear at all times, I've taken a few guesses.

If the works reaches a point of apparent stability, I'll merge it. For now, I've added 2 commits in the unit-delay branch.

Filter couldn't operate according to the structure shown, because the unit delay variable was not kept, meaning either algorithm or structure has been incorrect. Most probably algorithm, so I reformulated the filter with use of unit delays (as Xout_mem).

Capture du 2019-07-20 20-05-18

You have the processing made of 2 loops: (as denoted in the paper)

n-loop implements as: https://github.com/jpcima/bbd-delay-experimental/blob/543cd7ea6e4452ce421451f07262a6bc4eb70889/bbd_line.cc#L72

Further, this n-loop does two different jobs on even/odd ticks:

You will observe odd-case will perform this add: https://github.com/jpcima/bbd-delay-experimental/blob/543cd7ea6e4452ce421451f07262a6bc4eb70889/bbd_line.cc#L89

I have found suspicious that BBD's output Xout can be allowed to grow unlimited, as this sum is never reset.

From the paper, I was able to read this information. Capture du 2019-07-20 20-19-14

If I correctly interpret "over one interval of the k side", it suggest to reset the Xout sum at each frame. Hence I did. https://github.com/jpcima/bbd-delay-experimental/blob/e7f197e8c26f3c60f963006d955b84dd82d7844e/bbd_line.cc#L70-L71

Thinking of it, it means at some frames, accumulation may remain 0 if BBD doesn't tick (the case of a slow clock). But also intuitively, it makes sense: in analog world, when BBD doesn't tick, the bucket will empty itself (capacitor discharge), and it's the similar idea when zero are fed into the lowpass following.


Hope this little comment makes sense. I didn't try this action, you can tell me if this has worked or not.

b3ll commented 5 years ago

Ok, so I tried both of these, both together and individually and I didn't notice anything different, unfortunately. These things do seem like they're correct regardless (especially zeroing out the ever increasing output filter, and the emptying bucket when there's no clock ticks).

However, I did some more analysis and it looks like the issue is even worse when the delay time is being modulated.

Assuming that what I'm doing is correct (this should be a single delay line, with step size of 256, modulating a delay value from 0.00166 to 0.00535 over a period of 1hz with a simple sine wave:

  const float modFrequency = 1.0f;
  const float minDelay = 0.00166f;
  const float maxDelay = 0.00535f;

  float modDelayBase = (maxDelay + minDelay) / 2.0f;
  float modDelayDepth = (maxDelay - minDelay) / 2.0f;
  bucketSize = 256;
  auto *leftChannel = buffer.getWritePointer(0);
  float leftClock[numSamples];
  for (int i = 0; i < numSamples; i++) {
    float modDelayValue = sin(2 * M_PI * _leftTime);

    float delayTime = modDelayBase + (modDelayDepth * modDelayValue);

    leftClock[i] = ((bucketSize * 2.0) / delayTime) / sampleRate;

    _leftTime += (modFrequency / sampleRate);
    _leftTime -= (int)_leftTime;
  }
  _leftDelayLine->process(numSamples, leftChannel, leftClock);

It ends up still shelving a lot of artifacts in the top end of the spectrum (regardless of what note is being played). I attached a sound file with before / after, when you hear the delay being enabled you'll hear a buzzing noise (alongside the delayed signal). There's also a video showing the filter output after playing without the delay enabled, playing the same note again without the delay enabled, and then the delay being enabled. You'll see the top end is oscillating.

Part of me wants to think that somewhere along (perhaps the start of) the BBD implementation, some noise gets introduced into the delay and then that noise is persisted indefinitely. When I initially read this I thought that the Xout growing unlimited could've been the issue, but I don't think that's the cause. The buzzing is 100% related to the modulation since if I set a slower or faster modulation rate, the buzzing changes to match.

filter-output.zip

jpcima commented 5 years ago

For now I don't know, I'll just note a couple things.

b3ll commented 5 years ago

Yeah, I was about to try oversampling to see what happens. I imagine changing the sample rate to 88200hz isn't oversampling right? To oversample, I should be taking the input, and interpolating the input to 88200hz (by doubling and averaging, right?)

i.e. [sample1, sample2, sample3, …] becomes [sample1, (sample1+sample2) / 2, sample2, (sample2+sample3) / 2, …]

Also, in the technical writeup they talk about having the ƒBBD at 50000hz (for their test sources), how would I go about doing that in your implementation? Would I be supplying the delay as 0.000020s without modulation?

jpcima commented 5 years ago

Yeah, I was about to try oversampling to see what happens. I imagine changing the sample rate to 88200hz isn't oversampling right? To oversample, I should be taking the input, and interpolating the input to 88200hz (by doubling and averaging, right?)

No, it's a little more elaborate. More like:

  1. upsample input: fill gaps with zero samples, and smooth it with lowpass
  2. setup the bbd element with N× times scaled frequencies (filters rate divided by N, BBD clocks divided by N also)
  3. run then lowpass output
  4. decimate by N (keep every Nth output)

Anyway, better to fix the ordinary case first, before doing the oversampling.

Also, in the technical writeup they talk about having the ƒBBD at 50000hz (for their test sources), how would I go about doing that in your implementation?

It's a simple matter of passing (50000/Fs) as clock value. It handles the frequencies over Nyquist, it's not a problem to clock a line to be this high.

b3ll commented 5 years ago

Ok, now we're getting somewhere!

I setup the bbd as such:

  // setup
  _leftDelayLine = new BBD_Line(44100.0, 256, bbd_fin_j60, bbd_fout_j60);

  // inside process block
  auto *leftChannel = buffer.getWritePointer(0);
  float leftClock[numSamples];
  for (int i = 0; i < numSamples; i++) {
    leftClock[i] = 50000.0 / 44100.0;
  }
  _leftDelayLine->process(numSamples, leftChannel, leftClock);

  buffer.copyFrom(1, 0, buffer.getReadPointer(0), buffer.getNumSamples());

It turns out that when running this with the sample wav files supplied (they're attached), there's no sort of modulating / delaying effect heard at all. It just sounds super noisy and has lots of static added (which is probably what that static noise is when the BBD is modulated.

I looked at the sample files, and it looks like we're not making the same shape at all

Here's a photo of the original C chord input file:

orig chord no delay

Then this is what the research paper's output file looked like:

sample chord with delay

But this is what ours looks like (using what I wrote above):

bbd impl with delay

I'm currently trying to understand the whole algorithm to see if there's something missing at the beginning or if there's an error in how the signal is being delayed.

Here's a zip of the input chord, their output model chord, and what we're generating (as a before / after). chords.zip

jpcima commented 5 years ago

My next step in the plan is trying to implement "3.3 Real-valued systems".

The goal is to rewrite the complex-valued filter in a real-valued one. It takes a modification of the filter structure. (eg. the case M=5, it converts to 1 structure of first order and 2 structures of second order)

When it's real-valued, it's a return to a more conventional filter structure, and it makes the analysis easier.

Also now, the branch contains a test program which works on audio files. The buzzing distorsion is reproducible, even on the simple case, a sine wave at 1kHz.

b3ll commented 5 years ago

Okay, that sounds great. I tried of altering a lot of things last night, but I couldn’t seem to get the buzzing to go away… oversampling did help but only with absurdly high sample rates (i.e. 192000)

I was trying to see how the filter was done, but I genuinely don’t understand how the complex values map to a filter cutoff point. Filters I’ve normally written are usually something like (anything greater than some frequency is dropped, and others are dropped by a certain percent, depending on the curve).

Is there anything I can read to understand how the complex values / filter works? I’d very much like to help make the real value one

jpcima commented 5 years ago

You haven't picked the simplest of exercises for a start. it's for sure :D

This filter discretization is very unusual, but I think it exploits the fact of BBD being a sampling system. In this logic, it can treat signal input as a series of instantaneous pulses which occurs at the BBD ticks.

As I understand, it's because of this special case, that the filter can be simplified and expressed as such.

So yeah, considering it's so application-specific, you will have surely a hard time to get in the detail of this.

But if you search a generic explanation of filter theory, probably Miller Puckette's book has the best lesson I have found. (The Theory and Technique of Electronic Music)

b3ll commented 5 years ago

hah, yeah, it's fun to learn piece by piece though!

Upon checking out the book it has a section on delays which deals with complex numbers and it seems to share a lot of the same stuff as mentioned in the research paper, so I'll read that stuff too!

jpcima commented 5 years ago

Adam, I think the error is finally fixed at unit-delay, in the commit aadf660. As it seems, the output filter was just generated with a wrong formula. Please check with the new one.

b3ll commented 5 years ago

IT WORKS! 💯🔥🎹✨

And it sounds great!! Very nice catch! :D :D

b3ll commented 5 years ago

Ok I'm gonna get this working in the chorus effect plugin I'm working on and I'll post the results here :D

jpcima commented 5 years ago

Also, regarding the two fixes made regarding the algorithm.

(1) the reset of the sum must absolutely be made, or signal will go excentered around the DC, which is not the case in any of the website's output recordings

(2) about the filter restructuring and its unit delay node: it seems signal is outputting with more edge before the fix, and as more dull after applying. as it woulds seem, the website recordings appear to exhibit this tiny high frequency content.

So, who knows if it's algorithm which is right, or it's the structure pictured... It's testable with or without patch, and see what will happen...

b3ll commented 5 years ago

Is it possible to mix the filter output? Perhaps its too harsh and clipping more than it should…

Here's ours (with the filter using 1024 as an interp_size)

image

vs. theirs (theirs has additional modulation that I'm not seeing in our delay though so that's kinda weird)

image

These look rather similar and I don't think it's clipping too much, it still sounds great and significantly better than the buzzing beforehand

b3ll commented 5 years ago

In addition, I ran the test program against the input_chord.wav and compared it to the model_chord.wav and they look nearly identical. The model_chord.wav does have some oscillating frequencies in the very top end (around 20k), but that's the only difference I'm seeing… the article only mentions a clock of 50kHz and nothing about modulation, so I'm not sure about that

input_chord.wav ➡️ their file input_model.wav ➡️ input_chord.wav ran through bbd with 256 steps at 50kHz model_chord.wav ➡️ their file

chords2.zip

jpcima commented 5 years ago

It's quite possible that "theirs" has just misimplemented the output filter, and it's explaining the high frequency content. It's as I said above.

Did you compare with 971a7a8 modifications reverted ?

b3ll commented 5 years ago

Yeah, reverting 971a7a8 introduces the buzzing issue and isn't even close

image
jpcima commented 5 years ago

Yep, It's the same observation here. By the algorithm's way of implementing, output filter is inoperational.

Is it possible to mix the filter output?

Do you mean, to have a wet/dry control? you can't make this on chorus side?

b3ll commented 5 years ago

Oh sorry, I meant asking that in relation to the "dullness" of the filter. I think it sounds good but I was proposing the idea of having the output filter having a wet/dry control to add some of that harshness / edge back if desired

jpcima commented 5 years ago

Well it was just the implementation error, which makes the filter disabled. :clown_face:

The thing about these filters, it's they're designed from a fixed specification, it's not directly controllable. As you can see by BBD_Filter_Spec, and values are taken from article's "Table 1".

It's description of the "decomposed" filter form in M sections that the document's I/O filter sections are about. They must relate to filter's analog zeros and poles somewhat, but it's not told exactly how, or how are these values obtained from the Juno.

b3ll commented 5 years ago

LOL fair.

Yeah that part was definitely confusing to me, I had no idea how those values were computed (nor was it mentioned in the research paper).

Regardless, I'm happy with how this sounds! :D

b3ll commented 5 years ago

Here's some more good news, this is my plugin which aims to repro setting number II on the chorus unit vs. a Roland System 8 using its chorus II setting which uses hardware modeling to be nearly identical to the real unit:

Both are playing the note C3 endlessly, that's pretty close

image image
b3ll commented 5 years ago

Also here's a C Major Chord

Mine:

image

Theirs:

image
jpcima commented 5 years ago

I've tried plotting these I/O filters in numpy according to their R and P values. These describe analog transfer functions, R stands for residue, and P for pole.

But the bode plots came out as inverted, I don't know what's up with that. Here it is: filter.py.gz

b3ll commented 5 years ago

uh wait what, that's really weird… if this is inverted, that would explain the tail of the frequency spectrum getting magnified a lot, no?

jpcima commented 5 years ago

No, it means I don't have right interpretation of these values, or I don't know to use scipy software correctly :clown_face:

Anyway now's time for bed, and I don't have brainpower to check this

b3ll commented 5 years ago

I don't understand yet what the values mean, but from the docs it looks like you're using scipy correctly :P

Thank you for all your help though! I really appreciate it

jpcima commented 5 years ago

So, I've returned to this matter and figured out the filter. It's formulation known as Partial Fraction Expansion expressed as:

H(s) = sum for i=[1:M] : (Ri / (s - Pi)) having s=jω=j2πf

and added in the last commit. By plotting the response, I found a minor glitch which occurred on the input filter in the neighbourhood of 100kHz, and I fixed the coefficient.

Transfer function program: https://gist.github.com/jpcima/38c464509213620841385a3a09f96547

And a logarithmic plot: Capture du 2019-08-01 18-25-08

b3ll commented 5 years ago

oh that's awesome! That's super cool, I'll def read up on that… you've gotta see the plugin I've built with this, it sounds really great!

jpcima commented 5 years ago

I post here the Octave script creating custom filters. It's a 5th-order Butterworth 8kHz for this case.

Probably, the Juno one didn't use Butterworth, it used a steeper model for its cutoff. The TCA-350-Y chip from the Solina is known to have Butterworth output from its datasheet.

For other filter design kinds: https://fr.mathworks.com/help/signal/ug/comparison-of-analog-iir-lowpass-filters.html

The script output is R and P arrays, with M being filter's order.

format longG;

cutoff=8000;

[b,a]=butter(5,2*pi*cutoff,'low','s');
%[h,w]=freqs(b,a);
w=logspace(1,5);
h=freqs(b,a,w);
dB = mag2db(abs(h));

plot(w/(2*pi),dB);
xlabel('Frequency (Hz)');
ylabel('Magnitude (dB)');
grid on;
grid minor;

[r,p,k]=residue(b,a)
b3ll commented 5 years ago

Oh this is awesome, yeah the Juno had a -24/db lowpass filter, which seems right. Happy to finally understand where those magic numbers in the table came from hah