"Wave Collapse Function" algorithm in Processing

Dear all,

I would like to share with you a customized redesign of the Wave Collapse Function algorithm that works in Processing. My wish is that someone, someday, can port this sketch from Python to Java and make this algorithm available to the whole Processing community (more on that below).

What is the Wave Function Collapse algorithm ?

It is an algorithm written in 2016 by Maxim Gumin that can generate procedural patterns from a sample image or from a collection of tiles. You can see it in action here (2D “overlapping model”) and here (3D “tiled model”).

If you’re a game programmer / designer you’ve probably already heard of it as it has been used to generate terrain / buildings / cities proceduraly in a couple of recent games (see Bad North example).

Goal of this implementation:

To boil down the algorithm to its essence and avoid the redondancies and clumsiness of the original C# script (surprisingly long and difficult to read). This is an attempt to make a shorter, clearer and Processing compatible version of this algorithm, that can be understood and used by as many people as possible.

Characteristics of this implementation:

This is a rewritting of the 2D overlapping model in Processing Python mode. Unlike the original version, this implementation:

  • is not object oriented, making it easier to understand / closer to pseudo-code
  • is using 1D arrays instead of 2D arrays
  • is using array slicing for matrix manipulation
  • is populating the Wave with Python sets of indices whose size decrease as cells are “collapsed” (replacing large fixed size lists of booleans).
  • is storing data in dictionnaries whose keys are progressively deleted when possible.
  • is also storing data in tuples instead of lists
  • is skipping entropy computations (first entropy calculation not needed since equiprobable high level of uncertainty at start)
  • is displaying cells once (avoiding storing them in an array and redrawing at each frame)

For all these reasons, and also because it is Python (one-liner list comprehensions, minimalist syntax, …), this implementation is a lot, lot shorter than the original script (and its various ports) with only 70 lines of code against nearly 800 lines for the latter.
It is also faster in theory, but for obvious reasons (Jython under the hood) the following version won’t run as fast.

It should also be pointed that this implementation stores all the rotated and reflected matrices of each pattern and, like most versions, doesn’t use a backtracking system.

Full script

You’ll need to download a tiny bitmap picture (usually 8x8 or 16x16 pixels) beforehand to load it as an input. You can find a whole bunch on this dedicated github page.

from collections import Counter
from itertools import chain
from random import choice, sample

w, h, s = 96, 50, 9
N = 3

def setup():
    size(w*f, h*f, P2D)
    background('#FFFFFF')
    frameRate(1000)
    noStroke()
    
    global W, A, H, directions, patterns, freqs
    
    img = loadImage('Flowers.png') 
    iw, ih = img.width, img.height
    kernel = tuple(tuple(i + n*iw for i in xrange(N)) for n in xrange(N))
    directions = ((-1, 0), (1, 0), (0, -1), (0, 1))
    all = []
        
    for y in xrange(ih):
        for x in xrange(iw):
            cmat = tuple(tuple(img.pixels[((x+n)%iw)+(((a[0]+iw*y)/iw)%ih)*iw] for n in a) for a in kernel)
            for r in xrange(4):
                cmat = zip(*cmat[::-1])
                all.append(cmat)
                all.append(cmat[::-1])
                all.append([a[::-1] for a in cmat])

    all = [tuple(chain.from_iterable(p)) for p in all] 
    c = Counter(all)
    freqs = c.values()
    patterns = c.keys()
    npat = len(freqs) 
    
    W = dict(enumerate(tuple(set(range(npat)) for i in xrange(w*h))))
    A = dict(enumerate(tuple(set() for dir in xrange(len(directions))) for i in xrange(npat)))
    H = dict(enumerate(sample(tuple(npat if i > 0 else npat-1 for i in xrange(w*h)), w*h)))
    
    for i1 in xrange(npat):
        for i2 in xrange(npat):
            if [n for i, n in enumerate(patterns[i1]) if i%N!=(N-1)] == [n for i, n in enumerate(patterns[i2]) if i%N!=0]:
                A[i1][0].add(i2)
                A[i2][1].add(i1)
            if patterns[i1][:(N*N)-N] == patterns[i2][N:]:
                A[i1][2].add(i2)
                A[i2][3].add(i1)
           
          
def draw():    
    global H, W

    if not H:
        print 'finished'
        noLoop()
        return
        
    emin = min(H, key = H.get) 
    id = choice([idP for idP in W[emin] for i in xrange(freqs[idP])])
    W[emin] = {id}
    del H[emin]
    
    stack = {emin}
    while stack:
        idC = stack.pop() 
        for dir, t in enumerate(directions):
            x = (idC%w + t[0])%w
            y = (idC/w + t[1])%h
            idN = x + y * w 
            if idN in H: 
                possible = {n for idP in W[idC] for n in A[idP][dir]}
                if not W[idN].issubset(possible):
                    intersection = possible & W[idN] 
                
                    if not intersection:
                        print 'contradiction'
                        noLoop()
                        return
                        
                    W[idN] = intersection
                    H[idN] = len(W[idN]) - random(.1)
                    stack.add(idN)

    fill(patterns[id][0])
    rect((emin%w) * s, (emin/w) * s, s, s)
Fully annotated version
from collections import Counter
from itertools import chain
from random import choice, sample

w, h = 96, 50  # dimensions of output (array of wxh cells)
f = 9 # size factor 
N = 3 # dimensions of a pattern (NxN matrix)

def setup():
    size(w*f, h*f, P2D)
    background('#FFFFFF')
    frameRate(1000)
    noStroke()
    
    global W, A, H, directions, patterns, freqs, xs, ys
    
    img = loadImage('Flowers.png') # path to the input image
    iw, ih = img.width, img.height # dimensions of input image
    xs, ys = width//w, height//h # dimensions of cells (rect) in output
    kernel = tuple(tuple(i + n*iw for i in xrange(N)) for n in xrange(N)) # NxN matrix to read every patterns contained in input image
    directions = ((-1, 0), (1, 0), (0, -1), (0, 1)) # (x, y) tuples to access the 4 neighboring cells of a specific cell
    all = [] # array list to store all the patterns found in input
    
          
   
    #### Stores the different patterns found in input 
    
    for y in xrange(ih):
        for x in xrange(iw):
            
            ''' The one-liner below (cmat) creates a NxN matrix with (x, y) being its top left corner.
                This matrix will wrap around the edges of the input image.
                The whole snippet reads every NxN part of the input image and store the associated colors.
                Each NxN part is called a 'pattern' (of colors). Each pattern can be rotated or flipped (not mandatory). '''
                
            cmat = tuple(tuple(img.pixels[((x+n)%iw)+(((a[0]+iw*y)/iw)%ih)*iw] for n in a) for a in kernel)
            
            # Storing rotated patterns (90°, 180°, 270°, 360°)
            for r in xrange(4):
                cmat = zip(*cmat[::-1]) # +90° rotation
                all.append(cmat)
                all.append(cmat[::-1]) # vertical flip
                all.append([a[::-1] for a in cmat]) # horizontal flip 

    
    
    #### Flatten pattern matrices + count occurences 
    
    ''' Once every pattern has been stored,
        - we flatten them (convert to 1D) for convenience
        - count the number of occurences for each one of them (one pattern can be found multiple times in input)
        - select and store unique patterns only '''
    
    all = [tuple(chain.from_iterable(p)) for p in all]  # flattening arrays
    c = Counter(all) # Python counter
    freqs = c.values() # number of occurences for each unique pattern
    patterns = c.keys() # list of unique patterns
    npat = len(freqs) # number of unique patterns
    
    
    
    #### Initializes the 'wave' (W), entropy (H) and adjacencies (A) array lists
    
    ''' Array W (the Wave) keeps track of all the available patterns, for each cell.
        At start start, all patterns are valid anywhere in the Wave so each subarray
        is a list of indices of all the patterns'''
    
    W = dict(enumerate(tuple(set(range(npat)) for i in xrange(w*h)))) 
    
   
     ''' Array H should normally be populated with entropy values.
        Entropy is just a fancy way to represent the number of patterns 
        still available in a cell. We can skip this computation and populate 
        the array with the number of available patterns instead.
        
        At start all patterns are valid anywhere in the Wave, so all cells 
        share the same value (npat). We must however pick one cell at random and
        assign a lower value to it. Why ? Because the algorithm in draw() needs
        to find a cell with the minimum non-zero entropy value. '''

    H = dict(enumerate(sample(tuple(npat if i > 0 else npat-1 for i in xrange(w*h)), w*h))) 
    
    
    ''' Array A (for Adjacencies) is an index datastructure that describes the ways 
        that the patterns can be placed near one another. More explanations below '''
    
    A = dict(enumerate(tuple(set() for dir in xrange(len(directions))) for i in xrange(npat))) # explanations below
    
    
    
    #### Computation of patterns compatibilities (check if some patterns are adjacent, if so -> store them based on their location)
    
    ''' EXAMPLE:
        If pattern index 42 can placed to the right of pattern index 120,
        we will store this adjacency rule as follow:
    
                        A[120][1].add(42)
    
        Here '1' stands for 'right' or 'East'/'E'
    
        0 = left or West/W
        1 = right or East/E
        2 = up or North/N
        3 = down or South/S '''
        
    # Comparing patterns to each other
    for i1 in xrange(npat):
        for i2 in xrange(npat):
            
            ''' (in case when N = 3) If the first two columns of pattern 1 == the last two columns of pattern 2 
                 --> pattern 2 can be placed to the left (0) of pattern 1 '''
            
            if [n for i, n in enumerate(patterns[i1]) if i%N!=(N-1)] == [n for i, n in enumerate(patterns[i2]) if i%N!=0]:
                A[i1][0].add(i2)
                A[i2][1].add(i1)
                
            
            ''' (in case when N = 3) If the first two rows of pattern 1 == the last two rows of pattern 2
                --> pattern 2 can be placed on top (2) of pattern 1  '''
            
            if patterns[i1][:(N*N)-N] == patterns[i2][N:]:
                A[i1][2].add(i2)
                A[i2][3].add(i1)
           
          
def draw():    
    global H, W
    
    
    # Simple stopping mechanism
    
    ''' If the dict (or arraylist) of entropies is empty -> stop iterating.
        We'll see later that each time a cell is collapsed, its corresponding key
        in H is deleted '''
    
    if not H:
        print 'finished'
        noLoop()
        return
    
    
    
    #### OBSERVATION
    
    ''' Find cell with minimum non-zero entropy (not collapsed yet).'''
    
    emin = min(H, key = H.get) 
    
    
    
    #### COLLAPSE
    
    ''' Among the patterns available in the selected cell (the one with min entropy), 
        select one pattern randomly, weighted by the frequency that pattern appears 
        in the input image.'''
        
    id = choice([idP for idP in W[emin] for i in xrange(freqs[idP])]) # index of selected pattern 
    
    
    
    ''' The Wave's subarray corresponding to the cell with min entropy should now only contains 
        the id of the selected pattern '''
        
    W[emin] = {id} 
    
    
    
    ''' Its key can be deleted in the dict of entropies '''
        
    del H[emin] 
    
    
    
    #### PROPAGATION
    
    ''' Once a cell is collapsed, its index is put in a stack. 
        That stack is meant later to temporarily store indices of neighoring cells '''
        
    stack = {emin}
    
    
    
    ''' The propagation will last as long as that stack is filled with indices '''
    
    while stack:
        
    
    
    ''' First thing we do is pop() the last index contained in the stack (the only one for now) 
        and get the indices of its 4 neighboring cells (E, W, N, S). 
        We have to keep them withing bounds and make sure they wrap around. '''
        
        idC = stack.pop() # index of current cell
        for dir, t in enumerate(directions):
            x = (idC%w + t[0])%w
            y = (idC/w + t[1])%h
            idN = x + y * w # index of negihboring cell
            
            
            
    ''' We make sure the neighboring cell is not collapsed yet 
        (we don't want to update a cell that has only 1 pattern available) '''
            
            if idN in H: 
                
                
                
    ''' Then we check all the patterns that COULD be placed at that location. 
        EX: if the neighboring cell is on the left of the current cell (east side), 
        we look at all the patterns that can be placed on the left of each pattern 
        contained in the current cell. '''
        
                possible = {n for idP in W[idC] for n in A[idP][dir]}
                
              
                  
    ''' We also look at the patterns that ARE available in the neighboring cell '''
    
                available = W[idN]
                
                
                
    ''' Now we make sure that the neighboring cell really need to be updated. 
        If all its available patterns are already in the list of all the possible patterns:
         —> there’s no need to update it (the algorithm skip this neighbor and goes on to the next) '''
         
                if not available.issubset(possible):
                    
                    
                    
    ''' If it is not a subset of the possible list:
        —> we look at the intersection of the two sets (all the patterns that can be placed 
        at that location and that, "luckily", are available at that same location) '''             
                    
                    intersection = possible & available 
                
                
                
    ''' If they don't intersect (patterns that could have been placed there but are not available) 
        it means we ran into a "contradiction". We have to stop the whole WFC algorithm. '''
                   
                    if not intersection:
                        print 'contradiction'
                        noLoop()
                        return
                        
                        
                        
    ''' If, on the contrary, they do intersect -> we update the neighboring cell with that refined 
        list of pattern's indices '''
        
                    W[idN] = intersection
                    
                    
                    
    ''' Because that neighboring cell has been updated, its number of valid patterns has decreased
        and its entropy must be updated accordingly.
        Note that we're subtracting a small random value to mix things up: sometimes cells we'll
        end-up with the same minimum entropy value and this prevent to always select the first one of them.
        It's a cosmetic trick to break the monotony of the animation'''               
                    
                    H[idN] = len(W[idN]) - random(.1)
                    
                    
                    
    ''' Finally, and most importantly, we add the index of that neighboring cell to the stack 
        so it becomes the next current cell in turns (the one whose neighbors will be updated 
        during the next while loop) '''
        
                    stack.add(idN)



    #### RENDERING
    ''' The collapsed cell will always be filled with the first color (top left corner) of the
        selected pattern '''
        
    fill(patterns[id][0])
    rect((emin%w) * xs, (emin/w) * ys, xs, ys)

allopt

The need for a Processing Java version

The reasons are threefold:

A lot of Processing users like to create games. As mentionned above, the WFC algorithm can be very useful for game design: generating terrains, buildings or even entire cities, be it in 2D or 3D. This could be a welcomed contribution for any designer who want to speed-up and expand his creative process.

Long-awaited Coding Challenge. This algorithm has been submitted a couple of times on the Coding Train’s github (and mentionned multiple times on the YT channel) as an idea for a future Coding Challenge. It seems Daniel Shiffman really liked the suggestion but never had the chance to tackle it. Providing a clear and simple example sketch in Java could be a great opportunity to revive the submission and maybe have the chance to finally see the algorithm explained to a larger audience.

The generative nature of the algorithm is closely related to the philosophy of Processing.
In my opinion, Processing is mostly about generating things: design, art, drawings, colors, shapes… It is a tool to explore the vast field of possibilites that a simple idea can lead to.
So far, the WFC algorithm has primarly been used in the field of game design but I cannot help but to think there’s so much more to do with it. Some have started generating music with it, others poetry or even wine descriptions… and I believe Processing is the perfect environement to explore the possibilities even further. It is in its DNA.

Original Algorithm

The following explanations are here to help understanding the orginal algorithm, as it was designed by Maxim Gumin. It is not an exact representation of the Python implementation above but can be a good introductory step to the code.

1/ Read the input bitmap, store every NxN patterns and count their occurences. ( optional: Augment pattern data with rotations and reflections.)

For example, when N = 3:

enter image description here

2/ Precompute and store every possible adjacency relations between patterns. In the example below, patterns 207, 242, 182 and 125 can overlap the right side of pattern 246

enter image description here

3/ Create an array with the dimensions of the output (called W for wave). Each element of this array is an array holding the state ( True of False ) of each pattern.

For example, let’s say we count 326 unique patterns in input and we want our output to be of dimensions 20 by 20 (400 cells). Then the “Wave” array will contain 400 (20x20) arrays, each of them containing 326 boolan values.

At start, all booleans are set to True because every pattern is allowed at any position of the Wave.

W = [[True for pattern in xrange(len(patterns))] for cell in xrange(20*20)]

4/ Create another array with the dimensions of the output (called H ). Each element of this array is a float holding the “entropy” value of its corresponding cell in output.

Entropy here refers to Shannon Entropy and is computed based on the number of valid patterns at a specific location in the Wave. The more a cell has valid patterns (set to True in the Wave), the higher its entropy is.

For example, to compute the entropy of cell 22 we look at its corresponding index in the wave ( W[22] ) and count the number of booleans set to True . With that count we can now compute the entropy with the Shannon formula. The result of this calculation will be then stored in H at the same index H[22]

At start, all cells have the same entropy value (same float at every position in H ) since all patterns are set to True , for each cell.

H = [entropyValue for cell in xrange(20*20)]

These 4 steps are introductory steps, they are necessary to initalize the algorithm. Now starts the core of the algorithm:

5/ Observation:

Find the index of the cell with the minimum nonzero entropy (Note that at the very first iteration all entropies are equal so we need to pick the index of a cell randomly.)

Then, look at the still valid patterns at the corresponding index in the Wave and select one of them randomly, weighted by the frequency that pattern appears in the input image (weighted choice).

For example if the lowest value in H is at index 22 ( H[22] ), we look at all the patterns set to True at W[22] and pick one randomly based on the number of times it appears in the input. (Remember at step 1 we’ve counted the number of occurences for each pattern). This insures that patterns appear with a similar distribution in the output as are found in the input.

6/ Collapse:

We now assign the index of the selected pattern to the cell with the minimum entropy. Meaning that every pattern at the corresponding location in the Wave are set to False except for the one that has been chosen.

For example if pattern 246 in W[22] was set to True and has been selected, then all other patterns are set to False . Cell 22 is assigned pattern 246 . In output cell 22 will be filled with the first color (top left corner) of pattern 246. (blue in this example)

7/ Propagation:

Because of adjacency constraints, that pattern selection has consequences on the neighboring cells in the Wave. The arrays of booleans corresponding to the cells on the left and right, on top of and above the recently collapsed cell need to be updated accordingly.

For example if cell 22 has been collapsed and assigned with pattern 246 , then W[21] (left), W[23] (right), W[2] (up) and W[42] (down) have to be modified so as they only keep to True the patterns that are adjacent to pattern 246 .

For example, looking back at the picture of step 2, we can see that only patterns 207, 242, 182 and 125 can be placed on the right of pattern 246. That means that W[23] (right of cell 22 ) needs to keep patterns 207, 242, 182 and 125 as True and set all other patterns in the array as False . If these patterns are not valid anymore (already set to False because of a previous constraint) then the algorithm is facing a contradiction and the whole algorithm must be stopped.

This process doesn’t stop to the 4 direct neighbors of the collapsed cell and has to be extended recursively to the neighbors of the neighbors, and so on, until all constraints are propagated.

8/ Updating entropies

Because a cell has been collapsed (one pattern selected, set to True ) and all its surrounding cells updated accordingly (setting non adjacent patterns to False ), the entropy of these cells have changed and needs to be computed again. (Remember that the entropy of a cell is correlated to the number of valid pattern it holds in the Wave.)

In the example, the entropy of cell 22 is now 0, ( H[22] = 0 , because only pattern 246 is set to True at W[22] ) and the entropy of its neighboring cells have decreased (patterns that were not adjacent to pattern 246 have been set to False ).

By now the algorithm arrives at the end of the first iteration and will loop over steps 5 (find cell with minimum non zero entropy) to 8 (update entropies) until all cells are collapsed.

Useful additionnal ressources:

This detailed article from Stephen Sherratt and this explanatory paper from Karth & Smith. Also this comprehensive blogpost by Trasevol_Dog.

14 Likes

Please find below an implementation of the “tiled model”.

com-optimize(1)

The difference is that, unlike the “overlapping model” above, adjacency constraints are not inferred from a source image. As a result, constraints and frequencies (if needed) have to be entered manually (hard-coded).

A couple of improvements

Writing all the adjacencies of each tile can be cumbersome (especially with large tilesets) so I decided to add a connector-based system (inspired by this article by Marian42) to facilitate the task. Instead of storing all the possible adjacent tiles in each direction for each tile, it is now only required to associate each face of a tile with its corresponding symmetry index. All faces sharing the same symmetry index will be adjacent.

Also, because most tilesets don’t come with the rotated versions of each tile, I added a mechanism that rotates a tile and stores the transformated output based on its file name.

Here for example, tile “a” will be rotated 3 times at a 90-degree angle and each transformation will be added to the tiles array. “b” is not taken into account and “c” will be rotated only once. As a result the final array will contain 7 tiles:
Capture%20d%E2%80%99%C3%A9cran%20(134)

The mechanism also automatically creates and rotates the connectors (arrays of symmetry indices) of each new rotated tile.

Limits

In some cases you may want (for cosmetic reasons) to prevent a specific tile to be placed near another tile even though adjacent rules allow this placement. So far I’m writing these avoidance rules (as I call them) by hand but it would be really convenient (provided it is possible) to automate that part as much as possible.

Full script

from random import choice, sample
from collections import deque
from java.io import File

w, h, s = 80, 45, 10
directions = ((-1, 0), (1, 0), (0, -1), (0, 1))

def setup():
    size(w*s, h*s, P2D)
    background('#FFFFFF')
    frameRate(1000)
    
    global W, A, H, tiles
    
    path = File("\Knots").listFiles()
    connectors = [(0, 1, 1, 0), (0, 0, 0, 0), (1, 0, 1, 0)]
        
    tiles = []
    count = 0
    
    # Iterating through a folder of unique tiles
    for i, png in enumerate(path):
        
        # appending tile to 'tiles' array lit
        tiles.append(loadImage(png.getAbsolutePath()))
                
        # if rotation information included in tile name --> allow rotation
        name = png.getName().split('.')[0]
        if '-' in name:

            # rotate by 90° the number of times specified in tile name 
            for rot in xrange(int(name[-1])):
                cp = tiles[-1]
                pg = createGraphics(cp.width, cp.height)
                pg.beginDraw()
                pg.pushMatrix()
                pg.translate(cp.width, 0)
                pg.rotate(HALF_PI)
                pg.image(cp, 0, 0)
                pg.popMatrix()
                pg.endDraw()
                
                # add rotated tile to 'tiles' array list
                tiles.append(pg)
  
                # create a corresponding "connector" array by rotating the indices of the previous tile
                c = deque(connectors[i + count])
                c.rotate(1)
                count += 1
                connectors.insert(i + count, tuple(c))
                
            
    ntiles = len(tiles) 

    W = dict(enumerate(tuple(set(range(ntiles)) for i in xrange(w*h)))) 
    H = dict(enumerate(sample(tuple(ntiles if i > 0 else ntiles-1 for i in xrange(w*h)), w*h)))
    A = dict(enumerate([[set() for dir in xrange(4)] for i in xrange(ntiles)]))
    
    for i1 in xrange(ntiles):
        for i2 in xrange(ntiles):
            if connectors[i1][0] == connectors[i2][2]:
                A[i1][0].add(i2)
                A[i2][1].add(i1)
                
            if connectors[i1][1] == connectors[i2][3]: 
                A[i1][2].add(i2)
                A[i2][3].add(i1)
    
                
def draw():    
    global H, W
    
    if not H:
        print 'finished'
        noLoop()
        return

    emin = min(H, key = H.get) 
    id = choice(list(W[emin]))
    W[emin] = {id}
    del H[emin]
    
    stack = {emin}
    while stack:
        idC = stack.pop() 
        for dir, t in enumerate(directions):
            x = (idC%w + t[0])%w
            y = (idC/w + t[1])%h
            idN = x + y * w 
            if idN in H: 
                possible = {n for idP in W[idC] for n in A[idP][dir]}
                if not W[idN].issubset(possible):
                    intersection = possible & W[idN] 
                
                    if not intersection:
                        print 'contradiction'
                        noLoop()
                        break

                    W[idN] = intersection
                    H[idN] = len(W[idN]) - random(.1)
                    stack.add(idN)
   
    image(tiles[id], (emin%w) * s, (emin/w) * s, s, s)
6 Likes

Great approach, good explanation, interesting topic. One of the best posts I ever read here. I can’t believe no one answered yet. It’s a perfect match for processing. If you understand the algorithm, it should be easy for you to translate it into Java code, even if you don’t know the language yet. It derives from the C syntax family, a concept you can hardly avoid when dealing with programming in any way.

2 Likes

Well, somebody could translate but some is just gibberish…

:wink:

@needfulthing: Thanks for your feedback and the kind words !

If I had more time I could possibly, and laboriously, come up with a Java version but it would certainly be clumsy and inneficient. I would rather see someone with a good understanding of Java data types implement it in the best possible way.

An example: I chose Python sets() because they are efficient and easy to manipulate. I could simply look for the Java equivalent but:

  • Would they be as efficient ? (compared to other Java data types)
  • Would they be as easy to manipulate ? (it seems they’re not)
  • Is there a better alternative ?
  • Could I compute intersections with it ?

Only someone with a good experience of the Java language could answer these questions.

@Chrisir: Thank you for your interest !

Have you looked at the annotations I provided in the second script (just above the first gif) ?

''' Then we check all the patterns that COULD be placed at that location. 
    EX: if the neighboring cell is on the left of the current cell (east side), 
    we look at all the patterns that can be placed on the left of each pattern 
    contained in the current cell. '''

    possible = {n for idP in W[idC] for n in A[idP][dir]}

Breaking that line down:

  • brackets {} indicate we’re dealing with a Python set(). I chose this data type not only because it is efficient but also because we need later to compute the intersection of 2 sets (intersection only works with sets, not operable with lists)
  • idC is the index of the Current cell (the one whose neighbors are updated)
  • idP is the index of the available Pattern in that cell
  • n is the index of the pattern than can be adjacent to pattern idP in direction dir

Basically we look at all the indices of the patterns available in the current cell, then for each one of them we look at all the indices of patterns that can be placed next to it, in a specific direction.
Avoiding that one-liner, it would give something like this:

possible = []
for availablePattern in currentCell:
    for adjacentPattern in listOfAdjacentPatternsToavailablePatternInThisDirection:
        possible.append(adjacentPattern)

2 Likes

Hi solub – this is a wonderful contribution, and I’m a huge fan of this work – and your effort in implementing it and documenting it so clearly. Fantastic.

A huge quick question – is the entropy array named “H” by mathematical convention?

I’m tempted to post a slight revision of your annotated code, and give all the variables semantic names – e.g. W --> wave, A --> adjacencies, H --> entropy, etc. – to make the code more readable for students and potential adaptors. The extreme concision is very elegant, but it might be a barrier to comprehension for coders without a background in this kind of math.

2 Likes

Given these sketches got some image assets, it’d be convenient to host them (including their root project folder) on GitHub. :octopus:

1 Like

Unfortunately that annotated code won’t run in PDE due to the comment-indent style – it causes parse errors.

I’ve cleaned up solub’s code comments and then renamed most of the variables to make them semantic. So for example:

possible = {n for idP in W[idC] for n in A[idP][dir]}

is now written:

possible = {n for pattern_idx in wave[cell_idx] for n in adjacencies[pattern_idx][dir]}

…which is still nested loops written as a compact generator expression, but it is a bit more clear what is happening.

If someone is interested in doing a Java port, this would be one place to start. Code below:

from collections import Counter
from itertools import chain
from random import choice, sample

out_w, out_h = 96, 50  # dimensions of output (array of wxh cells)
f = 9 # size factor
N = 3 # dimensions of a pattern (NxN matrix)

