dithering using 1x2 blocks

Question

Is there any software, which can do "block error diffusion" dithering as described here and here? Often I end up with completely separated pixels of the same colour, when doing Floyd Steinberg dithering. I would like pixel blocks of 1x2 or 2x2 pixels of the same colour to be united. Is it possible with any existing software? I am happy to write the code myself. I just need to know what the approach is.

I only want the final image to contain pixel blocks of size 1x2 like these:

Answer 1

The study of the library to which you linked includes all of the python code necessary to generate all of the images examples contained in the document.

The attaced code is a copy of the relevant functions needed to create an error corrected dither of arbitrary tiles. It takes an image as an input(foo.png) and creates two PNG files as output (foo_GreyDither.png, foo_BWDither.png)

For example, python ditherCode.py foo.png

#!/usr/bin/env python

import math, gd, random, sys, os

class Image(gd.image):
    gd.gdMaxColors = 256 * 256 * 256
    def __init__(self, *args):
        if args[0].__class__ == str:
            print "[LOAD] %s" % (args[0],)
        gd.image.__init__(self, *args)
    def save(self, name):
        print "[PNG] %s" % (name,)
        self.writePng(name)
    def getGray(self, x, y):
        p = self.getPixel((x, y))
        c = self.colorComponents(p)[0] / 255.0
        return c
    def getRgb(self, x, y):
        p = self.getPixel((x, y))
        rgb = self.colorComponents(p)
        return [rgb[0] / 255.0, rgb[1] / 255.0, rgb[2] / 255.0]
    def setGray(self, x, y, t):
        p = (int)(t * 255.999)
        c = self.colorResolve((p, p, p))
        self.setPixel((x, y), c)
    def setRgb(self, x, y, r, g, b):
        r = (int)(r * 255.999)
        g = (int)(g * 255.999)
        b = (int)(b * 255.999)
        c = self.colorResolve((r, g, b))
        self.setPixel((x, y), c)
    def getRegion(self, x, y, w, h):
        dest = Image((w, h), True)
        self.copyTo(dest, (-x, -y))
        return dest
    def getZoom(self, z):
        (w, h) = self.size()
        dest = Image((w * z, h * z), True)
        for y in range(h):
            for x in range(w):
                rgb = self.getRgb(x, y)
                for j in range(z):
                    for i in range(z):
                        dest.setRgb(x * z + i, y * z + j, *rgb)
        return dest

# Manipulate gamma values
class Gamma:
    def CtoI(x):
        if x < 0:
            return - math.pow(-x, 2.2)
        return math.pow(x, 2.2)
    def ItoC(x):
        if x < 0:
            return - math.pow(-x, 1 / 2.2)
        return math.pow(x, 1 / 2.2)
    CtoI = staticmethod(CtoI)
    ItoC = staticmethod(ItoC)
    def CtoI3(x):
        return [Gamma.CtoI(x[0]), Gamma.CtoI(x[1]), Gamma.CtoI(x[2])]
    def ItoC3(x):
        return [Gamma.ItoC(x[0]), Gamma.ItoC(x[1]), Gamma.ItoC(x[2])]
    CtoI3 = staticmethod(CtoI3)
    ItoC3 = staticmethod(ItoC3)
    def Cto2(x):
        if x < Gamma.CtoI(0.50):
            return 0.
        return 1.
    def Cto3(x):
        if x < Gamma.CtoI(0.25):
            return 0.
        elif x < Gamma.CtoI(0.75):
            return Gamma.CtoI(0.5)
        return 1.
    def Cto4(x):
        if x < Gamma.CtoI(0.17):
            return 0.
        elif x < Gamma.CtoI(0.50):
            return Gamma.CtoI(0.3333)
        elif x < Gamma.CtoI(0.83):
            return Gamma.CtoI(0.6666)
        return 1.
    Cto2 = staticmethod(Cto2)
    Cto3 = staticmethod(Cto3)
    Cto4 = staticmethod(Cto4)

# Create matrices
def Matrix(w, h, val = 0):
    return [[val] * w for n in range(h)]

# Iterate in 2D space
def rangexy(w, h):
    for y in range(h):
        for x in range(w):
            yield (x, y)

inputImage = Image(sys.argv[1])
# grad256bw = Image((32, 256))
# for x, y in rangexy(32, 256):
#     grad256bw.setGray(x, 255 - y, y / 255.)
# grad256bw.save("gradient256bw.png")

def subblock(src, tiles, propagate, diff, gamma):
    (w, h) = src.size()
    # Gamma correction
    if gamma:
        ctoi = Gamma.CtoI
        itoc = Gamma.ItoC
    else:
        ctoi = itoc = lambda x : x
    # Propagating the error to a temporary buffer is becoming more and
    # more complicated. We decide to use an intermediate matrix instead.
    tmp = Matrix(w, h, 0.)
    for x, y in rangexy(w, h):
        tmp[y][x] = ctoi(src.getGray(x, y))
    dest = Image((w, h))
    # Analyse tile list
    ntiles = len(tiles)
    ty = len(tiles[0])
    tx = len(tiles[0][0])
    cur = Matrix(tx, ty, 0.)
    w, h = w / tx, h / ty
    # Analyse error propagate list
    for x, y in rangexy(w, h):
        # Get block value
        for i, j in rangexy(tx, ty):
            cur[j][i] = itoc(tmp[y * ty + j][x * tx + i])
        # Select closest block
        dist = tx * ty
        for n in range(ntiles):
            d = 0.
            e = 0.
            for i, j in rangexy(tx, ty):
                d += cur[j][i] - tiles[n][j][i]
                e += diff[j][i] * abs(cur[j][i] - tiles[n][j][i])
            if abs(d) / (tx * ty) + e < dist:
                dist = abs(d) / (tx * ty) + e
                best = n
        # Set pixel
        for i, j in rangexy(tx, ty):
            dest.setGray(x * tx + i, y * ty + j, tiles[best][j][i])
        # Propagate error
        for i, j in rangexy(tx, ty):
            e = ctoi(cur[j][i]) - ctoi(tiles[best][j][i])
            m = propagate[j][i]
            for px, py in rangexy(len(m[0]), len(m)):
                if m[py][px] == 0:
                    continue
                if m[py][px] == -1:
                    cx, cy = px, py
                    continue
                tmpx = x * tx + i + px - cx
                tmpy = y * ty + j + py - cy
                if tmpx > w * tx - 1 or tmpy > h * ty - 1:
                    continue
                tmp[tmpy][tmpx] += m[py][px] * e
    return dest

ERROR_SUBFS22 = \
    [[[[0, -1, 0, 8./64],
       [0, 0, 0, 10./64],
       [7./64, 22./64, 15./64, 2./64]],
      [[0, 0, -1, 20./64],
       [0, 0, 0, 14./64],
       [2./64, 11./64, 15./64, 2./64]]],
     [[[0, 0, 0, 0./64],
       [0, -1, 0, 6./64],
       [12./64, 32./64, 13./64, 1./64]],
      [[0, 0, 0, 0./64],
       [0, 0, -1, 20./64],
       [0./64, 12./64, 28./64, 4./64]]]]

DIFF_WEIGHTED22 = \
    [[51./128, 33./128],
     [25./128, 19./128]]

GREYLINES22 = []
for n in range(4*4*4*4):
    vals = [0., 0.333, 0.666, 1.]
    a, b, c, d = n & 3, (n >> 2) & 3, (n >> 4) & 3, (n >> 6) & 3
    if (a != b or c != d) and (a != c or b != d):
        continue
    GREYLINES22.append([[vals[a], vals[b]], [vals[c], vals[d]]])

LINES22 = \
    [[[0., 0.], [0., 0.]],
     [[0., 1.], [0., 1.]],
     [[1., 0.], [1., 0.]],
     [[1., 1.], [0., 0.]],
     [[0., 0.], [1., 1.]],
     [[1., 1.], [1., 1.]]]

foo = sys.argv[1]
foo = foo[:-4]

subblock(inputImage, GREYLINES22,
         ERROR_SUBFS22, DIFF_WEIGHTED22, False).save(foo+"_GreyDither.png")
subblock(inputImage, LINES22,
         ERROR_SUBFS22, DIFF_WEIGHTED22, False).save(foo+"_BWDither.png")