Expand noise section

10 files changed, 62 insertions, 12 deletions

diff --git a/default.nix b/default.nix
index b39ba49..fccf558 100644
--- a/default.nix
+++ b/default.nix
@@ -76,6 +76,7 @@ in pkgs.stdenv.mkDerivation rec {
       numpy
       Mako
       scipy
+      matplotlib
     ]);
 
   in with pkgs; [
diff --git a/lbm.org b/lbm.org
index 95fb004..1399d4f 100644
--- a/lbm.org
+++ b/lbm.org
@@ -4283,11 +4283,15 @@ less noticable for smaller step widths these are not desirable from a performanc
 Our renderer employs view-aligned slicing and random jittering to remove visible slicing.
 The choice of /randomness/ for jittering the ray origin is critical here as plain random numbers
 produce a ugly static-like pattern. A common choice in practice is to use so called /blue noise/
-instead.
+instead. Noise is called /blue/ if it contains only higher frequency components which makes it
+harder for the pattern recognizer that we call brain to find patterns where there should be none.
 
 For performance it makes sense to precompute such blue noise into tileable textures 
 than can then be easily used to fetch per-pixel-ray jitter offsets.
 
+The void-and-cluster algorithm cite:ulichneyVoidandclusterMethodDither1993
+provides a straight forward method for generating tileable blue noise textures.
+
 #+BEGIN_SRC python :session :results none
 import numpy as np
 import matplotlib.pyplot as plt
@@ -4295,21 +4299,37 @@ from numpy.ma import masked_array
 from scipy.ndimage import gaussian_filter
 #+END_SRC
 
+The first ingredient for this algorithm is a =filteredPattern= function that applies a
+plain Gaussian filter with given $\sigma$ to a cyclic 2d array. Using cyclic wrapping here is
+what makes the generated texture tileable.
+
 #+BEGIN_SRC python :session :results none
 def filteredPattern(pattern, sigma):
     return gaussian_filter(pattern.astype(float), sigma=sigma, mode='wrap', truncate=np.max(pattern.shape))
 #+END_SRC
 
+This function will be used to compute the locations of the largest void and tightest
+cluster in a binary pattern (i.e. a 2D array of 0s and 1s). In this context a /void/ describes
+an area with only zeros and a /cluster/ describes an area with only ones.
+
 #+BEGIN_SRC python :session :results none
 def largestVoidIndex(pattern, sigma):
     return np.argmin(masked_array(filteredPattern(pattern, sigma), mask=pattern))
 #+END_SRC
 
+These two functions work by considering the given binary pattern as a float array that is blurred by
+the Gaussian filter. The blurred pattern gives an implicit ordering of the /voidness/ of each pixel, the
+minimum of which we can determine by a simple search. It is important to exclude the initial binary
+pattern here as void-and-cluster depends on finding the largest areas where no pixel is set.
+
 #+BEGIN_SRC python :session :results none
 def tightestClusterIndex(pattern, sigma):
     return np.argmax(masked_array(filteredPattern(pattern, sigma), mask=np.logical_not(pattern)))
 #+END_SRC
 
+Computing the tightest cluster works in the same way with the exception of searching the largest array
+element and masking by the inverted pattern.
+
 #+BEGIN_SRC python :session :results none
 def initialPattern(shape, n_start, sigma):
     initial_pattern = np.zeros(shape, dtype=np.bool)
@@ -4326,6 +4346,11 @@ def initialPattern(shape, n_start, sigma):
     return initial_pattern
 #+END_SRC
 
+For the initial binary pattern we set =n_start= random locations to one and then repeatedly
+break up the largest void by setting its center to one. This is also done for the tightest cluster
+by setting its center to zero. We do this until the locations of the tightest cluster and largest
+void overlap.
+
 #+BEGIN_SRC python :session :results none
 def blueNoise(shape, sigma):
     n = np.prod(shape)
@@ -4354,8 +4379,18 @@ def blueNoise(shape, sigma):
     return noise / (n-1)
 #+END_SRC
 
+The actual algorithm utilizes these three helper functions in  four steps:
+1. Initial pattern generation
+2. Eliminiation of =n_start= tightest clusters
+3. Elimination of =n/2-n_start= largest voids
+4. Elimination of =n-n/2= tightest clusters of inverse pattern
+For each elimination the current =rank= is stored in the noise texture
+producing a 2D arrangement of the integers from 0 to =n=. As the last
+step the array is divided by =n-1= to yield a grayscale texture with values
+in $[0,1]$.
+
 #+BEGIN_SRC python :session :results file :cache yes
-test_noise = blueNoise((51,51), 1.5)
+test_noise = blueNoise((101,101), 1.9)
 test_noise_path = 'tangle/tmp/test_noise.png'
 
 fig, axs = plt.subplots(1,2, figsize=(8,4), sharey=True)
@@ -4370,7 +4405,7 @@ fig.savefig(test_noise_path, bbox_inches='tight', pad_inches=0)
 test_noise_path
 #+END_SRC
 
-#+RESULTS[0a8e452b001df3a7538378936788f3d942c7939c]:
+#+RESULTS[d74551fb6d155b85b3fce050ac0b5561ee217040]:
 [[file:tangle/tmp/test_noise.png]]
 
 #+NAME: blue-noise-values
@@ -4378,7 +4413,7 @@ test_noise_path
 fig, axs = plt.subplots(1,5, figsize=(8,2))
 
 for i in range(5):
-    noise = blueNoise((32,32), 1.5)
+    noise = blueNoise((32,32), 1.9)
     noise_rgb = np.dstack((noise,noise,noise))
     plt.imsave(f"tangle/asset/noise/blue_{ i }.png", noise_rgb)
     axs[i].imshow(noise_rgb)
@@ -4389,7 +4424,7 @@ fig.savefig(noise_overview_path, bbox_inches='tight', pad_inches=0)
 noise_overview_path
 #+END_SRC
 
-#+RESULTS[f7b10e2cb3a1e94192fb867be7572976f572cec5]: blue-noise-values
+#+RESULTS[124750b29bb542c10bf3a939d5226f502f0f98f9]: blue-noise-values
 [[file:tangle/tmp/noise_overview.png]]
 
 #+BEGIN_SRC cpp :tangle tangle/util/noise.h
@@ -4402,13 +4437,8 @@ noise_overview_path
 #include <SFML/Graphics.hpp>
 #+END_SRC
 
-#+BEGIN_SRC cpp :tangle tangle/LLBM/memory.h
-__device__ float noiseFromTexture(cudaSurfaceObject_t noisemap, int x, int y) {
-  uchar4 color{};
-  surf2Dread(&color, noisemap, x*sizeof(uchar4), y);
-  return color.x / 255.f;
-}
-#+END_SRC
+In order to manage the noise texture at runtime we implement a simple =NoiseSource=
+similar to the =ColorPalette= class.
 
 #+BEGIN_SRC cpp :tangle tangle/util/noise.h
 struct NoiseSource {
@@ -4425,6 +4455,19 @@ struct NoiseSource {
 };
 #+END_SRC
 
+The values of the selected noise texture are also accessed similarly.
+
+#+BEGIN_SRC cpp :tangle tangle/LLBM/memory.h
+__device__ float noiseFromTexture(cudaSurfaceObject_t noisemap, int x, int y) {
+  uchar4 color{};
+  surf2Dread(&color, noisemap, x*sizeof(uchar4), y);
+  return color.x / 255.f;
+}
+#+END_SRC
+
+The interaction UI only provides a image-displaying drop down box to choose
+between the five different noise textures generated in this section.
+
 #+BEGIN_SRC cpp :tangle tangle/util/noise.h
 void NoiseSource::interact() {
   ImGui::Image(texture, sf::Vector2f(32,20));
diff --git a/sources.bib b/sources.bib
index 17d8c4c..b9ab851 100644
--- a/sources.bib
+++ b/sources.bib
@@ -52,3 +52,9 @@
   journal = {IEEE Transactions on Visualization and Computer Graphics},
   number = {6},
 }
+@inproceedings{ulichneyVoidandclusterMethodDither1993,
+  title={Void-and-cluster method for dither array generation},
+  author={Ulichney, Robert},
+  booktitle={Electronic Imaging},
+  year={1993}
+}
\ No newline at end of file
diff --git a/tangle/asset/noise/blue_0.png b/tangle/asset/noise/blue_0.png
index d1c0534..1a033c7 100644
--- a/tangle/asset/noise/blue_0.png
+++ b/tangle/asset/noise/blue_0.png
diff --git a/tangle/asset/noise/blue_1.png b/tangle/asset/noise/blue_1.png
index c32ee28..bd3656a 100644
--- a/tangle/asset/noise/blue_1.png
+++ b/tangle/asset/noise/blue_1.png
diff --git a/tangle/asset/noise/blue_2.png b/tangle/asset/noise/blue_2.png
index 0db58de..2cc996d 100644
--- a/tangle/asset/noise/blue_2.png
+++ b/tangle/asset/noise/blue_2.png
diff --git a/tangle/asset/noise/blue_3.png b/tangle/asset/noise/blue_3.png
index e0d06cf..eed995a 100644
--- a/tangle/asset/noise/blue_3.png
+++ b/tangle/asset/noise/blue_3.png
diff --git a/tangle/asset/noise/blue_4.png b/tangle/asset/noise/blue_4.png
index 28d1b80..ef3f607 100644
--- a/tangle/asset/noise/blue_4.png
+++ b/tangle/asset/noise/blue_4.png
diff --git a/tangle/tmp/noise_overview.png b/tangle/tmp/noise_overview.png
index 2e657ac..8a3ea1e 100644
--- a/tangle/tmp/noise_overview.png
+++ b/tangle/tmp/noise_overview.png
diff --git a/tangle/tmp/test_noise.png b/tangle/tmp/test_noise.png
index e9056f7..b0490af 100644
--- a/tangle/tmp/test_noise.png
+++ b/tangle/tmp/test_noise.png