def setup():
    size(out_w*f, out_h*f, P2D)
    background('#FFFFFF')
    frameRate(1000)
    noStroke()

    global wave, adjacencies, entropy, directions, patterns, freqs, cell_w, cell_h

    img = loadImage('Flowers.png') # path to the input image
    img_w, img_h = img.width, img.height # dimensions of input image
    cell_w, cell_h = width//out_w, height//out_h # dimensions of cells (rect) in output
    kernel = tuple(tuple(i + n*img_w for i in xrange(N)) for n in xrange(N)) # NxN matrix to read every patterns contained in input image
    directions = ((-1, 0), (1, 0), (0, -1), (0, 1)) # (x, y) tuples to access the 4 neighboring cells of a specific cell
    all = [] # array list to store all the patterns found in input

    #### Stores the different patterns found in input

    for y in xrange(img_h):
        for x in xrange(img_w):
            
            ''' The one-liner below (cmat) creates a NxN matrix with (x, y)
            being its top left corner. This matrix will wrap around the edges
            of the input image. The whole snippet reads every NxN part of the
            input image and store the associated colors. Each NxN part is
            called a 'pattern' (of colors). Each pattern can be rotated or
            flipped (not mandatory).
            '''
            
            cmat = tuple(tuple(img.pixels[((x+n)%img_w)+(((a[0]+img_w*y)/img_w)%img_h)*img_w] for n in a) for a in kernel)
            
            # Storing rotated patterns (90°, 180°, 270°, 360°)
            for r in xrange(4):
                cmat = zip(*cmat[::-1]) # +90° rotation
                all.append(cmat)
                all.append(cmat[::-1]) # vertical flip
                all.append([a[::-1] for a in cmat]) # horizontal flip

    #### Flatten pattern matrices + count occurences

    ''' Once every pattern has been stored,
    - we flatten them (convert to 1D) for convenience
    - count the number of occurences for each one of them (one pattern can be 
      found multiple times in input)
    - select and store unique patterns only
    '''

    all = [tuple(chain.from_iterable(p)) for p in all]  # flattening arrays
    c = Counter(all) # Python counter
    freqs = c.values() # number of occurences for each unique pattern
    patterns = c.keys() # list of unique patterns
    npat = len(freqs) # number of unique patterns

    #### Initializes the 'wave' (wave), entropy (entropy) and adjacencies (adjacencies) array lists

    ''' Array wave (the Wave) keeps track of all the available patterns, for
    each cell. At start start, all patterns are valid anywhere in the Wave so 
    each subarray is a list of indices of all the patterns
    '''

    wave = dict(enumerate(tuple(set(range(npat)) for i in xrange(out_w*out_h))))

    ''' Array entropy should normally be populated with entropy values.
    Entropy is just a fancy way to represent the number of patterns
    still available in a cell. We can skip this computation and populate
    the array with the number of available patterns instead.

    At start all patterns are valid anywhere in the Wave, so all cells
    share the same value (npat). We must however pick one cell at random and
    assign a lower value to it. Why ? Because the algorithm in draw() needs
    to find a cell with the minimum non-zero entropy value.
    '''

    entropy = dict(enumerate(sample(tuple(npat if i > 0 else npat-1 for i in xrange(out_w*out_h)), out_w*out_h)))

    ''' Array adjacencies (for Adjacencies) is an index datastructure that
    describes the ways
    that the patterns can be placed near one another. More explanations below
    '''

    adjacencies = dict(enumerate(tuple(set() for dir in xrange(len(directions))) for i in xrange(npat))) # explanations below

    #### Computation of patterns compatibilities (check if some patterns are adjacent, if so -> store them based on their location)

    ''' EXAMPLE:
    If pattern index 42 can placed to the right of pattern index 120,
    we will store this adjacency rule as follow:

        adjacencies[120][1].add(42)

    Here '1' stands for 'right' or 'East'/'E'

    0 = left or West/W
    1 = right or East/E
    2 = up or North/N
    3 = down or South/S
    '''
    
    # Comparing patterns to each other
    for i1 in xrange(npat):
        for i2 in xrange(npat):
            
            '''(in case when N = 3) If the first two columns of
            pattern 1 == the last two columns of pattern 2
            --> pattern 2 can be placed to the left (0) of pattern 1
            '''

            if [n for i, n in enumerate(patterns[i1]) if i%N!=(N-1)] == [n for i, n in enumerate(patterns[i2]) if i%N!=0]:
                adjacencies[i1][0].add(i2)
                adjacencies[i2][1].add(i1)

            '''(in case when N = 3) If the first two rows of
            pattern 1 == the last two rows of pattern 2
            --> pattern 2 can be placed on top (2) of pattern 1
            '''

            if patterns[i1][:(N*N)-N] == patterns[i2][N:]:
                adjacencies[i1][2].add(i2)
                adjacencies[i2][3].add(i1)

def draw():
    global entropy, wave

    # Simple stopping mechanism

    '''If the dict (or arraylist) of entropies is empty -> stop iterating.
    We'll see later that each time a cell is collapsed, its corresponding key
    in entropy is deleted
    '''

    if not entropy:
        print 'finished'
        noLoop()
        return

    #### OBSERVATION

    ''' Find cell with minimum non-zero entropy (not collapsed yet).'''

    entropy_min = min(entropy, key = entropy.get)

    #### COLLAPSE

    ''' Among the patterns available in the selected cell (the one with min entropy),
        select one pattern randomly, weighted by the frequency that pattern appears
        in the input image.'''

    pattern_id = choice([pattern_idx for pattern_idx in wave[entropy_min] for i in xrange(freqs[pattern_idx])]) # index of selected pattern

    ''' The Wave's subarray corresponding to the cell with min entropy should now only contains
        the id of the selected pattern '''

    wave[entropy_min] = {pattern_id}

    ''' Its key can be deleted in the dict of entropies '''

    del entropy[entropy_min]

    #### PROPAGATION

    ''' Once a cell is collapsed, its index is put in a stack.
        That stack is meant later to temporarily store indices of neighoring cells '''

    stack = {entropy_min}

    ''' The propagation will last as long as that stack is filled with indices '''

    while stack:

        ''' First thing we do is pop() the last index contained in the stack (the only one for now)
        and get the indices of its 4 neighboring cells (E, W, N, S).
        We have to keep them withing bounds and make sure they wrap around. '''

        cell_idx = stack.pop() # index of current cell
        for dir, t in enumerate(directions):
            x = (cell_idx%out_w + t[0])%out_w
            y = (cell_idx/out_w + t[1])%out_h
            neighbor_idx = x + y * out_w # index of negihboring cell

            ''' We make sure the neighboring cell is not collapsed yet
            (we don't want to update a cell that has only 1 pattern available)
            '''

            if neighbor_idx in entropy:

                ''' Check all patterns that COULD be placed at that location.
                EX: if the neighboring cell is on the left of the current cell
                (east side), we look at all the patterns that can be placed on
                the left of each pattern contained in the current cell.
                '''

                possible = {n for pattern_idx in wave[cell_idx] for n in adjacencies[pattern_idx][dir]}

                ''' We also look at the patterns that ARE available in the
                neighboring cell
                '''

                available = wave[neighbor_idx]

                ''' Now we make sure that the neighboring cell really need to
                be updated. If all its available patterns are already in the
                list of all the possible patterns: —> there’s no need to update
                it (the algorithm skip this neighbor and goes on to the next)
                '''

                if not available.issubset(possible):

                    ''' If it is not a subset of the possible list:
                    —> we look at the intersection of the two sets
                    (all the patterns that can be placed at that location
                    and that, "luckily", are available at that same location)
                    '''

                    intersection = possible & available

                    ''' If they don't intersect (patterns that could have been
                    placed there but are not available) it means we ran into a
                    "contradiction". We have to stop the whole WFC algorithm.
                    '''

                    if not intersection:
                        print 'contradiction'
                        noLoop()
                        return

                    ''' If, on the contrary, they do intersect -> we update the
                    neighboring cell with that refined list of pattern's indices
                    '''

                    wave[neighbor_idx] = intersection

                    ''' Because that neighboring cell has been updated, its number of valid patterns has decreased
                    and its entropy must be updated accordingly.
                    Note that we're subtracting a small random value to mix
                    things up: sometimes cells we'll end-up with the same
                    minimum entropy value and this prevent to always select the
                    first one of them. It's a cosmetic trick to break the
                    monotony of the animation'''

                    entropy[neighbor_idx] = len(wave[neighbor_idx]) - random(.1)

                    ''' Finally, and most importantly, we add the index of that
                    neighboring cell to the stack so it becomes the next current
                    cell in turns (the one whose neighbors will be updated
                    during the next while loop)
                    '''

                    stack.add(neighbor_idx)

    #### RENDERING
    ''' The collapsed cell will always be filled with the first color (top left
    corner) of the selected pattern '''
    
    fill(patterns[pattern_id][0])
    rect((entropy_min%out_w) * cell_w, (entropy_min/out_w) * cell_h, cell_w, cell_h)
