Blue Noise Dithering

Initial seed pattern

This is a short description on blue noise dithering. For more detail see Robert Ulichney’s 1993 paper about blue noise dithering.

Initially a binary pattern BP is created where N pixels (about 10 percent) are set to 1 and the rest to 0.

The binary pattern then is convolved with the following filter function to generate the density array DA:

You can use a sigma value of 1.5.

The convolution is wrapped around to facilitate a tileable result:

Maxima of the density array are called clusters and minima are called voids.

The 1 value in BP with the highest density value DA (tightest cluster) is set to 0 and DA is updated accordingly. Now the 0 value in BP with the lowest density value DA (largest void) is set to 1 (and DA is updated).

This is repeated until disolving the tightest cluster creates the largest void. This is done to spread the 1 values evenly.

Phase 1

In phase 1 of the dithering algorithm the 1 values of a copy of the seed pattern are removed one by one starting where the density DA is the highest. A copy of the density array DA is updated accordingly. The corresponding positions in the resulting dither array are set to N-1, N-2, …, 0.

Phase 2

In phase 2 starting with the seed pattern a mask is filled with 1 values where the density DA is the lowest. The density array DA is updated while filling in 1 values. Phase 2 stops when half of the values in the mask are 1. The corresponding positions in the dither array are set to N, N+1, …, (M * M) / 2 - 1

Phase 3

In phase 3 the density array DA is recomputed using the boolean negated mask from the previous phase (0 becomes 1 and 1 becomes 0). Now the mask is filled with 1 values where the density DA is the highest (clusters of 0s) always updating DA. Phase 3 stops when all the values in the mask are 1. The corresponding positions in the dither array are set to (M * M) / 2, …, M * M - 1.

Result

The result can be normalised to 0 to 255 in order to inspect it. The blue noise dither array looks as follows:

Here is an example with constant offsets when sampling 3D clouds without dithering.

Here is the same scene using dithering to set the sampling offsets.

One can apply a blur filter to reduce the noise.

Note how the blurred image shows more detail than the image with constant offsets even though the sampling rate is the same.

Let me know any comments/suggestions in the comments below.

Test Driven Development with OpenGL

Test Driven Development

Test driven development (TDD) undoubtedly helps a great deal in preventing development grinding to a halt once a project’s size surpasses a few lines of code.

1. The reason for first writing a failing test is to ensure that the test is actually failing and testing the next code change.
2. A minimal change to the code is performed to pass the new test while also still passing all previously written tests.
3. If necessary the code is refactored/simplified. The reason to do this after passing the test is so that one does not have to worry about passing the test and writing clean code at the same time.

Testing rendering output

One can test OpenGL programs by rendering test images and comparing them with a saved image (a test fixture). In order to automate this, one can perform offscreen rendering and do a pixel-wise image comparison with the saved image.

Using the Clojure programming language and the Lightweight Java Game Library (LWJGL) one can perform offscreen rendering with a Pbuffer object using the following macro (of course this approach is not limited to Clojure and LWJGL):

; Set up OpenGL modes
(defn setup-rendering
  "Common code for setting up rendering"
  [width height]
  (GL11/glViewport 0 0 width height)
  (GL11/glEnable GL11/GL_DEPTH_TEST)
  (GL11/glEnable GL11/GL_CULL_FACE)
  (GL11/glCullFace GL11/GL_BACK)
  (GL11/glDepthFunc GL11/GL_GEQUAL)
  (GL45/glClipControl GL20/GL_LOWER_LEFT GL45/GL_ZERO_TO_ONE))

; Macro accepting a block of code to be rendered to a Pbuffer
(defmacro offscreen-render
  "Macro to use a pbuffer for offscreen rendering"
  [width height & body]
  `(let [pixels#  (BufferUtils/createIntBuffer (* ~width ~height))
         pbuffer# (Pbuffer. ~width ~height (PixelFormat. 24 8 24 0 0) nil nil)
         data#    (int-array (* ~width ~height))]
     (.makeCurrent pbuffer#)
     (setup-rendering ~width ~height)
     (try
       ~@body
       (GL11/glReadPixels 0 0 ~width ~height GL12/GL_BGRA GL11/GL_UNSIGNED_BYTE pixels#)
       (.get pixels# data#)
       {:width ~width :height ~height :data data#}
       (finally
         (.releaseContext pbuffer#)
         (.destroy pbuffer#)))))

Note that the code sets up reversed-z rendering as discussed in an earlier article

Using the Midje testing library one can for example test a command for clearing the color buffer as follows:

; Test rendering a red background
(fact "Render background color"
  (offscreen-render 160 120 (clear (matrix [1.0 0.0 0.0])))
  => (is-image "test/sfsim25/fixtures/render/red.png"))

The checker is-image is implemented using ImageJ:

; Implement equality test for RGBA images ignoring alpha channel.
(defn is-image
  "Compare RGB components of image and ignore alpha values."
  [filename]
  (fn [other]
      (let [img (slurp-image filename)]
        (and (= (:width img) (:width other))
             (= (:height img) (:height other))
             (= (map #(bit-and % 0x00ffffff) (:data img))
                (map #(bit-and % 0x00ffffff) (:data other)))))))

; Load an RGB image
(defn slurp-image
  "Load an RGB image"
  [^String file-name]
  (let [img (.openImage (Opener.) file-name)]
    (.convertToRGB (ImageConverter. img))
    {:width (.getWidth img)
     :height (.getHeight img)
     :data (.getPixels (.getProcessor img))}))

The image is recorded initially by using the checker record-image instead of is-image and verifying the result manually.

; Test method writing image to file instead of comparing it.
(defn record-image
  "Use this test function to record the image the first time."
  [filename]
  (fn [other]
      (spit-image filename other)))

; Method to save image.
(defn spit-image
  "Save RGB image as PNG file"
  [^String file-name {:keys [width height data]}]
  (let [processor (ColorProcessor. width height data)
        img       (ImagePlus.)]
    (.setProcessor img processor)
    (.saveAsPng (FileSaver. img) file-name)))

One can use this approach (and maybe only this approach) to test code for handling vertex array objects, textures, and for loading shaders.

Testing shader code

Above approach has the drawback that it can only test complete rendering programs. Also the output is limited to 24-bit RGB images. The tests are therefore more like integration tests and they are not suitable for unit testing shader functions.

However it is possible to use a Pbuffer just as a rendering context and perform rendering to a floating-point texture. One can use a texture with a single pixel as a framebuffer. A single pixel of a uniformly colored quad is drawn. The floating point channels of the texture’s RGB pixel then can be compared with the expected value.

; Use shaders to render to a single pixel texture and read out result.
(defn shader-test [setup probe & shaders]
  (fn [uniforms args]
      (let [result (promise)]
        (offscreen-render 1 1
          (let [indices  [0 1 3 2]
                vertices [-1.0 -1.0 0.5, 1.0 -1.0 0.5, -1.0 1.0 0.5, 1.0 1.0 0.5]
                program  (make-program :vertex [vertex-passthrough]
                                       :fragment (conj shaders (apply probe args)))
                vao      (make-vertex-array-object program indices vertices [:point 3])
                tex      (texture-render 1 1 true
                                         (use-program program)
                                         (apply setup program uniforms)
                                         (render-quads vao))
                img      (texture->vectors tex 1 1)]
            (deliver result (get-vector img 0 0))
            (destroy-texture tex)
            (destroy-vertex-array-object vao)
            (destroy-program program)))
        @result)))

; Minimal vertex shader.
(def vertex-passthrough "#version 410 core
in highp vec3 point;
void main()
{
  gl_Position = vec4(point, 1);
}")

; Macro to render output of a code block to a texture.
(defmacro texture-render
  "Macro to render to a texture"
  [width height floating-point & body]
  `(let [fbo# (GL45/glCreateFramebuffers)
         tex# (GL11/glGenTextures)]
     (try
       (GL30/glBindFramebuffer GL30/GL_FRAMEBUFFER fbo#)
       (GL11/glBindTexture GL11/GL_TEXTURE_2D tex#)
       (GL42/glTexStorage2D GL11/GL_TEXTURE_2D 1
         (if ~floating-point GL30/GL_RGB32F GL11/GL_RGBA8) ~width ~height)
       (GL32/glFramebufferTexture GL30/GL_FRAMEBUFFER GL30/GL_COLOR_ATTACHMENT0 tex# 0)
       (GL20/glDrawBuffers (make-int-buffer (int-array [GL30/GL_COLOR_ATTACHMENT0])))
       (GL11/glViewport 0 0 ~width ~height)
       ~@body
       {:texture tex# :target GL11/GL_TEXTURE_2D}
       (finally
         (GL30/glBindFramebuffer GL30/GL_FRAMEBUFFER 0)
         (GL30/glDeleteFramebuffers fbo#)))))

; Read out texture and return floating point RGB image.
(defn texture->vectors
  "Extract floating-point vectors from texture"
  [texture width height]
  (with-2d-texture (:texture texture)
    (let [buf  (BufferUtils/createFloatBuffer (* width height 3))
          data (float-array (* width height 3))]
      (GL11/glGetTexImage GL11/GL_TEXTURE_2D 0 GL12/GL_BGR GL11/GL_FLOAT buf)
      (.get buf data)
      {:width width :height height :data data})))

Furthermore it is possible to compose the fragment shader by linking the shader function under test with a main function. I.e. it is possible to link the shader function under test with a main function implemented just for probing the shader.

The shader-test function defines a test function using the probing shader and the shader under test. The new test function then can be used using the Midje tabular environment. In the following example the GLSL function phase is tested. Note that parameters in the probing shaders are set using the weavejester/comb templating library.

; Probing shader generated using templating.
(def phase-probe
  (template/fn [g mu] "#version 410 core
out lowp vec3 fragColor;
float phase(float g, float mu);
void main()
{
  float result = phase(<%= g %>, <%= mu %>);
  fragColor = vec3(result, 0, 0);
}"))

; Use shader test function and probing shader to test a shadder.
(def phase-test (shader-test (fn [program]) phase-probe phase-function))

; Midje tests using test function to test a shader.
(tabular "Shader function for scattering phase function"
         (fact (mget (phase-test [] [?g ?mu]) 0) => (roughly ?result))
         ?g  ?mu ?result
         0   0   (/ 3 (* 16 PI))
         0   1   (/ 6 (* 16 PI))
         0  -1   (/ 6 (* 16 PI))
         0.5 0   (/ (* 3 0.75) (* 8 PI 2.25 (pow 1.25 1.5)))
         0.5 1   (/ (* 6 0.75) (* 8 PI 2.25 (pow 0.25 1.5))))

Note that using mget the red channel of the pixel is extracted. Sometimes it might be more desirable to check all channels of the RGB pixel.

Here is the actual implementation of the tested function:

#version 410 core

float M_PI = 3.14159265358;

float phase(float g, float mu)
{
  return 3 * (1 - g * g) * (1 + mu * mu) /
    (8 * M_PI * (2 + g * g) * pow(1 + g * g - 2 * g * mu, 1.5));
}

The empty function (fn [program]) is specified as a setup function. In general the setup function is used to initialise uniforms used in the shader under test.

Here is an example of tests using uniform values:

; Probing shader generated using templating.
(def transmittance-track-probe
  (template/fn [px py pz qx qy qz] "#version 410 core
out lowp vec3 fragColor;
vec3 transmittance_track(vec3 p, vec3 q);
void main()
{
  vec3 p = vec3(<%= px %>, <%= py %>, <%= pz %>);
  vec3 q = vec3(<%= qx %>, <%= qy %>, <%= qz %>);
  fragColor = transmittance_track(p, q);
}"))

; Shader test involving some uniform variables.
(def transmittance-track-test
  (transmittance-shader-test
    (fn [program height-size elevation-size elevation-power radius max-height]
        (uniform-int program :height_size height-size)
        (uniform-int program :elevation_size elevation-size)
        (uniform-float program :elevation_power elevation-power)
        (uniform-float program :radius radius)
        (uniform-float program :max_height max-height))
    transmittance-track-probe transmittance-track
    shaders/transmittance-forward shaders/horizon-angle
    shaders/elevation-to-index shaders/interpolate-2d
    shaders/convert-2d-index shaders/is-above-horizon))

; Midje tests to test a shader.
(tabular "Shader function to compute transmittance between two points in the atmosphere"
         (fact (mget (transmittance-track-test [17 17 1 6378000.0 100000.0]
                                               [?px ?py ?pz ?qx ?qy ?qz]) 0)
               => (roughly ?result 1e-6))
         ?px ?py ?pz     ?qx ?qy ?qz     ?result
         0   0   6478000 0   0   6478000 1
         0   0   6428000 0   0   6478000 0.5
         0   0   6453000 0   0   6478000 0.75
         0   0   6428000 0   0   6453000 (/ 0.5 0.75))

Here a setup function initialising 5 uniform values is specified.

Mocking shader functions

If each shader function is implemented as a separate string (loaded from a separate file), one can easily link with mock functions when testing shaders. Here is an example of a probing shader which also contains mocks to allow the shader to be unit tested in isolation:

; Probing shader generated using templating and defining mock shaders.
(def cloud-track-base-probe
  (template/fn [px qx n decay scatter density ir ig ib]
"#version 410 core
out lowp vec3 fragColor;
vec3 transmittance_forward(vec3 point, vec3 direction)
{
  float distance = 10 - point.x;
  float transmittance = exp(-<%= decay %> * distance);
  return vec3(transmittance, transmittance, transmittance);
}
vec3 ray_scatter_forward(vec3 point, vec3 direction, vec3 light)
{
  float distance = 10 - point.x;
  float amount = <%= scatter %> * (1 - pow(2, -distance));
  return vec3(0, 0, amount);
}
float cloud_density(vec3 point)
{
  return <%= density %>;
}
float phase(float g, float mu)
{
  return 1.0 + 0.5 * mu;
}
vec3 cloud_track_base(vec3 p, vec3 q, int n, vec3 incoming);
void main()
{
  vec3 p = vec3(<%= px %>, 0, 0);
  vec3 q = vec3(<%= qx %>, 0, 0);
  vec3 incoming = vec3(<%= ir %>, <%= ig %>, <%= ib %>);
  fragColor = cloud_track_base(p, q, <%= n %>, incoming);
}
"))

Let me know if you have any comments or suggestions.

Enjoy!

Updates:

• submitted to Reddit
• submitted to Hackernews
• Rhawk187 pointed out that exact image comparisons are also problematic because updates to graphics drivers can cause subtle changes. This can be adressed by allowing a small average difference between the expected and actual image.
• In order to be independent of the graphics driver one can run the tests using xvfb-run to use software rendering instead.

The code of make-vertex-array-object and render-quads is added here for reference.

(defn make-vertex-array-object
  "Create vertex array object and vertex buffer objects"
  [program indices vertices attributes]
  (let [vertex-array-object (GL30/glGenVertexArrays)]
    (GL30/glBindVertexArray vertex-array-object)
    (let [array-buffer (GL15/glGenBuffers)
          index-buffer (GL15/glGenBuffers)]
      (GL15/glBindBuffer GL15/GL_ARRAY_BUFFER array-buffer)
      (GL15/glBufferData GL15/GL_ARRAY_BUFFER
                         (make-float-buffer (float-array vertices)) GL15/GL_STATIC_DRAW)
      (GL15/glBindBuffer GL15/GL_ELEMENT_ARRAY_BUFFER index-buffer)
      (GL15/glBufferData GL15/GL_ELEMENT_ARRAY_BUFFER
                         (make-int-buffer (int-array indices)) GL15/GL_STATIC_DRAW)
      (let [attribute-pairs (partition 2 attributes)
            sizes           (map second attribute-pairs)
            stride          (apply + sizes)
            offsets         (reductions + (cons 0 (butlast sizes)))]
        (doseq [[i [attribute size] offset] (map list (range) attribute-pairs offsets)]
          (GL20/glVertexAttribPointer (GL20/glGetAttribLocation
                                        ^int program (name attribute))
                                        ^int size
                                      GL11/GL_FLOAT false
                                        ^int (* stride Float/BYTES)
                                        ^int (* offset Float/BYTES))
          (GL20/glEnableVertexAttribArray i))
        {:vertex-array-object vertex-array-object
         :array-buffer        array-buffer
         :index-buffer        index-buffer
         :nrows               (count indices)
         :ncols               (count attribute-pairs)}))))

(defn render-quads
  "Render one or more quads"
  [vertex-array-object]
  (setup-vertex-array-object vertex-array-object)
  (GL11/glDrawElements GL11/GL_QUADS
    ^int (:nrows vertex-array-object) GL11/GL_UNSIGNED_INT 0))

Debian 11 on the Tuxedo Aura 15 Gen1 laptop

Buying a laptop from Tuxedo Computers

I upgraded to a laptop from Tuxedo Computers (german page). The nice thing about Tuxedo Computers is that

• they provide GNU/Linux support for new hardware
• they let you customise the laptop (e.g. CPU, SSD size, RAM configuration, keyboard layout)
• the laptops are repairable/upgradeable
• the laptops come with 2 years warranty

Also note that Tuxedo Computers are now a KDE patron.

I chose the Tuxedo Aura 15 Gen1 (german page) with an AMD Ryzen 7 4700U CPU and an AMD Radeon RX Vega 7 GPU. I selected the laptop because of

• a keyboard with a number pad and backlight
• 15.6 inch display
• 8 core, 8 threads AMD CPU using only 15W power
• no proprietary NVidia graphics drivers

It is recommended when choosing the RAM to use a multi-channel configuration for better performance. Also make sure to order enough RAM because some of it is used by the integrated graphics card. Finally doing parallelized builds on many cores requires more memory.

Note that there also is the Tuxedo Pulse 15 Gen1 which comes with a 8 core 16 threads AMD Ryzen 7 4800H. Also the Tuxedo Aura 15 Gen2 is due to arrive in April and also comes with 16 threads (see Phoronix review).

Here is an in-depth review in German language (and here is a shorter review). Also there is a Youtube video review in English below:

Installing Debian 11

On Reddit I got helpful information from Ferdinand from Tuxedo Computers on how to install Debian 11 on the Aura 15 Gen1. Also see Debian Wiki for some information. Debian is not supported by Tuxedo Computers but it works nonetheless on a Tuxedo Aura 15 Gen1. I followed the following steps (no warranty) to get Debian 11 running:

• Download the latest Debian 11 unofficial image including non-free firmware and create a bootable USB stick using Rufus for example. The unofficial image comes with proprietary firmware for the AMD GPU and for the Wifi adapter.
• Boot with the USB stick plugged in, press F7 and select the USB stick.
• You need to plug in the network cable for installing the updates. The Wifi adapter is recognized but doesn’t work during installation.
• The Debian installer lets you choose which desktop (I chose KDE here).
• After installation you need to comment out the cdrom entries in /etc/apt/sources.list.

Installing Tuxedo packages and other packages

• You then need to add the Tuxedo package repository (see here for more information):

sudo nano /etc/apt/sources.list.d/tuxedocomputers.list
# Enter "deb https://deb.tuxedocomputers.com/ubuntu focal main" and save.
wget -O - https://deb.tuxedocomputers.com/0x54840598.pub.asc | sudo apt-key add -
sudo apt-get update

• Install firmware-linux.
• Manually install old libappindicator3-1 and libindicator3-7 package from Debian 10 (update: this is not necessary any more with the new version of Tuxedo Control Center thanks to Vinz)
• You can then install the following three Tuxedo packages:
  • tuxedo-control-center (user interface for power management) (source)
  • tuxedo-keyboard (dkms kernel module for keyboard lighting) (source)
  • tuxedo-restore-audio-fix (options snd-hda-intel probe_mask=1) (see package repository)
• Install packages from deb-multimedia.org.
• Install Plymouth boot splash.

Then a few more tweaks:

• I installed plasma-browser-integration and the corresponding browser addon (available for Firefox and Chrome/Chromium)
• I installed ffmpegthumbs and kdegraphics-thumbnailers and then enabled the thumbnail generators in the Dolphin settings
• I installed bleachbit and run it to clean up unused disk space

Finally you can list the installed packages on an old system to install them on the new machine. I.e. you can run the following command on the old machine:

echo aptitude install `aptitude search '!~M ~i' -F '%p'` > install.sh

Then you copy the resulting file install.sh to the target machine and run

sh ./install.sh

Enjoy!

Updates

• Had to edit 40-libinput.conf to get touchpad tapping to work
• Tuxedocomputers also gives you 10GByte of space on their Nextcloud server which you can mount using davfs2
• Press F2 to enter the BIOS; press F7 to select the boot device
• If you want to set up a Windows/Linux dual-boot (using Grub) you need to disable Bitlocker encryption
• If you want to mount an NTFS partition you need to disable fast boot on Windows
• The BIOS supports Flexicharger which lets you set charge limits to extend battery life
• You can also set the keyboard backlight color using webketje/tuxedo-backlight-control

You can now use the following content for /etc/apt/sources.list.d/tuxedo.sources:

Types: deb
URIs: https://deb.tuxedocomputers.com/debian
Suites: trixie
Components: main
Signed-By: /usr/share/keyrings/tuxedo-archive-keyring.gpg

You can install the keyring as follows:

wget https://deb.tuxedocomputers.com/ubuntu/pool/main/t/tuxedo-archive-keyring/tuxedo-archive-keyring_2022.04.01~tux_all.deb
sudo dpkg -i tuxedo-archive-keyring_2022.04.01~tux_all.deb

SpaceX Crew Dragon paper model

Sometime ago I discovered the AXM Paper Space Scale Models website. On the website Alfonso X Moreno published various paper scale models of the Space Shuttle, the ISS, the Tiangong, SpaceX Starship, and others.

I downloaded the SpaceX Crew Dragon paper model and I donated a bit of money for it. I used pdftk to split the PDF document into single- and double-sided jobs. Then I used doxdirect.com to get it printed, because they offer printing on 160gsm paper (normal office paper is only 80gsm and is not strong enough).

I couldn’t find instructions on the website on how to build the docking version of the paper model. Finally I found an instructional video on Youtube:

Here is the finished Crew Dragon capsule (I will do the service module later).

I can really recommend to check out AXM’s website if you are into paper models.

Enjoy and Merry Christmas!

Creating harmonica tabs with Lilypond

I play a bit of music using a Hohner Special 20 bluesharp. However I am not quick enough at reading sheet music, so I need to add tabs. Tabs basically tell you the number of the hole to play and whether to blow or draw (-). Lilypond is a concise programming language to generate sheet music. It turns out, that it is very easy to add lyrics using Lilypond. By putting the tabs in double quotes it is possible to add the tabs to the music. I found the tune Hard Times on the Internet and I used Lilypond to create a version with tabs. Here is the Lilypond code:

\version "2.22.0"
\header {
  title = "Hard Times"
  composer = "Stephen Foster"
  instrument = "Bluesharp"
}
\relative {
  r2 c'4 d
  e2 e4 d
  e8 g4.~ g4 e
  d4 c c d8 c
  e2 c'4 a
  g2 e8 c4.
  d4. c8 e4 d
  c1~
  c2 c4 d
  e2 e4 d
  e8 g4.~ g4 e
  d4 c c d8 c
  e2 c'4 a
  g2 e8 c4.
  d4. c8 e4 d
  c1~
  c2 e4 f4
  g2. g4
  g2 f4 g4
  a1
  g2 r2
  c2 a8 g4.
  e2 d8 c4.
  d4. c8 d4 e
  d2 c4 d
  e2 e4 d
  e8 g4.~ g4 e
  d4 c c d8 c
  e2 c'4 a
  g2 e8 c4.
  d4. c8 e4 d
  c1~
  c4 r4 r2
}

\addlyrics {
  "4" "-4" "5" "5" "-4" "5" "6" "5" "-4" "4" "4" "-4" "4" "5" "7" "-6" "6" "5" "4" "-4" "4" "5" "-4" "4"
  "4" "-4" "5" "5" "-4" "5" "6" "5" "-4" "4" "4" "-4" "4" "5" "7" "-6" "6" "5" "4" "-4" "4" "5" "-4" "4"
  "5" "-5" "6" "6" "6" "-5" "6" "-6" "6" "7" "-6" "6" "5" "-4" "4" "-4" "4" "-4" "5" "-4"
  "4" "-4" "5" "5" "-4" "5" "6" "5" "-4" "4" "4" "-4" "4" "5" "7" "-6" "6" "5" "4" "-4" "4" "5" "-4" "4"
}

Here is the output generated by Lilypond. Notice how the tabs and the notes are aligned to each other. You can click on the image to open the PDF file.

By the way, I can really recommend to sign up for David Barrett’s harmonica videos to learn tongue blocking which is essential to play single notes as in this example.

Enjoy!

Update:

I added harmonica tabs for Silent Night: