TL;DR: 2*1LSB triangular-pdf dithering breaks in edgecases at 0 and 1 due to clamping. A solution is to lerp to a 1bit uniform dither in those edgecases.
I am adding a second answer, seeing as this turned out a bit more complicated than I originally thought. It appears this issue has been a "TODO: needs clamping?" in my code since I switched from normalized to triangular dithering... in 2012. Feels good to finally look at it :) Full code for solution / images used throughout the post: https://www.shadertoy.com/view/llXfzS
First of all, here is the problem we are looking at, when quantizing a signal to 3bits with 2*1LSB triangular-pdf dithering:
- essentially what hotmultimedia showed.
Increasing contrast, the effect described in the question becomes apparent: The output does not average to black/white in the edgecases (and actually extends well beyond 0/1 before doing so).
Looking at a graph provides a bit more insight:
(grey lines mark 0/1, also in grey is the signal we are trying to output, yellow-line is average of dithered/quantized output, red is the error (signal-average)).
Interestingly, not only is the average output not 0/1 at the limits, it is also not linear (likely due to the triangular pdf of the noise). Looking at the lower end, it makes intuitive sense why the output diverges: As the dithered signal starts to include negative values, the clamping-on-output changes the value of the lower dithered parts of the output (ie the negative values), thereby increasing the value of the average. An illustration appears to be in order (uniform, symmetric 2LSB dither, average still in yellow):
Now, if we just use a 1LSB normalized dither, there are no issues at the edges-cases at all, but then of course we lose the nice properties of triangular dithering (see e.g. this presentation).
A (pragmatic, empirical) solution (hack) then, is to revert to [-0.5;0.5[ uniform dithering for the edgecase:
float dithertri = (rnd.x + rnd.y - 1.0); //note: symmetric, triangular dither, [-1;1[
float dithernorm = rnd.x - 0.5; //note: symmetric, uniform dither [-0.5;0.5[
float sizt_lo = clamp( v/(0.5/7.0), 0.0, 1.0 );
float sizt_hi = 1.0 - clamp( (v-6.5/7.0)/(1.0-6.5/7.0), 0.0, 1.0 );
dither = lerp( dithernorm, dithertri, min(sizt_lo, sizt_hi) );
Which fixes the edgecases while keeping the triangular dithering intact for the remaining range:
So to not answer your question: I do not know if there is a more mathematically solid solution, and am equally keen to know what Masters of Past have done :) Until then, at least we have this horrific hack to keep our code functioning.
I should probably cover the workaround-suggestion given in The Question, on simply compressing the signal. Because the average is not linear in the edgecases, simply compressing the input-signal does not produce a perfect result - though it does fix the endpoints: