Newton Fractal Renderer

newton.py

import numpy as np
import colorsys
import cmath
import math
from multiprocessing import Process, Queue

# max: maximum coord to render to
# min: minimum coord to render to
# step: how much 1 pixel counts as
# chunkFactor: how many images to split it up into - allows for multithreading
# with the default settings, there are 16 chunks (so 16 resulting images)

iMax = 10 # 10
iMin = -10 # -10
iStep = 0.05 # 0.01
iChunkFactor = 4 # 4

jMax = 10 # 10
jMin = -10 # -10
jStep = 0.05 # 0.01
jChunkFactor = 4 # 4

darknessFactor = 10 # 10
maxIterations = 100 # 1000
precision = 0.001 # 0.001

def f(z):
	return z**3 - 1

def deltaf(z):
	return (f(z + precision) - f(z)) / precision

def newton(n):
	i = 0
	while abs(f(n)) > 0.001:
		i += 1
		n = n - f(n) / deltaf(n)
		if i > maxIterations:
			return (0, float("NaN"))
	return (i, n)

def colorize(angle, intensity, imDepth):
	angle = angle * (1/math.pi) #fix the angle, used to have a range of -pi to pi, now has a range -1 to 1
	if angle < 0:
		angle = 1 + angle #make angle loop, making the range 0 to 1

	if math.isnan(angle):
		r, g, b = colorsys.hsv_to_rgb(0, 0, 0)
		return (r, g, b)
	else:
		r, g, b = colorsys.hsv_to_rgb(angle, 1, 1)
		r = round(intensity * r)
		g = round(intensity * g)
		b = round(intensity * b)
		return (r, g, b)

def renderonechunk(iChunkMax, iChunkMin, jChunkMax, jChunkMin, chunkIdentifier, outputQueue):
	data = []
	for j in np.arange(jChunkMin, jChunkMax, jStep):
		#print(str(round(j/jStep)))
		for i in np.arange(iChunkMin, iChunkMax, iStep):
			data.append(newton(complex(j, i)))
	print("Chunk {0}x{1} rendered".format(chunkIdentifier[0], chunkIdentifier[1]))
	outputQueue.put({"chunkIdentifier": chunkIdentifier, "data": data})
	return {}

if __name__ == "__main__":
	processes = []
	outputQueue = Queue()
	for iChunk in range(iChunkFactor):
		for jChunk in range(jChunkFactor):
			#strange math to evenly partition each chunk
			iChunkMin = ((iChunk + 0) * (iMax - iMin) / iChunkFactor) + iMin
			iChunkMax = ((iChunk + 1) * (iMax - iMin) / iChunkFactor) + iMin
			jChunkMin = ((jChunk + 0) * (jMax - jMin) / jChunkFactor) + jMin
			jChunkMax = ((jChunk + 1) * (jMax - jMin) / jChunkFactor) + jMin
			chunkIdentifier = [jChunk, iChunk]
			renderArgs = (iChunkMax, iChunkMin, jChunkMax, jChunkMin, chunkIdentifier, outputQueue)
			processes.append(Process(target=renderonechunk, args=renderArgs))
			processes[-1].start()

	# block until all the processes finish
	# populate the list of chunks
	chunks = []
	for processIndex in range(len(processes)):
		chunks.append(outputQueue.get())
	print("Rendered {} chunks...".format(iChunkFactor * jChunkFactor))
	imDepthList = []
	# iterate through all the chunks
	for iChunk in range(iChunkFactor):
		for jChunk in range(jChunkFactor):
			# and find the depth of each of them
			imDepthList.append(sorted(chunks[(iChunk*iChunkFactor) + jChunk]["data"], key = lambda d: d[0])[-1][0])

	# colorize and combine chunks
	# find the highest image depth and use it
	imDepth = max(imDepthList)
	if imDepth == 0: #bound imDepth at 1, most image viewers don't open an image with 0 color depth
		imDepth = 1

	# iterate through all the chunks
	for iChunk in range(iChunkFactor):
		for jChunk in range(jChunkFactor):
			# perform transformations on the data
			chunks[(iChunk*iChunkFactor) + jChunk]["data"] = list(map(lambda d: (imDepth - d[0]*darknessFactor, d[1]), chunks[(iChunk*iChunkFactor) + jChunk]["data"])) #invert the iteration count
			chunks[(iChunk*iChunkFactor) + jChunk]["data"] = list(map(lambda d: (0 if d[0] < 0 else d[0], d[1]), chunks[(iChunk*iChunkFactor) + jChunk]["data"])) #turn all the negatives into 0

	print("Processed {} chunks...".format(iChunkFactor * jChunkFactor))
	iChunkMin = (0 * (iMax - iMin) / iChunkFactor) + iMin
	iChunkMax = (1 * (iMax - iMin) / iChunkFactor) + iMin
	jChunkMin = (0 * (jMax - jMin) / jChunkFactor) + jMin
	jChunkMax = (1 * (jMax - jMin) / jChunkFactor) + jMin
	yRes = round((iChunkMax - iChunkMin) / iStep)
	xRes = round((jChunkMax - jChunkMin) / jStep)
	imHeader = "P3\n{0} {1}\n{2}\n".format(xRes, yRes, imDepth)
	print("Created image header with depth {0}, X resolution {1}, and Y resolution {2}".format(imDepth, xRes, yRes))

	#finally, save all the chunks as images
	for iChunk in range(iChunkFactor):
		for jChunk in range(jChunkFactor):
			im = ""
			correctChunk = list(filter(lambda c: c["chunkIdentifier"][0] == jChunk and c["chunkIdentifier"][1] == iChunk, chunks))
			for d in correctChunk[0]["data"]:
				#print("ANGLE: {:.2f}".format(abs(cmath.phase(d[1])) * (1/math.pi)))
				#print("INTENSITY: {}".format(d[0]))
				r, g, b = colorize(cmath.phase(d[1]), d[0], imDepth)
				im += "{0} {1} {2} ".format(r, g, b)
			file = open("fractal{0}x{1}.ppm".format(jChunk, iChunk), "w")
			file.write(imHeader + im)
			file.close
	print("Saved {} chunks...".format(iChunkFactor * jChunkFactor))

zCompose.sh

#!/bin/bash
montage -tile 4x4 -geometry +0+0 -background none fractal* final.png