1 Like

@jeremydouglass Thank you, really glad to hear your feedback !

Yes, I named the entropy array “H” by convention even though it is the entropy itself that should named “H”, not the array. Thinking about it now, “E” would be fine as well, at least it would ensure greater consistency between the arrays: W --> wave, A --> adjacency, E --> entropy.

Totally agree. Thank you so much for your revised version. It’s beyond me why the annotated script doesn’t run on your side, it works perfectly fine on my computer…

@GoToLoop

In this regard I would also add that the overlapping model is more suitable for students, mainly for 3 reasons:

  • tilesets are not easy to get a hold of / create
  • it is easier to manage a single bitmap image than a whole set of tiles (no rotations to compute and store, no adjacencies or frequencies to enter by hand)
  • bitmap images can easily be created in Processing (downloading an example image is not mandatory)

Actually, the output that looks like an archipelago (first gif) was based on a input image made in Processing. Alternatively, this could be a great incentive for students / hobbysts to create their own bitmaps while learning the basics of Processing

Example

N = 16 # 16 x 16 pixels

def setup():
    size(N, N)
    
    # colors
    w = color(242, 242, 240) # white
    b = color(124, 215, 235) # blue
    g = color(66, 217, 134) # green
    B = color(0) # black
    c = color(198, 149, 99) # chestnut brown
    
    
    plant = [w, w, w, w, w, w, B, B, B, w, w, w, w, w, w, w,
             w, w, w, w, w, B, g, g, g, B, w, w, w, w, w, w,
             w, w, B, B, w, B, g, B, g, B, w, w, B, B, w, w,
             w, B, g, g, B, B, g, g, g, B, w, B, g, g, B, w,
             w, B, g, B, g, B, g, B, g, B, B, g, B, g, B, w,
             w, B, g, g, g, B, g, g, g, B, g, g, g, g, B, w,
             w, w, B, g, g, g, B, g, B, g, B, g, g, B, w, w,
             w, w, w, B, g, g, B, g, B, g, g, g, B, w, w, w,
             w, B, B, B, B, B, B, B, B, B, B, B, B, B, B, w,
             B, b, b, b, b, b, b, b, b, b, b, b, b, b, b, B,
             B, b, b, b, b, b, b, b, b, b, b, b, b, b, b, B,
             B, b, b, b, b, b, b, b, b, b, b, b, b, b, b, B,
             w, B, B, b, b, b, b, b, b, b, b, b, b, B, B, w,
             w, w, w, B, B, B, B, B, B, B, B, B, B, w, w, w,
             w, w, B, b, b, b, b, b, b, b, b, b, b, B, w, w,
             w, w, B, B, B, B, B, B, B, B, B, B, B, B, w, w]
    
    
    
    island = [b, b, b, b, b, b, b, b, b, b, b, b, b, b, b, b,
              b, b, b, b, b, b, b, b, b, b, g, g, b, b, b, b,
              b, b, b, b, b, b, b, b, g, g, g, g, g, b, b, b,
              b, b, b, b, b, b, b, b, g, g, g, g, g, g, b, b,
              b, b, b, b, b, g, g, g, g, g, g, g, g, g, b, b,
              b, b, b, b, g, g, g, g, g, g, g, g, g, g, b, b,
              b, b, g, g, g, g, g, g, g, g, g, g, g, g, g, b,
              b, b, g, g, g, g, g, g, g, g, g, g, g, g, g, b,
              b, b, c, g, g, g, g, g, g, g, g, g, g, g, g, b,
              b, b, b, c, g, g, g, g, g, g, g, g, g, g, c, b,
              b, b, b, b, g, g, g, g, g, g, g, g, c, c, b, b,
              b, b, b, g, g, g, g, g, g, g, g, g, b, b, b, b,
              b, b, b, c, g, g, c, c, g, g, g, g, b, b, b, b,
              b, b, b, b, c, c, b, b, c, c, c, c, b, b, b, b,
              b, b, b, b, b, b, b, b, b, b, b, b, b, b, b, b,
              b, b, b, b, b, b, b, b, b, b, b, b, b, b, b, b,]
    
    
    
    loadPixels()
    for i in xrange(width*height):
        pixels[i] = plant[i]
    updatePixels()
    
    
    save('C:\Path\to\Folder\name.png')

1 Like

Very interesting stuff, @solub. Thank you for sharing this.

This fit exactly with a personal project and I finally sat down long enough to work it out.

Notes:

  • Compared to Python, Java is… verbose.
  • This implementation is object oriented. I did so to keep things straight in my mind.
  • It also does not handle a “contradiction” gracefully. The sequence will flood the remaining area with “null” values, which in this case, turn into yellow pixels.

I am (:calendar: :astonished: quite) late to the party but I hope someone will find this helpful.

While thinking about this algorithm, I realized it can be used to create word-search tables.
Being able to specify a few words and have the program fill in the rest is basically built-in.

A very simple word-search-builder is on Github as well.



import java.util.Arrays;
import java.util.HashSet;
import java.util.Set;
import java.util.TreeMap;
import java.util.Stack;

// Patterns will be square
int S = 3;
int CUT_SIZE = S-1;
// Helper names
int SOUTH=0, NORTH=1, WEST=2, EAST=3;
// Keep the program updating even with
//   a lot of stack usage
int MAX_ITERATIONS = 100;

PGraphics orig;
PGraphics output;
PatternLibrary archive;
Wave wave;

boolean running = true;
boolean keepRecording = true;
int saveCounter = 0;

class Pattern
{
  int[][] contents;
  
  Pattern(int[][] in)
  {
    contents = in;
  }
  
  public int hashCode()
  {
    int key = 0;
    for (int y=0; y<S; y++)
    {
      for (int x=0; x<S; x++)
      {
        int pix = contents[y][x];
        key = (31*key) + ((pix>>24) & 0xFF);
        key = (31*key) + ((pix>>8) & 0xFF);
        key = (31*key) + (pix & 0xFF);
      }
    }
    return key;
  }
  
  Pattern rotNxNinetyDeg(int reps)
  {
    int[][] rota = new int[S][S];
    float offset = (S-1) / 2.0;
    for (int y=0; y<S; y++)
    {
      for (int x=0; x<S; x++)
      {
        float nx = x - offset;
        float ny = y - offset;
        for (int i=0; i<reps; i++)
        {
          float tx = nx;
          nx = -ny;
          ny =  tx;
        }
        nx = nx + offset;
        ny = ny + offset;
        rota[y][x] = contents[(int)nx][(int)ny];
      }
    }
    return new Pattern(rota);
  }
  
  Pattern flipH()
  {
    int[][] f = new int[S][S];
    for (int y=0; y<S; y++)
    {
      for (int x=0; x<S; x++)
      {
        f[y][x] = contents[y][S-x-1];
      }
    }
    return new Pattern(f);
  }
  
  Pattern flipV()
  {
    int[][] f = new int[S][S];
    for (int y=0; y<S; y++)
    {
      for (int x=0; x<S; x++)
      {
        f[y][x] = contents[S-y-1][x];
      }
    }
    return new Pattern(f);
  }
  
  boolean canBeAbove(Pattern other) { return compareTopToBottom(contents, other.contents); }
  boolean canBeBelow(Pattern other) { return compareTopToBottom(other.contents, contents); }
  boolean canBeToTheRightOf(Pattern other) { return compareLeftToRight(other.contents, contents); }
  boolean canBeToTheLeftOf(Pattern other) { return compareLeftToRight(contents, other.contents); }
  
  boolean compareTopToBottom(int[][] top, int[][] bottom)
  {
    boolean same = true;
    for (int y=0; y<CUT_SIZE; y++)
    {
      for (int x=0; x<S; x++)
      {
        int tp = top[y+1][x];
        int bt = bottom[y][x];
        
        if (tp != bt)
        {
          same = false;
          break; // get out early when possible
        }
      }
    }
    return same;
  }
  
  boolean compareLeftToRight(int[][] left, int[][] right)
  {
    boolean same = true;
    for (int y=0; y<S; y++)
    {
      for (int x=0; x<CUT_SIZE; x++)
      {
        int lf = left[y][x+1];
        int rt = right[y][x];
        
        if (lf != rt)
        {
          same = false;
          break; // get out early when possible
        }
      }
    }
    return same;
  }
}

class PatternLibrary
{
  int[][][] adjacencies = new int[4][][]; // four cardinal directions
  Pattern[] patterns;
  int[] elements;
  int[] frequencies;
  int totalTiles = 1;
  
  PatternLibrary(PImage template)
  {
    totalTiles = gatherPatterns(template);
    for (int dir=0; dir<4; dir++)
    {
      adjacencies[dir] = new int[totalTiles][totalTiles];
    }
    findAdjacencies();
    println("Total tiles: " + totalTiles);
  }
  
  void placePattern(Pattern p, HashMap<Integer,Pattern> tiles, HashMap<Integer,Integer> freqs)
  {
    int k = p.hashCode();
    tiles.put(k, p);
    freqs.put(k, freqs.getOrDefault(k,0) + 1);
  }
  
  int gatherPatterns(PImage src)
  {
    HashMap<Integer,Integer> freqs = new HashMap();
    HashMap<Integer,Pattern> tiles = new HashMap();
    println("gathering patterns...");
    for (int y=0; y<src.height+S; y++)
    {
      for (int x=0; x<src.width+S; x++)
      {
        int[][] ptrn = new int[S][S];
        for (int dy=0; dy<S; dy++)
        {
          for (int dx=0; dx<S; dx++)
          {
            int tx = (x+dx) % src.width;
            int ty = (y+dy) % src.height;
            int idx = ty*src.width + tx;
            ptrn[dy][dx] = src.pixels[idx];
          }
        }
        
        Pattern block = new Pattern(ptrn);
        Pattern rota = block.rotNxNinetyDeg(1);
        Pattern rotb = block.rotNxNinetyDeg(2);
        Pattern rotc = block.rotNxNinetyDeg(3);
        placePattern(block, tiles, freqs);
        placePattern(block.flipH(), tiles, freqs);
        placePattern(block.flipV(), tiles, freqs);
        placePattern(rota, tiles, freqs);
        placePattern(rotb, tiles, freqs);
        placePattern(rotc, tiles, freqs);
        
        // Optional extra patterns
        placePattern(rota.flipH(), tiles, freqs);
        placePattern(rota.flipV(), tiles, freqs);
        placePattern(rotb.flipH(), tiles, freqs);
        placePattern(rotb.flipV(), tiles, freqs);
        placePattern(rotc.flipH(), tiles, freqs);
        placePattern(rotc.flipV(), tiles, freqs);
      }
    }
    
    println("Collecting elements...");
    Set<Integer> keys = tiles.keySet();
    patterns = new Pattern[keys.size()];
    elements = new int[patterns.length];
    frequencies = new int[patterns.length];
    int counter = 0;
    for (Integer key : keys)
    {
      Pattern selected = tiles.get(key);
      patterns[counter] = selected;
      elements[counter] = selected.contents[0][0];
      frequencies[counter] = freqs.get(key);
      counter++;
    }
    return patterns.length;
  }
  
  void markAllowed(int keyA, int dirA, int keyB, int dirB)
  {
    adjacencies[dirA][keyA][keyB] = 1;
    adjacencies[dirB][keyB][keyA] = 1;
  }
  
  void findAdjacencies()
  {
    println("building adjacency tables...");
    for (int key=0; key<patterns.length; key++)
    {
      Pattern selected = patterns[key];
      for (int otherIndex=0; otherIndex<patterns.length; otherIndex++)
      {
        Pattern other = patterns[otherIndex];
        if ( selected.canBeAbove(other) ) markAllowed(key, SOUTH, otherIndex, NORTH);
        if ( selected.canBeBelow(other) ) markAllowed(key, NORTH, otherIndex, SOUTH);
        if ( selected.canBeToTheRightOf(other) ) markAllowed(key, WEST, otherIndex, EAST);
        if ( selected.canBeToTheLeftOf(other) ) markAllowed(key, EAST, otherIndex, WEST);
      }
    }
  }
}

class Field
{
  int[] state;
  int[] allowBuffer; // Used to prevent object creation during main algorithm
  int tileType = -1;
  float entropy = 1.0;
  int count = 0;
  boolean needsUpdate = true;
  
  Field(int[] counts)
  {
    state = Arrays.copyOf(counts, counts.length);
    allowBuffer = new int[state.length];
    count = state.length;
    entropy = count; // all are possible initially
    needsUpdate = true;
  }
  
  void restrict(int[] allowed)
  {
    count = 0;
    for (int i=0; i<state.length; i++)
    {
      state[i] *= allowed[i];
      if (state[i] != 0) count++;
    }
    entropy = count - random(0.2);
  }
  
  void gatherNeighborhood(int[][] lookingDirection)
  {
    Arrays.fill(allowBuffer, 0);
    for (int i=0; i<state.length; i++)
    {
      if (state[i] != 0)
      {
        int[] mult = lookingDirection[i];
        for (int j=0; j<state.length; j++)
        {
          allowBuffer[j] |= mult[j];
        }
      }
    }
  }
  
