Some examples of dithering


In the following examples, c is the number of tones of gray allowed in dithering.

First of all, let us see the need of dithering. Consider for instance the following Matlab/Octave program that draws a spherical gradient and quantize it:
% GRADIENT
N = 400;
[X,Y]=meshgrid( linspace(-1,1.01,N), linspace(-1,1.01,N) );
Z = 0.99*sqrt( max(1-(X.^2+Y.^2),0) );
I = mat2gray(Z);
figure(1)
imshow(I)
imwrite(I, 'gradient.png')

% QUANTIZATION
c = 14
figure(2)
A = floor( c*Z );
I = mat2gray(A);
imshow(I)
imwrite(I, 'output.png')

The gradient is:
ogdith
and the corresponding quantized images for c=6 and c=14 gray tones are quite deffective:
ndith6 ndith14
c = 6
c = 14


The most pedestrian way of performing dithering is to apply a random noise before quantization. These are the results with uniformly distributed noise at each pixel:
gdith2 gdith6
c = 2
c = 6

gdith10 gdith14
c = 10
c = 14


Now we consider the following photo (500x500) as input to appreciate finer details:
tcol500
We work with the B/W version.

In ordered dithering, Bayer matrices are employed instead of random noise to assure the differences of near pixels in the middle tones. The Bayer matrices for dimension 2k, say Mk, follow a simple recurrence. Starting by M0=(0), the first block 2k-1x2k-1 is Mk-1, the second (to the right) is Mk-1 + 2·2-2k, the third is Mk-1 + 3·2-2k and the fourth Mk-1 + 2-2k.

These are examples with c=4 taking M2 and M4. The last one is very common in practice.
dith4_4
dith16_4
dim = 4, c = 4
dim = 16, c = 4

Finally, these are some examples with diffusion dithering, namely with Floyd-Steinberg algorithm. It distributes the quantization error among the neighboring pixels with a simple formula. The results are nevertheless impressive, especially taking into account how easy is to code it.

dif_dith2
dif_dith4
c = 2
c = 4


dif_dith2 dif_dith8
c = 6
c = 8

For more examples and theory on dithering see http://caca.zoy.org/study/index.html.

The code

This is the matlab/octave code employed to generate the images:


% note: octave: pkg load image

% oct = 1 -> octave (matlab otherwise)

oct = 1;

if oct == 1
  pkg load image
end


%----------------------------------------------------
%  GRADIENT
%----------------------------------------------------


disp(' ')
disp('============================================================')


N = 400;
[X,Y]=meshgrid( linspace(-1,1.01,N), linspace(-1,1.01,N) );
Z = sqrt( max(1-(X.^2+Y.^2),0) );

% Z = matrix with the gradient

I = mat2gray(Z);
imshow(I)
imwrite(I,'../images/ogdith.png');

disp('Gradient saved as "ogdith.png" (400x400)')
disp('============================================================')


%----------------------------------------------------
%  TRIVIAL QUANTIZATION
%----------------------------------------------------

for c = 2:4:14
  figure(1)
  A = floor( (c-1)*Z+0.5 )/(c-1);
  I = mat2gray(A);
  imshow(I)
  namefig = '../images/gndith%d.png';
  namefig = sprintf(namefig,c);
  imwrite(I,namefig);
end
disp('Fig 1. Trivial quantization saved as "gndith*.png"')
disp('============================================================')


%----------------------------------------------------
%  RANDOM NOISE
%----------------------------------------------------

for c = 2:4:14
  figure(2)
  nois = rand(N)-0.5;
  A = floor( (c-1)*Z+nois +0.5)/(c-1);
  I = mat2gray(A);
  imshow(I)
  namefig = '../images/gdith%d.png';
  namefig = sprintf(namefig,c);
  imwrite(I,namefig);
end
disp('Fig 2. Random noise saved as "gdith*.png"')
disp('============================================================')




%----------------------------------------------------
%  PHOTO
%----------------------------------------------------


disp(' ')
disp('============================================================')


imag = imread('../images/zorig_tcol500.jpg');
imag = rgb2gray(imag);
Z = im2double(imag)-0.00000001;
N = size(Z,1);

% Z = matrix with photo

disp('Read photo "zorig_tcol500.jpg" (500x500)')
disp('============================================================')



%----------------------------------------------------
%  TRIVIAL QUANTIZATION (P)
%----------------------------------------------------

for c = 2:4:6
  figure(3)
  A = floor( (c-1)*Z+0.5 )/(c-1);
  I = mat2gray(A);
  imshow(I)
  namefig = '../images/tndith%d.png';
  namefig = sprintf(namefig,c);
  imwrite(I,namefig);
end
disp('Fig 3. Trivial quantization (P) saved as "tndith*.png"')
disp('============================================================')


%----------------------------------------------------
%  RANDOM NOISE (P)
%----------------------------------------------------

for c = 2:4:6
  figure(4)
  nois = rand(N)-0.5;
  A = floor( (c-1)*Z+nois+0.5)/(c-1);
  I = mat2gray(A);
  imshow(I)
  namefig = '../images/tdith%d.png';
  namefig = sprintf(namefig,c);
  imwrite(I,namefig);
end
disp('Fig 4. Random noise (P) saved as "tdith*.png"')
disp('============================================================')



%----------------------------------------------------
%  ORDERED DITHERING
%----------------------------------------------------


M1 = [0,2;3,1]/4;
M2 = [M1,M1+2/16;M1+3/16,M1+1/16];
M3 = [M2,M2+2/64;M2+3/64,M2+1/64];
M4 = [M3,M3+2/256;M3+3/256,M3+1/256];
M5 = [M4,M4+2/1024;M4+3/1024,M4+1/1024];

% only using M = M2 (4x4) and M = M4 (16x16)

for mat = 4:12:16
  if mat==4
    M = M2;
  elseif mat==16
    M = M4;
  end

  Nsub = floor(N/size(M,1)) +1;
  DT = repmat(M, 1, Nsub);
  DT = repmat(DT, Nsub, 1);
  DT = DT(1:N,1:N);
  % DT = full size dithering matrix

  for c = 2:2:8
    figure(5)
    A = floor( (c-1)*Z+DT)/(c-1);
    I = mat2gray(A);
    imshow(I)
    namefig = '../images/ord_dith%d_%d.png';
    namefig = sprintf(namefig,size(M,1),c);
    imwrite(I,namefig)
  end
end
disp('Fig 5. Ordered dithering saved as "ord_dith*_*.png"')
disp('============================================================')



%----------------------------------------------------
%  FLOYD-STEINBERG DITHERING
%----------------------------------------------------

for c = 2:2:8
  figure(6)
  A = (c-1)*Z;
  for y = 1:N
    for x = 1:N
      oldpixel = A(x,y);
      newpixel = round(oldpixel);
      A(x,y) = newpixel;
      quant_error = oldpixel - newpixel;
      if x<N
        A(x + 1,y)= A(x + 1,y) + quant_error * 7 / 16;
      end
      if ( x>1 ) && ( y<N )
        A(x - 1,y + 1) = A(x - 1,y + 1) + quant_error * 3 / 16;
      end
      if y<N
        A(x ,y + 1) = A(x ,y + 1) + quant_error * 5 / 16;
      end
      if ( x<N ) && ( y<N )
        A(x + 1,y + 1) = A(x + 1,y + 1) + quant_error * 1 / 16;
      end
    end
  end

  I = mat2gray(A);
  imshow(I)
  namefig = '../images/dif_dith%d.png';
  namefig = sprintf(namefig,c);
  imwrite(I,namefig)
end

disp('Fig 6. Floyd-Steinberg dithering saved as "dif_dith*.png"')
disp('============================================================')