Correct indexing bugs due to incorrect refactor in robustIteration()

[AlexHarker]

Fixes the indexing bug discussed today with @weefuzzy.

It could probably do with a quick check over to see that this looks similar to the code here:

https://github.com/flucoma/flucoma-core/blob/2011fc0b21d4ec597d973bc64657542768269452/include/algorithms/ARModel.hpp

Which I think is the last good version.

[AlexHarker]

Impossible for me to tell if the tests should fail given the change or not. The tests that fail are all to do with checking the output of transient processes, but if they are based on old results then it would make sense for them to now be different.

[weefuzzy]

Yes, those are all unsurprising test failures – we'll have to quite literally play it by ear.

[weefuzzy]

This feels qualitatively easier to use with the real-time slicer and extractor.

[AlexHarker]

Thanks - I had no doubt it would be better as all my original audio testing was done on code that was functional. The question for me is how far the improvements go. Is it still missing obvious transients (as per your tests sent separately last month with graphical accompaniment)?

[weefuzzy]

Here's the windowed forward and backward error, plotted against their source material

[AlexHarker]

To my mind that's vastly improved from the previous plots and for good reason as most of the previous ARModel would have been almost entirely based on the raw data, which means that it is likely that if there's enough transient in the window, or it's a bit softer that it likely won't detect (as it did not before).

Levels of noisiness are about what I'd expect - I presume they are about the same on both (it's just that the spikes are higher for backwards), but it's a bit hard to tell due to the scaling differences on display. but it seems to be mostly just under 1, which for a normalised error (assuming the factor is around 3) would be correct/expected.

I see nothing here that points to anything other than merging this and then considering other improvements either to the variance estimate or possibly the manner of prediction (although I still have reservations about whether that route should be taken, or at least how to validate that it is really an improvement from either an empirical or theoretical standpoint) as future enhancements. As far as I can tell this represents what I intended the code to do and performs better so I would vote merge and then we can look at the block edge effects discussed and then leave you to make whatever further changes you wish after that.

It might be a good time to update the tests and pull in the updated autocorrelation code and possibly way of calculating the AR parameters here as there would then be a good basis for checking those changes. That can also be done without me, but I'd be happy to be involved in reviewing that as it removes some other non-optimal aspects of the code to date. Likewise the safety changes to ARModel.hpp

[weefuzzy]

I've adjusted the slicer tests and will run them by @tremblap in supercollider form before fixing up here.

Bafflingly, however, the transients decomposition objects no longer null sum when using the reported latency. I'm going to go and sit in a dark cave to try and figure out why this single change should / could do that (yes, I've been back to the dev branch and rebuilt / rechecked 😉 )

(less sure of this now...)

Well, it 'nulls' to about -100dB and a mean around -200 dB, whereas the existing test in SC has a threshold of -120 dB. Personally, I'm not hugely inclined to lose sleep over that, if it cancels to better than 16-bits.

[AlexHarker]

Hmmmm. If this were an FFT process I'd agree and consider that acceptable because there are loads of points of accuracy to consider. However, given that the output is literally formed from time domain subtraction then the question is what bit depth is the null sum test performed at - I'm guessing given that it's SC we are down to 32 bit? (my thinking is that even if the interpolated output was formed of generated noise, it should still null sum with the input given how the transient output is formed). On the other hand, floating point dB limits are complex to reason about (because floating and also potential catastrophic subtraction relative error). -100dB seems a bit high for me, even for a single sample, but I'd have to look at the whole of the processing chain in the test to have any decent thoughts about what acceptable performance would be. It would also be worth having a relative error test to indicate what that -100dB equates to (more like how many ULP is that?).

In theory 32 bit is obviously 24bits of precision, so you'd expect a worst case 1ULP error in relation to the loudest value to be -144dB - however, this assumes that all signals are under 0dB - if that assumption is not the case then things change, but you'd have to be quite a few dB up and I see no reason that would be the case in the test unless there's some funky gain somewhere.

[weefuzzy]