  int[] getAllowed(int[][] lookingDirection)
  {
    if (count == 1)
    {
      if (tileType != -1) return lookingDirection[tileType]; // attempt shortcut
      else gatherNeighborhood(lookingDirection); // fall back to regular check
    }
    else
    {
      gatherNeighborhood(lookingDirection);
    }
    return allowBuffer;
  }
  
  int weightedChoice()
  {
    int len = 0;
    TreeMap tm = new TreeMap<Float,Integer>();
    float total = 0;
    for (int i=0; i<state.length; i++)
    {
      if (state[i] > 0)
      {
        total += state[i];
        tm.put(total, i);
        len++;
      }
    }
    if (len > 0)
    {
      float select = random(total);
      if (tm.higherEntry(select) != null)
      {
        int v = (Integer)tm.higherEntry(select).getValue();
        Arrays.fill(state, 0);
        state[v] = 1;
        count = 1;
        return v;
      }
      return -1;
    }
    return -1;
  }
  
  boolean collapse()
  {
    if (count < 1) return false; // contradiction
    tileType = weightedChoice();
    if (-1 == this.tileType) return false; // contradiction
    entropy = 0.0;
    count = 1;
    needsUpdate = true;
    return true;
  }
}

class Wave
{
  int superPosition = color(0,88,201);
  Stack<Integer> todo;
  
  int w = 1;
  int h = 1;
  boolean initial = true;
  boolean isStable = true;
  
  Field[] area; // flattened 2D array
  PatternLibrary lib;
  
  Wave(int wWidth, int wHeight, PatternLibrary pl)
  {
    this.w = wWidth;
    this.h = wHeight;
    this.area = new Field[this.w*this.h];
    this.lib = pl;
    this.reset();
    this.todo = new Stack();
    println("WxH: " + this.w + "x" + this.h);
    println("Patterns: " + lib.totalTiles);
  }
  
  int findLowestEntropyNotEqualToZero()
  {
    int idx = -1;
    float mx = 9999999;
    for (int key=0; key<area.length; key++)
    {
      Field f = area[key];
      float e = f.entropy;
      if ((e < mx) && (e > 0) && (-1 == f.tileType))
      {
        idx = key;
        mx = e;
      }
    }
    return idx;
  }
  
  void step()
  {
    if (isStable)
    {
      observe();
    }
    isStable = propagate();
  }
  
  void observe()
  {
    int idx = -1;
    if (initial) // random selection
    {
      int rx = (int)random(w);
      int ry = (int)random(h);
      idx = (ry * w) + rx;
      initial = false;
    }
    else
    {
      idx = findLowestEntropyNotEqualToZero();
    }
    if (idx != -1)
    {
      Field f = area[idx];
      if (f.collapse()) todo.push(idx); // do not use contradiction for remainder
    }
  }
  
  void adjustNeighbor(int n, int dir, Field ref)
  {
    Field nei = area[n];
    int prevCount = nei.count;
    if (-1 != nei.tileType) return; // prevent an error from wiping out the entire design
    nei.restrict( ref.getAllowed(lib.adjacencies[dir]) );
    if (nei.count != prevCount)
    {
      nei.needsUpdate = true;
      todo.push(n);
    }
  }
  
  boolean propagate()
  {
    int workLoad = 0;
    boolean done = todo.empty();
    while (!done && (workLoad < MAX_ITERATIONS))
    {
      int i = todo.pop();
      int x = i % w;
      int y = i / w;
      int dpx = (x + w + 1) % w;
      int dnx = (x + w - 1) % w;
      int dpy = (y + h + 1) % h;
      int dny = (y + h - 1) % h;
      int id_W = (y * w) + dnx;
      int id_E = (y * w) + dpx;
      int id_N = (dny * w) + x;
      int id_S = (dpy * w) + x;
      
      Field f = area[i];
      // Method:
      //   update each neighbor based on the contents of selected
      //   if neighbor changed, push it onto the stack
      adjustNeighbor(id_N, NORTH, f);
      adjustNeighbor(id_S, SOUTH, f);
      adjustNeighbor(id_E, EAST, f);
      adjustNeighbor(id_W, WEST, f);
      
      workLoad++;
      done = todo.empty();
    }
    return done;
  }
  
  void reset()
  {
    for (int x=0; x<w; x++)
    {
      for (int y=0; y<h; y++)
      {
        int idx = (y * w) + x;
        area[idx] = new Field(lib.frequencies);
      }
    }
    initial = true;
  }
  
  void displayTo(PGraphics pg)
  {
    pg.beginDraw();
    for (int y=0; y<h; y++)
    {
      for (int x=0; x<w; x++)
      {
        int k = (y * w) + x;
        Field select = area[k];
        if (select.needsUpdate)
        {
          int len = select.count;
          if ((1 == len) && (-1 != select.tileType))
          {
            int col = lib.elements[select.tileType];
            pg.set(x,y,col);
          }
          else if (lib.totalTiles == len)
          {
            int col = color(103,62,15);
            pg.set(x,y,col);
          }
          else
          {
            float ratio = float(len) / lib.totalTiles;
            int orange = color(214,75,0);
            int yellow = color(255,255,0);
            int col = lerpColor(orange, yellow, 1.0-ratio);
            pg.set(x,y,col);
          }
          select.needsUpdate = false;
        }
      }
    }
    pg.endDraw();
  }
}

void mouseClicked()
{
  wave.reset();
}

void setup()
{
  size(768,512,P2D);
  
  int outWidth = 60;
  int outHeight = 60;
  
  // Blocky rendering
  ((PGraphicsOpenGL)g).textureSampling(POINT);
  
  output = createGraphics(outWidth, outHeight, P2D);
  ((PGraphicsOpenGL)output).textureSampling(POINT);
  output.beginDraw();
  output.noStroke();
  output.background(0);
  output.endDraw();
  
  // Load image
  //PImage orig = loadImage("Hogs.png");
  //orig.loadPixels();
  
  // To me, looks like a contour map
  orig = createGraphics(32,32,P2D);
  orig.beginDraw();
  orig.background(0);
  orig.stroke(255);
  orig.noFill();
  orig.strokeWeight(2);
  orig.ellipse(orig.width/2,orig.height/2, 14,14);
  orig.loadPixels(); // this must be between begin- and endDraw()
  orig.endDraw();
  
  archive = new PatternLibrary(orig);
  wave = new Wave(outWidth, outHeight, archive);
  
  frameRate(30);
}

void draw()
{
  background(0);
  
  int SOME_NUMBER = 10; // repetitions per frame
  for (int i=0; i<SOME_NUMBER; i++)
  {
    wave.step();
  }
  wave.displayTo(output);
  
  image(output, 0, 0, 512, 512);
  image(orig, 600, 0);
  
  fill(255);
  text("frame: " + frameCount, 512, 20);
}

2 Likes

Thanks so much for sharing this! A wonderful contribution.

Great – you should add a comment about it to the repo README.md

@jeremydouglass Right you are. Added.

Beautiful ! Thank you so much for your contribution @noahbuddy.