Shannon Exact64 — Shannon interpolation taken perhaps a little too seriously

[Moderator edit: Project moved to Tinkering]

My challenge to @jussi_laako author of HQP.

Shannon Exact64 is an offline 44.1 → 88.2 kHz sample-rate converter built around one deliberately uncompromising idea:

Take Shannon interpolation literally.

No windowed-sinc approximation. No arbitrary finite tap count. No deliberately introduced transition band. The interpolation is evaluated for the complete finite input file, with zero extension beyond its boundaries.

And because merely calculating the result wasn’t pedantic enough, Exact64 also tries to prove the correctly rounded IEEE-754 binary64 result.

It starts with a high-precision P96 certification pass. If even a single output sample cannot be unambiguously certified to its binary64 value, the complete certified pass is escalated to P128. There is no approximate fallback.

The accompanying telemetry log reports when this happens, which samples required escalation, the result after P128, and the resulting sample-peak/overshoot information.

In short:

No magic. No new sampling theory. Just textbook Shannon interpolation pursued with a slightly unreasonable level of numerical pedantry.

Shannon Exact64 — Quick Start

No installation is required.
No Visual Studio, runtime, or additional software is needed.

Download ShannonExact64.exe.
Drag and drop a 44.1 kHz / 16-bit stereo WAV file onto the EXE.
Do not press Enter or click anything while SRC is running. Processing can take some time — this is normal and does not mean the application has frozen.
When finished, the 88.2 kHz / 64-bit float WAV and _LOG.txt will appear in the same folder as the original WAV.

That’s it. Do not double-click the EXE — drag the WAV onto it instead.

https://drive.google.com/drive/folders/11l4urcGNRxS_3lqxpRJt7nOcFVn89gs0?usp=sharing

This sounds all good in theory, but there are some significant drawbacks to this approach:

  • It’s offline only. I’m not sure why this is challenging Jussi and HQP, since HQP is designed to do this in real time. The functionality is not overlapping.
  • I assume it’s not using FFT? FFT algorithms not only run faster, they also reduce the errors by performing fewer operations. Is this why you need to go to 128 bits (I assume that’s what P128 means)?
  • If you pad all tracks with zeros, you break gapless playback. Some albums (e.g. Pink Floyd 's) are continuous, with no real track boundaries. How long would it take to convert a 1hr continuous album, and what kind of errors would it entail?
  • Shannon’s reconstruction is supposed to retain all frequencies, all the way up to the Nyquist limit. That is not feasible practically, since frequencies close to the limit are problematic. For example, a track that consists of alternating +1/-1 samples is aliased, so it cannot be uniquely up-sampled. It’s best to filter out those frequencies.

Exact64 is deliberately not intended to compete with real-time SRCs on throughput or latency. Offline-only is a design choice. The goal is different: given a complete finite 44.1 kHz PCM sequence, evaluate its zero-extended Shannon reconstruction at the 88.2 kHz sample positions and determine the correctly rounded IEEE-754 binary64 result.

P96/P128 is therefore not simply “using more precision because the algorithm accumulates too much numerical error.” P96 is used as a certification precision. For each newly interpolated sample, the implementation attempts to prove which binary64 rounding cell contains the exact result. If P96 cannot establish that unambiguously, the certified computation is repeated at P128. If P128 still cannot certify it, there is no approximate fallback.

And that is an important part of the philosophy behind Exact64: it is deliberately unforgiving. It does not try to “improve” the source or hide its native imperfections. Whatever is present in the 44.1 kHz discrete sequence — including awkward content near Nyquist — is part of the input it is asked to reconstruct. It simply applies the specified Shannon convolution as faithfully as possible. No protective roll-off, no psychoacoustic optimization, no attempt to make the result prettier.

P128 is not there to compensate for accumulated numerical error from "not using an FFT.

Exact64 uses transform-based acceleration (NTT/CRT) rather than evaluating every sinc term by naive direct floating-point summation. This dramatically reduces computational complexity while using exact integer modular arithmetic for the convolution stage, avoiding the accumulation of floating-point rounding error that a huge direct summation would introduce.

Just think of it this way: Shannon Exact64 is exact down to the very last bit that 64-bit floating point can represent. It isn’t merely very close to the mathematically exact result — it is the correctly rounded binary64 result.

Can you explain what that means exactly? Unless there is some rounding done on the partial results, you’d need a precision with more bits than there are bits in the whole audio track for each output sample. So, the length of the computation sequence does have an effect on the final accuracy.

Also, zero-extending a track is fine if that track is independent. If however the track is part of a gapless album, a zero extension will simply apply a rectangular window, so you’ll end up with a windowed sinc with the worst kind of window.

P96 and P128 are simply the internal precision used to make sure that final rounding is correct. They are not output bit depths and they are not different quality modes — the output is always 64-bit floating point.

Think of it like calculating an answer with extra decimal places before rounding it. Exact64 first calculates with roughly 96 bits of precision. If that is enough to determine the correct final binary64 value with certainty, it is done.

In the rare case where the exact result lies so close to the boundary between two possible binary64 numbers that 96 bits are not enough to prove which one is correct, Exact64 recalculates the certification at 128-bit precision.

So P128 does not mean that 96-bit calculation was “bad” or inaccurate. It simply means that more precision was required to prove which of two extremely close 64-bit numbers is the correctly rounded answer.

The easiest analogy is calculating 1.23456… before writing the answer as 1.235. Usually a few extra digits are enough to know how to round. If the value happens to be extremely close to 1.2345…, you calculate more digits before making the decision.

That is essentially what P96 → P128 does — just at binary64 precision and for every interpolated sample.

You keep repeating the same thing, but not getting my point. You cannot just round only once at the end, you have to round intermediate results. Otherwise, you’d have to accumulate a huge number of bits until the end.

Say the input track uses 16 bits per sample. An IEEE 754 binary64 number has a significand with 53 bits. When you multiply the two, you get a result with 16+53=69 bits. If you don’t want to round it, you have to keep all of them. Now, when you multiply the next sample, you’ll get a result with 69+16=85 bits etc. Process the next 3 samples, and you end up with 133 bits, which is already more than 128. If you want to skip rounding until the end, you’ll have to accumulate 16 additional bits every sample, so you’ll end up with more bits in the result that there are bits in the original audio channel. For a 1-minute track at 44.1kHz, you’ll end up with more than 42 million bits of precision before you can finally round. For 24bit tracks, it’s over 65 million bits.

@Marian

I agree that you cannot simply ignore rounding of the intermediate operations. Of course P96/P128 arithmetic itself has finite precision and intermediate operations are rounded. My point is that these rounding errors are accounted for in the certification bounds. The goal is not to perform the complete calculation exactly and round only once at the very end, but to prove that despite the intermediate rounding, the exact Shannon result can only round to one particular binary64 value.

What I disagree with is the idea that this requires another 16 bits for every input sample. Shannon interpolation is a sum of products, not a chain of multiplications.

Answer was made with help from Chat gpt because I am a German

The central part is odd64_pass(). There, the Shannon coefficients are represented using integer bounds, the convolution is performed using NTT + CRT with integer arithmetic, and afterwards two endpoints are constructed for each sample. In simplified form, it works like this:

C = exact integer convolution result

lower = C + sum(negative input samples)
upper = C + sum(positive input samples)

RN64(lower) = lo
RN64(upper) = hi

if lo == hi:
certified
else:
unresolved

You can see this directly here:

uint64_t lo=round_binary64_bits(addSmall(C,neg),D),
hi=round_binary64_bits(addSmall(C,pos),D);

if(lo!=hi) {
++unresolved;
…
}

o[k]=lo;

And round_binary64_bits() explicitly implements round-to-nearest, ties-to-even.

The actual convolution is specifically not performed as a long chain of rounded floating-point multiplications. The NTT operates modularly with integers, and CRT then reconstructs the integer convolution result.

Instead, the uncertainty comes from the approximation of 1/\pi. For this, Exact64 uses certified floors:

\left\lfloor\frac{2^{97}}{\pi}\right\rfloor

and

\left\lfloor\frac{2^{129}}{\pi}\right\rfloor.

The direction of the possible error is then accounted for according to the sign of the individual PCM samples. This is why there are neg and pos terms. From these, an enclosure of the actual Shannon result is constructed.

And only if both endpoints round to exactly the same binary64 bit pattern is the sample accepted as certified.

Ok, so the more roundings you do, the more errors you accumulate. The result may be the closest binary64 approximation you can get, but it’s the closest for that particular set of operations, not the closest to the theoretical, infinite precision result, which doesn’t depend on the method. Given that an FFT greatly reduces the number of operations, its best approximation will probably be closer to the theoretical than the best approximation of convolution.

You are absolutely right, my bad. Yes, all multiplications are independent. Still, for 16bit samples, you’d have to add as many 65bit numbers as there are samples, and each such addition has the potential of adding one extra carry bit that you would have to retain. Again, for a 1-minute 44.1kHz track, there are about 2.6M samples, which could, in theory, add as many bits. Not as many as I initially estimated, but still a lot. Yes, samples usually alternate between positive and negative, so you’ll not only add, but also subtract; however, the sinc function also alternates, and you can get into a situation where many products are positive. The obvious example is doing a 2x up-sampling of the alternating +1/-1 signal.

Also, when you multiply two numbers that are orders of magnitude apart, the result may still be exact. For example, multiplying 1e100 by 2 will only change the exponent and remain exact. On the other hand, if you add 1e100 and 2, you’ll get 1e100, which is not exact anymore. Since the sinc function attenuates with the distance from the central peak, at one point, the products may be small enough to fall off the partial sum.
And finally, the accuracy of the result likely depends on the track length. With FFT, which has a fixed buffer length, accuracy is predictable.

But I think we’re rat-holing a bit too much on this. It’s unlikely any precision gains will lead to any audible gains when using binary64. Two issues remain:

  • It’s offline, so it requires storing up-sampled tracks and thus wasting storage.
  • There is the implicit rectangular windowing of gapless tracks.

I think there is still a misunderstanding about how Exact64 actually works.

Exact64 does not calculate the Shannon convolution by adding millions of floating-point products one after another. The convolution uses NTT/CRT with integer arithmetic, so small contributions cannot simply disappear because they are added to a much larger floating-point accumulator.

The P96/P128 certification then constructs lower and upper bounds which are intended to contain the theoretical Shannon result.

Both bounds are rounded to binary64 using round-to-nearest, ties-to-even. If both give exactly the same binary64 bit pattern, the result is certified. If they differ, Exact64 repeats the certification at P128. If the result still cannot be certified, it fails instead of returning an approximation.

So the goal is not to find the closest binary64 value to the result of a particular sequence of floating-point operations. The goal is to find the correctly rounded binary64 value of the theoretical Shannon reconstruction itself.

Of course, this depends on the mathematical correctness of the bounds and the NTT/CRT implementation. These are exactly the parts I still want to validate independently.

Your point about floating-point accumulation is correct in general. I just don’t think it applies to the way Exact64 performs its convolution.

My AI has checked the Code

I can send you if you want

Thanks. I looked up NTT and it seems it is a form of FFT after all over an integer modulo field. Regardless, the convolution can be computed directly in modulo arithmetic. The modulo keeps the integers bounded, so they require a relatively small number of bits. Doing this for 2 or 3 smaller primes can further reduce the number of bits, and the CRT will recover the exact integers within the two bounds computed from samples. The sums of samples are also small enough to fit in a relatively small number of bits. I’m not sure if there are pathological cases that can break these assumptions, but they can be detected and flagged. I learned something new, and I will concede that the number of additions doesn’t matter. However, the number of multiplications does matter, because the number of “kernel coefficients”, which are not exact, depends on the number of samples. A fixed-buffer FFT wouldn’t have this problem.

Now, regarding the graphs: They show that for both algorithms, the stop band attenuation is 200dB or better (which is better than the theoretical SNR of 32bit samples) and a pass band ripple of about 1E-9 decibels or better. Are there any humans on earth who will be able to distinguish the two?

No, I don’t expect anyone to hear a difference at these levels. :sweat_smile: But audibility isn’t really the point of Exact64.

What I find interesting in those graphs is that Exact64 achieves its passband behaviour without windowing the Shannon kernel and without deliberately introducing a transition band below Nyquist. It simply follows the specified Shannon reconstruction.

So it hasn’t traded away bandwidth in order to obtain prettier filter numbers. In fact, the measurement shows Exact64 remaining essentially flat closer to the 22.05 kHz Nyquist boundary.

ARDFTsrc is obviously already far beyond anything required for audibility. Exact64 is simply pursuing a different goal: how faithfully can the unwindowed Shannon reconstruction itself be evaluated and represented in binary64?

At this point the remaining differences are interesting mathematically, not because I expect anyone to hear them.

AI generated Text sent by Julian

So the purpose of exact64 is to extend the bandwidth by a few hertz, when it’s already above 22kHz? While at the same time doing worse in terms of passband ripple and stopband attenuation? Again, if this is not about audibility, what is it about? Why would one want to convert the files offline and use a few times more storage for no practical purpose? You can of course use NTT/CRT with fixed-buffer FFT if you want, and still use it in real-time, but again, I don’t think it will result in any audible improvement.

Also, you still haven’t addressed the issue of windowing in gapless tracks. That’s both a practical and a theoretical problem.

Yes, that was a deliberate design decision. Exact64 treats each input file as one finite discrete sequence and defines the samples outside that sequence as zero.

The reason was to make the reconstruction completely deterministic and mathematically well-defined without making assumptions about audio that is not contained in the file. At the beginning and end of a file, Exact64 has no knowledge of what came before or what comes after. Rather than inventing samples, extrapolating them, or applying a tapering window, I chose zero extension.

I agree that this becomes an issue with gapless albums. If two tracks are actually parts of one continuous recording, processing them independently introduces an artificial boundary that does not exist in the original continuous sequence. In that case, the mathematically cleaner solution would be to concatenate the gapless tracks, perform the Shannon reconstruction on the complete sequence, and split them again afterwards at the corresponding sample positions.

So I wouldn’t call the current behaviour an accidental windowing effect. It is an explicit boundary condition of Exact64: finite file + zero extension. But I agree that for truly gapless material, the album rather than the individual track should ideally be treated as the finite sequence.

AI generated text

Thanks. It’s been a nice theoretical debate. I’ll conclude it here: I don’t see this as a challenge to HQP or to established real-time algorithms. I would be interested in seeing the code; if that’s possible, please PM me a link.

From first post …