Well, yes, being SC there's quite a lot of gubbins in the the way between language, server buffers, process, and back again.

Probably worth checking this as close to source as possible before getting too excited.

---

ULPs are new things for me, but with some hasty googling, the answer would be appear to be 'many ULPs':

The peak difference between the sum of the transients + residual and the source (done in sclang) is
-0.015426635742188 - -0.015437006950378 = 1.0371208e-05

Now, SC doesn't have the facilities for working out the ULP (as far as can see), but with numpy

(np.float32(-0.015426635742188) - np.float32(-0.015437006950378)) / np.spacing(np.float32(-0.015426635742188))

comes out at a chunky 11136. Did I do that right though?

[AlexHarker]

Well - we don't need to get as technical as ULP probably, although over a small number of simple calculations you are looking for a few ULP at most lost accuracy, unless there is some significant reason that accuracy is lost somewhere.

If the numbers above are correct that is at the 4th significant figures or around the 10th significant bit. That seems far too much - something funky is going on, but it is hard to see how it could be in the transient separator given the way it forms it output.

I suspect this might be an investigation of the testing system, rather than the algo - but I would need to see the whole of the test process to really understand it. For now the next step is to see what the two components are for the step that has that "large" error.

[weefuzzy]

So, for the afflicted area, it seems like something goes batshit. Here is the source waveform over 1000 samples around the glitch

Looks innocent enough. But here is the transient channel output. -250 seems excessive as a value.

Here is the residual

Source and residual together, on a sensible scale,

and some detritus when we do (trans + residual) - source

[AlexHarker]

Well - that explains the null summing error perfectly - so this is definitely something to do with what is coming out of the transient separator and the lack of null summing accuracy is the fact that one component has gone well about 0dB and we are losing accuracy in the bit we care about (so the -250 is cancelled but not perfectly because it obliterates some of the low bits in sum).

Next thing to investigate is why this is the output and what this means for the transient separator. I'd very much hope that that section is a detected transient (not that it necessarily should be) and thus the interpolator is going a bit haywire - then we need to figure out why.

[weefuzzy]

Source is Tremblay-AaS-SynthTwoVoices-M.wav

Params are default:

order	20
blockSize	256
padSize	128
threshFwd	2
threshBack	1.1
windowSize	14
clumpLength	25

[AlexHarker]

OK. See (original code from 2018):

flucoma-core/include/algorithms/ARModel.hpp

Line 151 in 2011fc0

for (int i = 0; i < mOrder + size; i++)

And now see (code now):

flucoma-core/include/algorithms/util/ARModel.hpp

Line 166 in e01e237

for (index i = mParameters.size(); i < mParameters.size() + size; i++)

Which are not the same loop size/points.

I see why that's been done (to avoid reading from input that is prior to the pointer), but as far as I can see the first part of the vector isn't explicitly initialised (it may default to zeros - not sure) and in fact this pre-reading was intentional (it's another example of the sort of issue already discussed about the prediction methods and should probably be addressed similarly). This may be unrelated to the issues of accuracy, but it should be addressed.

I didn't make that change.

[weefuzzy]

Good find. That fixes it, thanks. I can see you didn't make the change 😉

[AlexHarker]

Good. I hoped it might, but I wasn't sure. Happy to discuss why I do this, but it relates to the vagueness in the Kleiner/Martin/Thomson paper about how to deal with the start of time series (they don't tell you if I recall).

It's also important in the context of the transient code (if I recall correctly) that the input that gets passed in at the beginning of the blocks as you move through time is the residual (in other words because the beginning of the estimates array never updates, it's important that we feed it something that has a high probility of already being "correct" - which in this case means that the transients/outliers have already been removed).

That's what happens here (might be worth correcting the comment that should say "into" but has been messed up by an "int"->"index" find replace situation:

flucoma-core/include/algorithms/public/TransientExtraction.hpp

Line 292 in cdbe03c

// Copy the residual indexo the correct place

It's all coming back to me now...

[weefuzzy]

superseded by #117