Skip to content

Instantly share code, notes, and snippets.

@x-or
Created July 11, 2013 08:07
Show Gist options
  • Save x-or/5973480 to your computer and use it in GitHub Desktop.
Save x-or/5973480 to your computer and use it in GitHub Desktop.
Multi-Layered XOR Pattern Motion Generator
# Multi-Layered XOR Pattern Motion Generator
# Copyright (c) 2013, Sergey Shishmintzev
# License: MPL/GPL
# Interpreter: Python 2.7.3 + numpy 1.6.2 + PIL 1.1.7
# P.S. Sorry for my English. It isn't my mother tongue.
import numpy as np
def bitreverse(n, l):
# Reversing a binary number
r = 0
for _ in range(l):
r = 2*r + (n&1)
n /= 2
return r
# Self-check
assert bin(bitreverse(1, 4))=='0b1000'
assert bin(bitreverse(2, 4))=='0b100'
assert bin(bitreverse(5, 4))=='0b1010'
def modulo_solve4(x, y, z, a, modulo):
# Solve: x + y + z + ? == a (mod modulo).
#return x^y^z^a # bitwise, modulo 2
return (a-x-y-z) % modulo
# Self-check
assert modulo_solve4(1, 2, 3, (1+2+3+4)%10, modulo=10)==4
assert modulo_solve4(1, 2, 3, (1+2+3+9)%10, modulo=10)==9
assert modulo_solve4(9, 2, 3, (9+2+3+1)%10, modulo=10)==1
assert modulo_solve4(4, 5, 6, (4+5+6+7)%10, modulo=10)==7
def undo_rule225(unk_state, next_state, i0, j0, modulo):
"""
Find the a fourth value, when three values are known,
satisfying a condition:
unk_state[i0, j0] + unk_state[i0, j0+1] + unk_state[i0+1, j0] + unk_state[i0+1, j0+1] == next_state[i0, j0] (mod modulo).
In the other words, find the state of the 2D additive cellular
automata using a known next state after the applying rule 225
(rule numbering from [Khan1997]). State to find is a partially known.
Input values:
unk_state - partially known state;
next_state - known next state, after applying rule 225;
i0, j0 - working position index;
modulo - modulo of the operators.
Changed values:
unk_state - state with zero or one solution added (solution will be added
at the some adjacent position near working position index).
Reference:
[Khan1997] A.R. Khan, et.al., "VLSI architecture of a cellular automata machine", 1997,
http://dx.doi.org/10.1016/S0898-1221(97)00021-7
"""
assert unk_state.ndim == 2
assert next_state.shape == unk_state.shape
# additional indices for the adjacent values
i1 = (i0 + 1) % unk_state.shape[0]
j1 = (j0 + 1) % unk_state.shape[1]
# prefetch a partially known state
s00 = unk_state[i0, j0]
s01 = unk_state[i0, j1]
s10 = unk_state[i1, j0]
s11 = unk_state[i1, j1]
# prefetch a next state
n00 = next_state[i0, j0]
# solve (s00 + s01 + s10 + s11) % modulo == n00
# Walk over adjacent values to find an unknown one
if s00 is not None and s01 is not None and s10 is not None and s11 is None:
unk_state[i1, j1] = modulo_solve4(s00, s01, s10, n00, modulo)
elif s00 is not None and s01 is not None and s10 is not None and s11 is not None:
pass # all adjacent values are known
elif s00 is not None and s01 is not None and s10 is None and s11 is not None:
unk_state[i1, j0] = modulo_solve4(s11, s00, s01, n00, modulo)
elif s00 is not None and s01 is None and s10 is not None and s11 is not None:
unk_state[i0, j1] = modulo_solve4(s10, s11, s00, n00, modulo)
elif s00 is None and s01 is not None and s10 is not None and s11 is not None:
unk_state[i0, j0] = modulo_solve4(s01, s10, s11, n00, modulo)
else:
assert False, "less than three adjacent values are known"
def bitwise_rule60(x):
"""
Perform rule 60 ( http://mathworld.wolfram.com/Rule60.html )
for every row using bitwise operators (only modulo 2.)
"""
#return x^np.roll(x, 1, axis=1) # considerably slow
assert x.ndim == 2
X = np.empty(x.shape, dtype=int)
I = x.shape[0]
J = x.shape[1]
for i in xrange(I):
for j in xrange(J):
X[i, j] = x[i, (j - 1) % J] ^ x[i, j]
return X
def normalize(x):
# Scale all array values into range [0, 255]
mn = np.min(x)
mx = np.max(x)
if mx - mn != 0:
X = 255*(x - mn) / (mx - mn)
else:
X = np.zeros(x.shape, dtype=x.dtype)
return X
_last_mean = 0
def mean_split(x):
"""
When modulo is large (like >3000), values of the state after
rule 60 are splitted into two noticeable close parts.
When modulo is small - split by mean as failback.
"""
global _last_mean
l = list(set(x.flat))
mean = max(np.mean(l), _last_mean)
_last_mean = mean # hack to preventing the oscilation
#l.sort(); print len(l), np.mean(l), (l)
top = x.copy()
bot = x.copy()
#print top
#print bot
top[x<mean] = 0
bot[x>=mean] = 0
return top, bot
def makecolor2(a, b):
# Convert two number arrays into the array of the bi-color RGB tuples.
assert a.ndim == 2
assert a.shape == b.shape
color = np.empty(a.shape, dtype=tuple)
for i in xrange(a.shape[0]):
for j in xrange(a.shape[1]):
color[i, j] = (a[i, j], (a[i, j]+b[i, j])/2, b[i, j])
return color
if __name__ == '__main__':
from PIL import Image
# Image size
power = 9
size = 2**power
half = 2**(power-1)
image = Image.new("RGB", (size, size), 0)
pix = image.load()
# known state after applying rule 225
br_pattern_state = np.zeros((size, size), int)
for i in xrange(size):
ri = bitreverse(i, power)
br_pattern_state[ri, i] = 1
for modulo in xrange(32, 1024*3+1): # image motion tiny cycle is 1024 (when power is 9)
print "modulo =", modulo
# partially known state
state = np.empty((size, size), object)
# fill initial (known) values
for i in xrange(size):
state[i, 0] = i % modulo
state[0, i] = i % modulo
# find rest values of the state
for i in xrange(size):
for j in xrange(size):
undo_rule225(state, br_pattern_state, i, j, modulo)
# the magic pass
state2 = bitwise_rule60(state)
# color magic
a, b = mean_split(state2 % modulo)
#a, b = mean_split(state2) # alternative view
c = makecolor2(normalize(a), normalize(b))
# produce image
for x in xrange(size):
for y in xrange(size):
pix[x, y] = c[y, x]
image.save("img-p%d-m%05d.png"%(power, modulo), "png")
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment