Cartesian Functions as Voxels

  • 28 November 2014
Voxel posts

  1. Intro
  2. Cartesian Functions (this post)
  3. Spherical Functions
  4. Curve functions

Of the three, the easiest class of functions to turn into voxels was the Cartesian functions. In order to understand the process the key thing I eventually realised was that voxelisation is almost the reverse of triangulation.

For the latter, we step through the x and y co-ordinates at a given resolution and calculate the consequent z-values to define the vertex co-ordinates of the triangle.

In contrast, for the former we step through not just the x and y co-ordinates, but also the z co-ordinate. The question we then have to answer is whether this (x, y, z) point is inside or outside the volume.

In practice, we can’t answer this question without defining the boundaries of the volume. Functy has an option to specify the thickness for a function, but the simplest case is to say anything on or below the function that defines the surface is inside the volume, and anything above is outside.


Inside and outside a Cartesian surface volume

The z-value generated by the function for given x and y co-ordinates represents the height of the function above the x-y plane. For a given voxel position (x, y, z) we therefore calculate the height of the function z‘ = f(x, y) at the point x, y and compare this to our actual position (x, y, z). If zz‘ then we’re inside the volume, so should fill in the voxel. If z > z‘ on the other hand, we’re outside the volume, so should leave the voxel empty.

And that’s it. To voxelise the function we check each of the points that make up the voxel space, perform the comparison, and either fill in or leave empty the voxels as we go. There are a lot of voxels to go through because we’re working with cubic dimensions (so even a low resolution of 100 steps per dimension gives us a million separate voxels to consider), but the comparison to perform is pretty straightforward in itself, as is clear from the code that performs the check.

fZSlice = psFuncData->fZMin + (nSlice * fZStep);

for (nX = 0; nX < nResolution; nX++) {
	fXFunc = psFuncData->fXMin + (nX * fXStep);
	for (nY = 0; nY < nResolution; nY++) {
		fYFunc = psFuncData->fYMin + (nY * fYStep);

		if (psCartesianData->psVariableX) {
			SetVariable (psCartesianData->psVariableX, fXFunc);
		if (psCartesianData->psVariableY) {
			SetVariable (psCartesianData->psVariableY, fYFunc);

		fZFunc = ApproximateOperation (psFuncData->psFunction);

		if ((fZFunc < fZSlice) && ((psFuncData->boMaterialFill == TRUE) || (fZFunc > (fZSlice - psFuncData->fMaterialThickness)))) {
			ucFill = 255u;
		else {
			ucFill = 0u;

Spherical functions are a little more complicated, but not much. I’ll write about them in the next post.

Sorry, the comment form is now closed.