4D Quaternion Julia Set Ray Tracer

Published on 20 September 2009

A quaternion Julia set is a four-dimensional equivalent of the standard two-dimensional fractal. By taking a 3D 'slice' through the 4D space it is possible to visualise a solid fractal.

This quaternion Julia set ray tracer is a port of Keenan Crane's Cg GPU implementation for Adobe Pixel Bender and a more general GLSL shader patch for Quartz Composer. The Pixel Bender version lets you render out high resolution images in Photoshop CS4 (or animate the parameters in After Effects). The Quartz Composer patch lets you animate any parameter in real-time.

The maths bit

Complex numbers have two components that define a point in a plane and are the core build blocks of fractal calculations. In 1843 the Irish mathematician Sir William Rowan Hamilton developed as a way of describing complex points in three-dimensional space. The only snag was he had to add a fourth dimension to make the maths work (hence the quad part of the name). Whilst complex numbers are described as the sum of a real and imaginary component: z = a + bi, quaternions are similar but have three imaginary components: z = a + bi + cj + dk

Fractals can be calculated using quaternions with the usual equation z n+1 = z n 2 + c where we track which points have a magnitude greater than the bailout threshold. The problem is trying to visualise the set of 4D points in 3D space. At this point I'll direct you to an excellent write-up Keenan has posted that describes the process in detail.

The ray tracing bit

The ray tracing is actually more like the ray-marching process I used in the procedural terrain renderer; a ray vector is fired into the scene and stepped forward until it gets within a collision threshold of the fractal. This brute-force approach is very inefficient but fortunately there are couple of optimisations that can be taken advantage of.

Firstly, a bounding sphere is defined as the limit of the fractal and used as the starting point for the initial ray vector for each pixel. This removes any unnecessary stepping at the start.

Secondly, there is a distance estimator function that will return the distance to the closest point on the Julia set for any z point in quaternion space. The distance estimator accelerates the ray tracing by a method called unbounding volumes, described in the paper Ray Tracing Deterministic 3-D Fractals [John Hart et al, 1989]. At each point of the rays journey through quaternion space the distance estimator returns the distance to the closest point in the Julia set. It means the next step the ray takes can be of this amount, which greatly reduces the overall number of steps required to intersect the fractal surface.

When the ray comes within a defined collision distance, epsilon, of the fractal surface the ray marching is stopped. Note, a fractal has an infinite level of detail so the epsilon factor acts to smooth out the surface and make it renderable.

The normal vector for each point is generated from the gradient of the approximate fractal surface. The full technique is described in the paper mentioned above.

Shading and ambient occlusion

My contribution with this shader is the addition of a basic ambient occlusion (AO) parameter. AO darkens areas of a surface that are in close proximity and less likely to receive light from its surroundings. When combined with standard shading, ambient occlusion gives a far more realistic looking image.

AO has been applied to quaternion Julia set fractals before by Inigo Quilez and animated here. However, these implementations are quite involved requiring multiple passes and extra sampling. The implementation used in this ray tracer is just an approximation and not physically accurate, but does add a visual improvement at very little additional computational cost.

Simply darken the surface by a factor proportional to the number of ray tracing steps taken for each point and the result is comparable to ambient occlusion. It works thanks to a side effect of the unbounded volume ray tracing approach; points that are least occluded are likely to be reached in fewer steps as the distance estimator function will return a large step size for the ray. Whereas rays to points nested in creases or dips will require more steps because the neighbouring surfaces will cause the distance estimator to return a smaller step size.

The image below illustrates the different components that makes up the final surface. Top row, left to right: pure ambient occlusion, pure phong. Bottom row: ambient occlusion and phong combined, and finally everything; a base colour with specular highlights and ray traced shadows.

Download and installation

Download the Quaternion Julia set ray tracer

Updated 14/12/2009 to fix an After Effects rendering bug.

For Pixel Bender open the QuaternionJulia.pbk file with the Adobe Pixel Bender Toolkit or copy it into the Pixel Bender Files folder in your Photoshop CS4 installation directory (you will need to have installed the PB plugin for Photoshop first).

The Quartz Composer patch was created with Quartz Composer 4, which is part of the new XCode. It will also work with QC3 as the patch is just a GLSL shader. I've had it running without antialising at up to 30fps @ 600x400 on my iMac with an NVIDIA GeForce 8800GS.

Just after porting this to Quartz Composer I discovered that Apple have released an OpenCL implementation for Snow Leopard. With native support now built into Snow Leopard, OpenCL opens up some very exciting possibilities for future projects.

Older: Artforms of Nature Newer: 3D Mandelbulb Ray Tracer