Coding Curves 13: Superellipses and Superformulas

coding curves

Chapter 13 of Coding Curves

In this chapter we’ll be talking about some interesting shapes. Superellipses are quite useful in design and UI work, especially a specific superellipse we’ve come to know affectionately as the “squircle”. These are basically rounded rectangles, but with some neat properties. Superformulas are an extension of superellipses. They are more complicated, probably less useful, but interesting in their own right. Let’s dive in.

Superellipses

Sometimes you need to draw a rectangle for something. We’ve all been there. Drawing a rectangle is easy enough in almost any graphics package. So you draw a rectangle.

But that’s kind of boring. Those square corners are just so… square. So you come up with a rounded rectangle. Some drawing apis have this built in, but it’s not too hard to hand code. You draw a 90-degree arc at each corner and a line between each arc. This is ok.

In fact, you can make those corners have a relatively small radius, or you can make it quite round:

But you might not love this. There’s kind of a break in continuity in the edges. It’s going along nice and straight and then it abruptly turns into an arc. Maybe you want something that smoothly transitions into that round corner. Enter the superellipse.

You can think of a superellipse as a blend between an ellipse and a rectangle. Like an ellipse, and many of the other curves we’ve drawn in this series, we’ll consider this shape as a center point and a varying radius drawn from a t of 0 to 2 * PI. This is different from a lot of rectangle functions which take an x, y location (usually top, left of the shape) and a width and height. But you can adjust the function as you want.

While there are different ways to express a superellipse mathematically, we’ll dig around until we find one that is parameterized with that 0 to 2 PI t value. And here is is:

x = pow(abs(cos(t)), 2 / n) * a * sign(cos(t))
y = pow(abs(sin(t)), 2/ n) * b * sign(sin(t))

Hmm… more complicated than you might have expected?

There is a somewhat simpler formula that looks like this:

x = pow(cos(t), 2 / n) * a
y = pow(sin(t), 2 / n) * b

But the problem with this is that it only handles one quadrant: 0 to PI / 2. Rather than running a for loop four times, once for each corner, and figuring out the signs of x and y, we do some fancy math with abs and sign.

Your language probably has an abs function, but it just returns the absolute value of a value. You could easily write your own like:

function abs(val) {
  if (val < 0) {
    return -val
  }
  return val
}

This will always return zero or a positive value.

You might not have a sign function though. This just returns -1 if the given value is negative, and +1 if it’s positive. Returns 0 for 0. You can write this like so:

function sign(val) {
  if (val < 0) {
    return -1
  }
  if (val > 0) {
    return 1
  }
  return 0
}

You might be tempted to get fancy and do something like:

function sign(val) {
  return val / abs(val)
}

This works wonderfully. Until val is 0. Then it crashes. So you’ll still need a conditional in there to catch that case.

Now in the formula I gave above, there’s some duplication and some not very clear variables, and there’s no way to position the shape. So here’s an actual function that’s more usable:

function superellipse(xc, yc, rx, ry, n) {
  for (t = 0; t < PI * 2; t += 0.01) {
    c = cos(t)
    s = sin(t)
    x = pow(abs(c), 2 / n) * rx * sign(c)
    y = pow(abs(s), 2 / n) * ry * sign(s)
    lineTo(xc + x, yc + y)
  }
  closePath()
}

Most path drawing apis have some sort of closePath function that draws a final line back to the starting point.

And there you go. Let’s recreate our above rounded rectangle using a superellipse.

width = 600
height = 400
canvas(width, height)

// set color to orange however you do that...
superellipse(300, 200, 250, 150, 10)
fill()

You might or might not like that better than the simple rounded rectangle, but let’s explore it some more. That last parameter, n, controls how curved the corners are. The higher the value, the closer you get to a rectangle. Here is an n of 20:

And here is an n of 4:

This looks like the screen shape of an old fashioned TV set.

And this is a good time to bring up squircles.

Squircle

A superellipse where the x and y radii equal is sometimes also known as a supercircle. And has also been called a squircle. A combination of a square and circle. Some definitions of squircle state that it must have an n of 4. It looks like this:

This shape looks very satisfying, and has become very popular recently in user interface design. It’s often used for icons, especially for those used to launch apps on mobile devices. Sometimes the n value might be not quite 4 as above, but close enough.

There is actually an alternate formula for a squircle which John Cook wrote about at https://www.johndcook.com/blog/2022/10/27/variant-squircle/. His site is a great resource by the way. I can’t say I always understand all of what he is talking about, but he’s sent me down many very interesting paths of exploration.

Back to Superellipses

Let’s explore that n parameter a bit more. We already saw that as we make it larger, the corners get tighter. When we get an n of 2, something interesting happens. We just get an ellipse (or a circle if the height and width are equal).

When we go from 2 down to 1, we see the straight edges turn to corners, and the corners become more straight. Here’s 1.5:

When we get down to n = 1, corners and edges have completely reversed and we get a diamond shape.

Then as we go below 1, the corners start curving inwards. Here’s 0.75:

At 0.2, the superellipse is almost disappearing.

It will be essentially invisible with an n just a bit lower than that.

And when you go negative, things get really weird. Here’s n at -4. I had to reduce the size a bit so you could see what’s going on – we get a filled rect and some inverted superellipses extending out of each corner.

Not sure there’s any use for that last part, but there it is.

And that’s about it for superellipses. A great shape to have in your graphics toolbox. But let’s move on.

Superformulas

A superformula is a generalization or extension of s superellipse. Again, we have a radius that differs as we go from 0 to 2 * PI. The formula for that radius is:

That’s straight from the Wikipedia article on superformulas. It was easier to copy and paste that than it would be to try and type out the formula in a way that made sense. We’ll convert that to pseudocode shortly, but you should see some similarities between this and the superellipse formula. You see that we take the cosine and sine of something, divide it by a value, then take a power of the absolute value of that. The squiggly symbol there is the Greek letter phi, and it’s what we will call t. So we’re saying the radius at angle t is …

Then we have some other variables, m, n1, n2 and n3 as well as a and b which influence the x and y radii. In my implementation I’m going simplify it by just having a single radius. So I’ll set a and b to 1, which means we can just ignore them in the code and multiply the result by the radius we want.

Since this is so complicated, I’ll break down some of these parts into separate variables and then combine them for the final result.

First we’re taking the cosine and sine of the same value, so we can precalculate that. And I’ll call m symmetry.

angle = symmetry * t / 4

Then in between the big parentheses we have two terms we are adding together. We’ll calculate each of those:

term1 = pow(abs(cos(angle), n2)
term2 = pow(abs(sin(angle), n3)

Again, we’re considering a and b to be 1, so we can ignore them here.

Next we can plug these terms into the rest of the formula to get the final radius. Remember, we multiply that by the overall radius we want:

r = pow(term1 + term2, -1/n1) * radius

And finally, use the radius and angle to get the next point to draw to:

x = xc + cos(t) * r
y = yc + sin(t) * r

Put all together, this is our function:

function superformula(xc, yc, radius, symmetry, n1, n2, n3) {
  for (t = 0; t < 2 * PI; t += 0.01) {
    angle = symmetry * t / 4
    term1 = pow(abs(cos(angle), n2)
    term2 = pow(abs(sin(angle), n3)
    r = pow(term1 + term2, -1/n1) * radius
    x = xc + cos(t) * r
    y = yc + sin(t) * r
    lineTo(x, y)
  }
  closePath()
}

OK, now what the heck can we do with this? And what are all those parameters supposed to control?

When I’m faced with something like this, I usually start by finding some parameters that create some relatively simple, stable configuration and then start tweaking just one of the parameters to see what it changes. Then move on to another, and eventually look at how various parameters interact. Fortunately, the Wikipedia page on superformulas gives us this nice chart as a starting point:

The numbers on top of each image are the parameters m (what we call symmetry), n1, n2 and n3.

So we can immediately see that the symmetry parameter controls how many “nodes” will be in the shape.

So here’s a symmetry of 3, with all the n parameters set to 1:

And symmetry 5:

And 8:

Now, keeping symmetry at 8, let’s change n1 to 10:

Then down to 3:

And then below 1, down to 0.2:

OK, so we get that n1 makes the sides of the shape bulge in or out.

Alright, now we’ll stick with symmetry 8 and put n1 back to 1. Then we’ll up n2 to 1.5:

We still have 8 nodes, but every other one is a bit rounder, and the in-between nodes are a bit sharper. This gets more obvious at n2 = 2:

At n2 = 5, the rounder nodes have started to double, and the sharper ones are receding:

And here we are at n2 = 10. I had to reduce the radius because the shape was outgrowing the canvas:

Finally, setting n1 and n2 back to 1, we can try increasing n3 to 3:

This has a similar effect to that of n2 but changes alternate nodes, as if the whole shape was rotated 45 degrees.

At this point, I have a pretty good idea of what’s going on and I can just start trying some random parameters:

Symmetry 16, n1=0.5, n2=0.75, n3=2

Symmetry 32, n1=0.9, n2=0.2, n3=-0.3

Anyway, you don’t need me to choose random numbers for you at this point. Play with it. Use big numbers, small numbers, fractions, negatives, primes, numbers that are divisible by each other or not. And see what happens.

Remember, too, that we removed the a and b values from the original formula and used a single radius. You might want to try putting those back in and create some elliptical superformulas. The Wikipedia page also gives an alternate formula where the m parameter is split into two separate values, y and z:

This can allow for even more complex shapes.

I played around with this a little and came up with shapes like this, which I quite like:

But I’ll leave the coding of that up to an exercise for the reader.

Coding Curves 12: Guilloche Patterns

coding curves, misc

Chapter 12 of Coding Curves

Guilloche patterns are very intricate and fascinating patterns. You’ll often find them on bank notes and other official documents as well as watches and other intricate machinery. Because they are so intricate and complex, and often engraved into metal, they are usually done by machines themselves. Imagine a high end spirograph machine with a metal etching tool rather than a ball point pen. The name “Guilloche pattern” is rather vague and can apply to all kinds of similar patterns. I’m going to explain how to make a pattern like the one you see here:

three ringed guilloche pattern

This is a lot like the kind of thing you’d see on some kind of certificate or bank note, and once you understand what’s going on, you can adapt the code to make other similar patterns.

Step One – A Simple Ring

We’ll start by making a single, simple ring that looks like this:

a circle composed of many looping waves

This is very much like a trochoid or a rose curve. In fact, you could probably make this with those formulas, but I’m going to do it a bit differently to set it up for the more complex stuff we’ll do later on.

This is essentially a sine wave that’s wrapped around a circle. Note that it has an inner radius and an outer radius. The sine wave has 80 nodes, but it overlaps itself. To simplify what’s happening here, I’ll increase the inner radius and make it so it doesn’t overlap:

a circle drawn with a sine wave instead of a single curve

Now you can see more clearly that this is, as I said, a sine wave wrapped around a circle. Now I’ll add a little bit of overlap:

a circle drawn with an overlapping sine wave

Here, you should still be able to see the sine wave even with the overlap. The first image is the same idea, but with a lower inner radius and more overlap. Now let’s look how to draw this.

We start with an inner and outer radius. With a bit of easy math we can figure out a “mid” radius. This will be the the zero-point of the sine wave. And we’ll need a range, which is how much the sine wave will have to extend to either side of this radius to hit the inner and outer limits. You might also call this the amplitude of the wave.

width = 600
height = 600
canvas(width, height)

translate(width / 2, height / 2)

inner = 50
outer = 250
range = (outer - inner) * 0.5
mid = inner + range

Next, we need values for how many cycles in the wave, and how much it overlaps. I’ll call these nodes and div. These should be whole numbers and should not be evenly divisible. This is very similar to the n and d parameters in the rose curves we made in the last chapter, but I’ll name them differently to keep things clear.

nodes = 80
div = 11

Now we can loop t from 0 to 2 * PI * div and draw some line segments. The angle of each new point for the next segment is simply t and the radius will be computed as shown. We need to multiply 2 * PI by div to make sure we go around the circle enough times to meet back up correctly.

for (t = 0; t < 2 * PI * div; t += 0.01) {
  radius = mid + sin(t * nodes / div) * range
  x = cos(t) * radius
  y = sin(t) * radius
  lineTo(x, y)
}
stroke()

If you look back at the roses chapter, you’ll see this is very similar to how we got the radius there, but rather than figuring the radius as a sine and single multiplication, we’re using mid and range to fit the radius between inner and outer.

You can play with this a bit. Try different inner and outer radii and number of nodes and values for div. For a good Guilloche pattern, you probably want to keep nodes rather high, and div lower. But Most importantly, they should not be evenly divisible. An easy formula for decent patterns is to make div a small prime number and make sure nodes is not a multiple of div.

For example, if div is 17, you wouldn’t want to use 170 for nodes, or you’ll get this:

a rounded star kind of shape with 10 nodes

But changing nodes to 171 gives you a much better pattern:

a circle drawn with many overlapping loops

If you want to play with an interactive version of this, I made this a few years back:

https://bit101.github.io/lab/dailies/170120.html

Step 2 – A Complex Ring

Next, we’ll add some complexity to the mix. Rather than having our sine wave go back and forth between fixed inner and outer radii, we’ll have these radii themselves vary with a separate sine wave each! The result will look like this:

a single ring of a guilloche pattern

We’ll need a few more parameters for this. To calculate the final outer radius at any point, we’ll need the base outer radius, how many cycles of that sine wave and how far the final radius will differ from the base outer radius. And we’ll need the same three values for the inner radius. And of course we’ll still need nodes and div. The values I used for the above image are:

inner = 100.0
n0 = 7.0
h0 = 10.0

outer = 250.0
n1 = 17.0
h1 = 20.0

nodes = 142.0
div = 89.0

So the inner radius will have 7 nodes, and the radius will vary from 90 to 110, which means a base of 100, plus or minus 10. Similarly the outer radius will have 17 nodes and vary between 230 and 270.

All the code we did early to find the “mid radius” and range for the sine wave will have to be redone for each iteration of the for loop now, using these values and the dynamic inner and outer radii.

for (t = 0; t < 2 * PI * div; t += 0.01) {
  r0 = inner + sin(t * n0) * h0
  r1 = outer + sin(t * n1) * h1

  range = (r1 - r0) * 0.5
  mid = r0 + range

  radius = mid + sin(t * nodes / div) * range
  x = cos(t) * radius
  y = sin(t) * radius
  lineTo(x, y)
}
stroke()

You can surely simplify this code, but I wanted to write it all out explicitly for the sake of clarity.

Again this is a good stopping point to play with some of the parameters and see the different types of shapes you can make. We have one more level of complexity to tackle next.

Here’s an interactive version of this step:

https://bit101.github.io/lab/dailies/170121.html

Step 3 – Multiple Rings

The next step is to create multiple rings that will exactly fit together, like the image that appears at the start of this chapter:

a three ringed guilloche pattern

This is easier than it might seem. If you draw one smaller ring, and then one larger one around that, you just have to make sure that the parameters for the larger ring’s inner radius match the parameters for the smaller ring’s outer radius. Since we’ll be calling the same code multiple times, this would be a good time to turn it into a reusable function. For the most part, this just means wrapping the for loop and stroke call in a function with all the parameters passed in. I also added x and y parameters, so you can center the pattern anywhere on the canvas.

function guilloche(x, y, ir, n0, h0, or, n1, h1, nodes, div) {
  for (t = 0; t < 2 * PI * div; t += 0.01) {
    r0 = ir + sin(t * n0) * h0
    r1 = or + sin(t * n1) * h1

    range = (r1 - r0) * 0.5
    mid = r0 + range

    radius = mid + sin(t * nodes / div) * range
    lineTo(x +cos(t) * radius, y +sin(t) * radius)
  }
  stroke()
}

I upped the canvas size to 800×800 and then used this set of function calls:

guillloche(400, 400, 50, 6, 10, 120, 12, 20, 137, 37)
guillloche(400, 400, 120, 12, 20, 220, 18, 30, 141, 41)
guillloche(400, 400, 220, 18, 30, 350, 24, 20, 164, 53)

This resulted in the following image:

another three ringed guilloche pattern

Note that in the first call, the outer radius params are 120, 12 and 20. These correspond to the inner radius params in the next call. And the outer radius params for the second call, 220, 18 and 30, correspond to the inner radius params on the final call. This way, each ring lines up perfectly.

As usual, play around with this. You can make as many ring as you want. You might want to set up a custom data type that encapsulates all of the radius parameters to make that part more reusable. I’ll leave it to you.

Another thing you might want to try is making each ring a different color:

a colored three ringed guilloche pattern

Here’s an interactive version of the final product:

https://bit101.github.io/lab/dailies/170122.html

Note that there’s no rule that says the rings have to match exactly. Try putting some space between them or letting the overlap and see what kind of interesting patterns you can create that way. This next example has the first and last rings the same as in the previous example, but the middle ring is defined quite differently:

guillloche(400, 400, 20, 4, 5, 120, 8, 10, 137, 37)
guillloche(400, 400, 160, 5, 24, 160, 11, 24, 80, 17)
guillloche(400, 400, 220, 18, 16, 350, 10, 20, 164, 53)

Here, the inner and outer radius is the same, 160, which makes some interesting blobs rather than a discrete ring. All is fair in Guilloche patterns!

a three ringed guilloche pattern. the inner ring does not tightly interlock with the others

And of course, don’t forget to animate!

an animated guilloche pattern

Summary

Do an image search for Guilloche patterns and you’ll find things that are much different, and often more complex than what I’ve presented here. This is just one take on them.

You might try altering the code do draw something other than circles. Ellipses would be a good place to start. And then move onto other shapes. But all that is a bit beyond the scope of this chapter.

Coding Curves 11: Roses

coding curves

Chapter 11 of Coding Curves

Now we come to another one of my favorite types of curves – roses or rose curves. To me, these look a lot like circular Lissajous curves, or very regular harmonographs. In fact, they are a special instance of hypotrochoids, but special enough to look at on their own. Just to give you some instant visuals, here’s a rose curve:

Like many other curves we’ve looked at, we can get a parameterized formula that will take a t value that goes from 0 to 2 * PI and give us back a value that will let us plot the curve. Let’s look back at the formula for a circle first…

function circle(x, y, r) {
  for (t = 0; t < 2 * PI; t += 0.01) {
    x1 = x + cos(t) * r
    y1 = y + sin(t) * r
    lineTo(x1, y1)
  }
}

You could simplify that into a single line within the for loop, but I wanted to spread it out for clarity.

A rose curve uses the same strategy, but instead of a fixed radius, that radius is constantly changing, also based on the t value, as well as another parameter. Here’s the formula for the radius:

r = a * cos(n * t)

So we have two new variables here. a is the overall radius of the rose, and n controls the number of petals in the rose (although the petal count gets a bit complicated, so we’ll come back to that shortly). So we can make a rose function like this:

function rose(x, y, a, n) {
  for (t = 0; t < 2 * PI; t += 0.01) {
    r = a * cos(n * t)
    x1 = x + cos(t) * r
    y1 = y + sin(t) * r
    lineTo(x1, y1)
  }
}

And now, if you want to, you can clean this up a bit:

function rose(x, y, a, n) {
  for (t = 0; t < 2 * PI; t += 0.01) {
    r = a * cos(n * t)
    lineTo(x + cos(t) * r, y + sin(t) * r)
  }
}

For now, let’s just say that n should be a positive whole number. But we’ll explore ranges beyond that of course.

Now we can draw our first rose like so:

width = 800
height = 800
canvas(800, 800)

rose(width / 2, height / 2, width * 0.45, 5)
stroke()

I’ll be using the width * 0.45 a lot here. It just makes the radius a bit less than half the size of the canvas, so the curve will go almost to the edge of the canvas, but never hit it.

And this gives us a 5-petal rose:

The first example at the top of this page used an n of 7. And here is a rose with an n of 11:

So far we’re seeing a good correlation between n and the number of petals. At least for odd values of n. But what if we use an n of 4?

Interesting. This gives us eight petals. This holds true for any value of n. Odd values create n petals. Even values create 2 * n petals. Just to go way out in one direction, here’s one with n = 40, which gives 80 petals. I had to up the resolution – incrementing t in the for loop by 0.001 to keep it from getting jagged.

In the opposite direction, going down to n = 1, gives you a single node:

A bit strange, but it works out mathematically. You’ll find that for negative values, the rose looks the same as for positive values of n. Here’s 5 on the left and -5 on the right:

Unsurprisingly, n = 0 gives us nothing. And so that covers all the whole number roses. If that’s all there was to roses, it would be nice, but there’s a lot more to go.

An Alternate Rose

Actually, before I move beyond whole numbers of n, I want to just mention an alternate rose formula. Instead of using cosine in the radius formula, you can use sine instead:

r = a * sin(n * t)

This gives you the same roses as the original, but rotated. Here’s a 5-petal rose using the original cosine on the left and sine on the right:

And the same for a 8-petal rose (n = 4):

The actual amount of rotation is PI / (2 * n) radians, or 90 / n degrees. For odd values of n, this always has the visual effect of rotating the rose by 90 degrees (the actual rotation may be different, but due to rotational symmetry, it appears to rotate 90 degrees). For even values of n, it rotates the rose so the petals will now be where the spaces between the petals were in the original version.

Fractional values of n

Things start to get more interesting when we start using fractional values for n. We can try it generating a rose with:

rose(width/2, height/2, width * 0.45, 5.0 / 4.0)
stroke()

But this gives us, rather disappointingly, the following:

The problem is that it’s going to have to go beyond 2 * PI to finish it’s cycle. How far beyond? Well, to figure that out programatically, we’ll need to first ensure that the n value is rational. If it’s an irrational number, the rose will continue forever without reaching its exact starting point. We’ll also need to know both the numerator and denominator of that fraction. We can adjust the rose function to take an extra value, so we have n and d for numerator and denominator.

rose(x, y, a, n, d) {
  for (t = 0; t < 2 * PI; t += 0.01) {
    r = a * cos(n / d * t)
    lineTo(x + cos(t) * r, y + sin(t) * r)
  }
}

This doesn’t solve the problem yet, but gets us the first step. If you want you can enforce n and d to be integers to make sure you’re getting a rational fraction, but make sure you convert them so the division in line 3 returns a floating point value.

Now we need to change the for loop limit from 2 * PI to the actual value we need. That limit value is:

limit = PI * d * m

But what is this new m value there? Well, m should be equal to 1 if d * n is odd. And m should be 2 if d * n is even. Woo! A bit complex. But we can simplify it.

We usually test for evenness by taking a number modulo 2. If the result is 0, that means the number is even. If the result is 1, the original number is odd. So we want:

m = 1 when d * n % 2 == 1

and

m = 2 when d * n % 2 == 0

So we can say:

m = 2 - d * n % 2

This gives us:

rose(x, y, a, n, d) {
  m = 2 - d * n % 2
  limit = PI * d * m
  for (t = 0; t < limit; t += 0.01) {
    r = a * cos(n / d * t)
    lineTo(x + cos(t) * r, y + sin(t) * r)
  }
}

Remember, if you are enforcing integers for n and d, you might need to do some casting or conversion to make everything work correctly. I’ll leave that to you. Now we can redo the fractional one like so:

rose(width/2, height/2, width * 0.45, 5, 4)
stroke()

And now we get something much nicer:

This time, the rose continued all the way around and completed itself.

Now you can go to town trying different fractions. I find that things get really interesting when you use higher numbers that are very close to each other. For example, n = 22, d = 21:

Or even 81 and 80:

Roses with fractions less than 1

Things become a whole different type of interesting when you get fractions that are less 1.0. For example, here are roses with n and d of 1,2 on the left, 1,3 in the middle, and 1,4 on the right.

A trick to find interesting patterns is to take a pair of numbers that would usually reduce down, like 17 / 51 will reduce to 1 / 3, giving us the middle figure above. But then shift one of the values a bit. Here’s 17 and 52:

A big difference for just a shift of 1.

Named Roses

Some of these rose curves have special names. I’ll share some of them.

Limaçon Trisectrix

This has a ratio of 1 / 3. We already saw this one above.

Dürer Folium

With a ratio of 1 / 2. Also seen previously.

Quadrifolium

Ratio is 2 / 1

Trifolium

Ratio of 3 / 1

Maurer Roses

If you thought we were almost done, wrong! There’s a whole other type of rose curve to explore – Maurer roses!

Maurer roses start with the basic rose function, but instead of just drawing the curve all the way around, it draws a series of line segments to points along the rose curve. Although it doesn’t have to be so, this is often done with 360 segments and the angles used are specified in degrees. We construct a rose, here using a ratio of 4 / 1, and then pick a degree value to step by. In this case, I chose 49. Then we loop t from 0 to 360 and multiply t by that degree value. So the degrees goes from 0, to 49, 98, 147, 196 and so on. We use that value in our rose (converting to radians of course) and use that at the next point. Here’s what it looks like in action for the first 30 iterations:

To put it a different way, in a normal rose curve, we are incrementing in very tiny increments, so we get a very smooth curve. Here, we are incrementing in gigantic jumps, so we get what looks like is going to be a chaotic mess. But, if we let it finish its full path through to 360 iterations, we get…

Aha! Not a chaotic mess after all! In fact, quite nice. Actually, above I’ve drawn the regular rose on top of the Maurer rose. Here is the Maurer all by itself:

I think the two combined look really nice.

So how do we do this?

Well, again, we start out with the basic rose function. But in this case, we’ll just stick to a single integer value. So just n rather than n and d. But we also want to specify how many degrees to jump on each iteration. To avoid confusion with the earlier d parameter, I’ll call this deg. So the signature is:

function maurer(x, y, a, n, deg) 

Again, we want to loop from 0 to 360 for our initial t value. And then we want to get that value that is t multiplied by deg. This is the degree value shown in the animation above. We’ll call it k but at this point we’re done with degrees, so we’ll convert it to radians by multiplying by PI and dividing by 180

function maurer(x, y, a, n, deg) {
  for (t = 0; t < 360; t++) {
    k = t * deg * PI / 180
    r = a * cos(n * k)
    lineTo(x + cos(k) * r, y + sin(k) * r)
  }
}

We’ll then just execute the rose algorithm, but using k instead of t.

Now we can set something up like the following.

width = 800
height = 800
canvas(800, 800)

maurer(width / 2, height / 2, width * 0.45, 5, 37)
stroke()

// drawing the regular rose is optional
rose(width / 2, height / 2, width * 0.45, 5, 1)
stroke()

And get this:

Play around with different values for n and deg. You’ll find that n works the same way it did for regular roses. But minor variations in deg can create radically different images. For example here is n = 7 and deg = 23:

But moving deg up to 24 gives you this:

Not nearly as nice.

Generally, you’ll find that even numbers for deg will have a lower chance of being interesting than odd numbers (with exceptions).

And anything that divides evenly into 360 is not going to be great. For example, here’s 4, 120:

I drew the full rose too, but the Maurer is just the triangle on the right hand side. Increase that to 121 though, and you get this beauty:

Also, lower prime numbers usually always work pretty well. I’ve noticed that the lower values of n let you get away with higher prime numbers for deg. But it’s something I haven’t tested very thoroughly. Something to play around with.

One more thing you might want to try is fractional Maurer roses. You don’t even have to alter the code at this point. You can just enter the fraction. Because we are always looping from 0 to 360, we don’t need to adjust for a different number of loops. Here’s one to start with. Make sure you put both fraction values into the rose function separately, if you are using that.

maurer(0, 0, width * 0.45, 5.0 / 4.0, 229)
stroke()
rose(0, 0, width * 0.45, 5, 4)
stroke

See what you can find among all the possible variations.

2022 in Review

misc

I’ll keep this light this year. Professionally, a strange and difficult year. But at least I still have a job, so I’m grateful for that.

Beyond the job though, I did do some fun coding projects.

Probably the biggest was my raytracing journey. I’ve never been super interested in that kind of photo-realistic rendering from an artistic viewpoint. But working through the books and creating the renderer itself was the most fun I’ve had coding in ages. Of course, once it was sufficiently “complete”, I kind of lost interest in it. It didn’t really change the basic fact that 3D stuff is not what I’m really interested in creatively. That said, I’ll probably play with it some more in the future. My original thinking was to do some kind of generative stuff that could then be rendered in 3D. But “generative” fairly often means lots of individually created units interacting in some way – at least the way I usually do generative. so it can get expensive rendering hundreds or thousands of objects. I might need to look more into optimizing the engine – I skipped that chapter in the book!

Other than that, I created a few vim plugins. One of which, bufkill, I use on a regular basis and really love. I made a few others, but none that were really worth talking about too much. During my plugin-writing phase, I also did a lot of upgrading of my neovim setup, switching my whole config to use Lua, and dropping Conquer of Completion in favor of neovim’s built-in LSP system for completion and syntax checking. That was a lot of work, but I’m happy with the result.

And of course, I started the Coding Curves series, which I’m still really excited about. As I said in the initial post, this was a book that I’ve wanted to write for years. I finally came to the conclusion that I was not going to sit down and write it from start to finish, but writing each chapter as a blog post was way more doable. I banged out ten “chapters” so far and they’ve been pretty well received. I’m taking a short break, but will be back to do the rest before too long. Most of what’s left I already have lots of code and examples written for, and a lot of the explanations are kind of half written in my head already.

I haven’t put out a lot of creative-code artwork this year, but in the last few weeks I’ve been working hard to upgrade my personal code libraries. I currently use a couple of libraries of my own in every creative coding project I do: blcairo which is the drawing api with a ton of custom rendering routines, and bitlib which has math, geometry, color, a PRNG, noise, and other various utilities I use all the time. I love this library and I’ve added a lot to it in recent months.

Oh, and I quit Twitter. I fully deleted my account and set up a new private account with the old user name just to keep anyone from taking it over and impersonating me. It was tough to walk away from 8,000 followers, but every time I see the toxic mess that the platform has become, I have zero regrets. Find me now at https://mstdn.social/@bit101. I remain very happy with Mastodon.

Plans for 2023

There are a few things I want to do in the new year.

  • Create my own interpreted language. A while back I got the book Writing an Interpreter in Go and worked partway through it. It’s quite good and I want to do more with this. My idea would be a simple interpreted language that would depend on my own Go libraries as described above. This would probably just be for my own education and personal use. I don’t have a lot of interest in supporting a language in the wild.
  • I want to finally do something with music. At this stage in life, I doubt I will ever learn any kind of instrument, but I have enjoyed doing stuff with electronic music on and off over the years. I’m currently messing around with MilkyTracker, which reminds me of playing with mod trackers back in my Amiga years. I know trackers are kind of toys, but that’s all I’m really looking for now. Something to play with and have fun and maybe learn a bit. Eventually I might like to dip my toes into musical hardware, synths, sequencers, etc. Anyway, one of my goals will be to release a song in 2023. “Release a song” sounds way more serious than I mean it to sound. I will create a song and put it on my web site and hopefully it won’t suck too much (but, fair warning, it probably will). If anyone has any advice or resources for learning this stuff, let me know.
  • I’d love to do another side project creating graphics for something. My last big projects were back in 2018. I’m open for projects!
  • Of course, I’ll finish the Coding Curves project. Probably in January or February. I’m toying with the idea of a similarly structured project with a different subject, tba.

Coding Curves 10: Spirals

coding curves

Chapter 10 of Coding Curves

OK, let’s talk spirals.

Spirals are a lot like circles, in that they are a set of points with a distance relationship to a fixed center point. But unlike circles, where that distance is fixed, with spirals, that distance varies. The distance from a given point to the center point is generally a function based on the angle between those two points. So you’ll usually have some function that takes in an angle and returns a radius. Then you can use the radius and angle to find the x, y position of the point at that angle. There are many different spiral formulas, which give you spirals that have a different look and feel. Let’s start with one of the most basic spirals.

Side note: a circle is sometimes called a “degenerate” spiral. No disrespect intended. It just means that a circle follows the “rules” of spirals but isn’t really what we would think of as a spiral. Like a triangle where one side has a length of zero. All the triangle math generally works fine, but it’s really just a line to our eyes.

The Archimedean Spiral

As you can see, on each cycle, the radius grows by a fixed amount. Here’s the formula for this spiral:

r = a * t

Here, t is the angle and a is some constant, in the case of the above image, it’s 5. The product of these two is the radius at the angle t. If we increase t to 10, we get this:

Now you can see that the a constant determines the spacing between each cycle of the spiral. In this case, it’s grown outside the bounds of the canvas.

Now when we are drawing spirals, we have to decide how many cycles we’ll be drawing. How many times are we going to go around? If t goes from 0 to 2 * PI, we have one cycle:

That’s one cycle. As you can see, for this spiral, the curve starts in the center and expands out from there. Three cycles means t goes from 0 to 2 * PI * 3.

Given all that, we can start to put together a spiral playground, like so:

width = 800
height = 800
canvas(width, height)
translate(width / 2, height / 2)

cycles = 10
res = 0.01

for (t = res; t < 2 * PI * cycles; t += res) {
  r = archimedean(5, t)
  x = cos(t) * r
  y = sin(t) * r
  lineTo(x, y)
}
stroke()

function archimedean(a, t) {
  r = a * t
  return r
}

We say we want 10 cycles, so the for loop goes to 2 * PI * cycles. We call the archimedean spiral function, passing in 5 as the constant, and t as the angle. That gives us a radius, which we use with the angle to find an x, y point and draw a line to it.

Throughout the chapter, I’ll be showing you other functions for different spirals. You can just replace the call to archimedean with the other functions.

The above code gives us this image:

In the drawing api I’m using, angles that increase positively go around in a clockwise direction. So this spiral is going around clockwise from the center. This may be different for you. It depends on how your drawing api handles angles. To reverse the direction, you have a few options:

A. Use scale to flip the canvas

scale(1, -1)

B. Change the code that creates the points to make one axis negative. For example:

  x = cos(t) * r
  y = -sin(t) * r

C. Change the for loop to go backwards:

for (t = -res; t > -2 * PI * cycles; t -= res) {

Doing any of those should get you going around the other way.

One last thing…

One thing to note in the code, is that I defined the for loop to start with res, not 0.

for (t = res; t < 2 * PI * cycles; t += res) {

To see why I did that, let’s move to our next spiral.

The Hyperbolic Spiral

This spiral looks quite different from the first one. Each cycle is not a fixed distance from the past one. You’d be tempted to say that the distance increases on each cycle, but hold that thought.

Here’s the function for the this spiral:

function hyperbolic(a, t) {
  r = a / t
  return r
}

Not very different at all. We’re just dividing a by t rather than multiplying. For the above image, I passed in an a value of 1000. If I bring a down to 10, we get a tiny little spiral like this:

It should also now be obvious why I started the for loop with res rather than 0. If t is 0, then we’d be in a divide-by-zero situation here, which would cause some problems. This might cause a crash, or might just give you a radius of NaN (not a number), which would not be very useful. So we start with a number we know will not be zero.

Another thing to note about this spiral is the direction it goes. I drew the above spiral with 20 cycles. If I move that down to 5 (moving a back up to 1000), we get:

Maybe you can now see that the spiral starts large and grows smaller as t increases. With 100 cycles, the spiral starts to jam up in the center.

So, as I said earlier, at first glance the radius seems to be growing on each cycle, but now you can see that it’s actually starting out large and getting smaller on each cycle.

Let’s move on.

The Fermat Spiral

This is a nice one. Here’s the formula:

function fermat(a, t)  {
	r = a * pow(t, 0.5)
	return r
}

Here we are multiplying a by t to the power of 0.5. This assumes you have a function called pow that raises it’s first argument to the power of the second argument. Depending on your language, you might also be able to say:

  r = a * (t^0.5)

or

  r = a * (t**0.5)

The above image was drawn with 20 cycles and an a of 20. Drawn with 10 cycles, we get:

So you can see this spiral draws from the center outwards. It starts out with a relatively large space between each cycle, but that space reduces slightly the further you go out.

Going back to 20 cycles and changing a to 40:

We see that the distance between cycles has increased.

Here we are with 100 cycles and an a of 15.

At the end, there is hardly any space at all between each cycle. And we have some interesting moiré patterns going on there.

The Lituus Spiral

This looks a lot like the Hyperbolic spiral, but the formula is a lot closer to the Fermat spiral:

function fermat(a, t)  {
	r = a * pow(t, -0.5)
	return r
}

We just used a power of -0.5 rather than 0.5. The above spiral was drawn with 20 cycles and an a of 500. Here it is with 10 cycles:

So you can see that this is another one that draws in towards the center. Going back to 20 cycles, if we lower a to 50, we get:

The word “lituus” originally meant a curved staff, wand or horn. The above image explains best how this spiral got its name.

If we raise a to 1000, we get:

So you could say a lower a causes the spiral to get sucked in to the center more quickly.

So many spirals! Let’s do some more!

The Logarithmic Spiral

This has a feel of being the opposite of the Fermat spiral. The spacing starts quite small, and increases the further it goes out. The formula:

function logarithmic(a, k, t)  {
	r = a * exp(k * t)
	return r
}

This one is a bit more complex than the others we’ve seen. First off, it has two parameters besides the t angle. And we’ll have to get into the exp function.

We actually touched on exp in the Harmonographs chapter. To recap, we have a mathematical constant, e known as Euler’s Number. It’s value is roughly 2.71828. When working with logarithms, it’s quite common to raise e to a power. So some math libraries have included a function to do that directly, often called exp. For a concrete example, let’s take a look at Javascript. It has e as a constant, Math.E. So to get e to the power of 2, you could say:

Math.pow(math.E, 2)

But it also has an exp function, so you can do the exact same thing by saying:

Math.exp(2)

So with this line,

r = a * exp(k * t)

we are multiplying a by e raised to the power of k * t. In the above image, I had a set to 0.5 and k at 0.05.

If we raise a to 1, we get a much larger spiral, though it looks pretty much the same.

Bringing a down to 0.25 gives us:

Again, similar look and feel, but smaller overall.

I’ll reset a back to 0.5 and move k up to 0.1.

Now you can see it expands much faster. Bringing k down to merely 0.04 has a much bigger effect than I expected:

The Golden Spiral

This spiral increases by a rate equal to the “golden ratio”, approximately 1.618. Here’s the formula:

function golden(t) {
	r = pow(PHI, 2 * t / PI)
	return r
}

First, note that this function has no parameters other then t. It’s hard-coded with the value PHI, which is the golden ratio. Many math libraries have that value built in as a constant. If yours doesn’t, you can approximate it by saying

PHI = 1.61803

Or, if you want to be exact, you can use the following to get the exact value:

PHI = (1 + sqrt(5) ) / 2

Then save that as a constant somewhere and use it as needed.

You’ve probably seen this spiral in images like this:

By Romain – Own work, CC BY-SA 4.0, https://commons.wikimedia.org/w/index.php?curid=114415511

This is actually the Fibonacci spiral. The size of each square the sum of the next two next smaller squares and each curve is made up by a 90-degree arc centered in one corner of each square. It isn’t precisely the same as the Golden spiral, but very close.

… and more

Peruse this list of spirals if you want to find some other interesting ones to try out:

https://en.wikipedia.org/wiki/List_of_spirals

Here’s another good one:

https://mathworld.wolfram.com/topics/Spirals.html

Spirangles

I ran across this term while writing this chapter. Essentially it’s a spiral made by straight line segments. By changing the angle between each segment, you can form different shapes.

These are very easy to make with our existing setup. The above one was made with the archimedean function, an a of 3, and 20 cycles. The trick is reducing the resolution. For this one, I set:

res = PI * 2 / 3

Now on each step of the for loop, t will increase by one-third of a circle. Here are some others, dividing by 4 and 5.

You get the idea. You might want to experiment with some of the other functions too. Some of them are quite satisfying.

Sunflowers

No discussion of spirals would be complete without talking about sunflowers. Just do a search on some combination of the terms, sunflowers, fibonacci, and spiral and you’ll get a ton of reading material and pretty pictures. What it all means, I’ll leave to others, but it’s fun to draw the kind of spiral you find looking at sunflowers. The golden ratio is built right into this pattern, and here’s what it looks like:

You can see not one, but many, many spirals going in and out at different angles. Here’s the code I used to make this (i.e. the pseudocode that represents the actual code I used):

width = 800
height = 800
canvas(width, height)
translate(width / 2, height / 2)

count = 1000
for (i = 0; i < count; i++) {
  percent = i / count
  size = 14.0 * percent
  r = 380.0 * percent
  t = i * PI * 2 * PHI
  x = cos(t) * r
  y = sin(t) * r
  circle(x, y, size)
  fill()
}

We choose a count for how many “sunflower seeds” we want to draw. Here it’s 1000.

Then we loop through using i and get a percent value that it i / count.

The size value is the radius of each “seed”. As i approaches count, percent will approach 1, so size will approach a maximum of 14.

Similarly r is the radius at which to place each seed. It will max out at 380, just shy of half the width of the canvas.

We calculate t using line 11 there, which is just a magical sunflower Fibonacci formula. Seriously though, read up on it if you want to know more. With an angle and radius, we can get an x and y and draw the seed there with the current size.

And that’s about all I have to say about spirals for now. Catch you next time!

Coding Curves 09: Roulette Curves

coding curves

Chapter 9 of Coding Curves

Initially I was going to title this chapter “Trochoids and Cycloids”. I thought they were two different, but related things. As I got into it, I realized I was very confused about what each thing was. Actually a cycloid is just a very specific type of trochoid. I’ll forgive myself though. Here are the definitions of each on Wikipedia:

In geometry, a trochoid (from Greek trochos ‘wheel’) is a roulette curve formed by a circle rolling along a line.

In geometry, a cycloid is the curve traced by a point on a circle as it rolls along a straight line without slipping.

Well, not only are they not two different things, they sound almost identical from those descriptions. The devil is in the details. So let’s explore.

There are three different types of trochoids:

  • Common trochoids (also called cycloids!)
  • Prolate trochoids
  • Curtate trochoids

Beyond those, there are a number of related curves, most of which we’ll cover here:

  • Epitrochoids
  • Hypotrochoids
  • Epicycloids
  • Hypocycloids
  • Involutes

All together, these make up the family of roulette curves. We’ll cover all but Involutes in this chapter.

Now that I’ve thrown a whole bunch of meaningless words at you, let’s figure out what all these things are. Starting with trochoids.

Trochoids

So, as described above, a trochoid is the curve formed by rolling a circle on a line. Let’s visualize that before we get into coding:

If we trace the path of that black dot on the edge of the circle…

…that new curve is a trochoid. In fact, because that drawing point lies exactly on the edge of the circle, it’s a “common” trochoid, also called a cycloid.

If we extend that point out beyond the edge of the circle, we get what’s called a prolate trochoid.

And if the point is inside, the circle, it’s a curtate trochoid.

To make these animations, I moved the circle along from left to right, figured out what its rotation should be at each location, and used sine and cosine to figure out where that point would be based on the position of the circle and its rotation. Then drew lines to that point. But there’s a somewhat more direct formula for trochoids:

x = a * t - b * sin(t)
y = a - b * cos(t)

Here, t is the angle of rotation of the circle, and can just increase infinitely, a is the radius of the circle, and b is the distance from the center of the circle to the drawing point – you could say the radius of that point. Let’s try it out.

width = 800
height = 300
canvas(width, height)

translate(0, height/2)
scale(1, -1)
moveTo(0, 0)
lineTo(width, 0)
stroke()

a = 20.0
b = 20.0
res = 0.05

for (t = 0.0; t < width; t += res) {
	x = a * t - b * sin(t)
	y = a - b * cos(t)
	lineTo(x, y)
}
stroke()

Since this formula is written for Cartesian coordinates, we’ll translate the y-axis to the center of the canvas and flip it.

Then we’ll draw a horizontal line through the center to represent the “floor” that the circle is rolling on. That part is optional.

We set a and b both to 20, so well get a cycloid here. We loop through from 0 to the width of the canvas as t, apply our formulas to that value, and draw a line to the resulting point.

If we raise b to 60, we get prolation.

And if we lower it to 10, we see significant curtation.

I don’t know if prolation and curtation are actual words, but they sound pretty cool.

Anyway, that’s about all there really is to trochoids. Try different values for a and b, or whatever else you want to mix it up with, but I don’t have much more to say about them. But there are several other roulette curves we have left to discuss here.

Centered Trochoids

The next four curves we’ll look at are called centered troichoids. The difference is that rather than the circle rolling along a line, it’s rolling along another circle, either inside or outside that other circle.

The four curves are:

  • Epicycloid – a curve formed by the path of a point exactly on the edge of a circle which is rolling around the outside another circle
  • Epitrochoid – same as above, but the point is within or outside of the moving circle, not on the edge, the same as a curtate or prolate trochoid.
  • Hypocycloid – a curve formed by the path of a point exactly on the edge of a circle which is rolling around the inside another circle
  • Hypotrochoid – same as above, but the point is within or outside of the moving circle, not on the edge

In addition to these there are a number of named curves that are just one of the above with a particular ratio between the radii of the two circles. We’ll take a look at a few of them.

I find these curves to be a lot more interesting than regular trochoids, with a lot more variety. So let’s dig in.

Epitrochoids

First, let’s visualize a circle rolling around the outside of another circle.

And let’s see what we get if we trace the path that black point is drawing.

There’s your epitrochoid!

And in fact, because that point is on the edge of the moving circles, it’s also an epicycloid!

If we move the point out a bit…

And we can move it in a bit.

Again, I made these animations by figuring out where the circle would be, drawing the circle, figuring out what it’s rotation would be at that stage of its journey and drawing the dot, and then connecting the paths of where that dot was on each frame into a curve. It was worth doing in order to show everything at work like that, but there’s a (relatively) simple formula you can just apply:

The formula for an epitrochoid:

x = (r0 + r1) * cos(t) - d * cos(((r0 + r1) * t) / r1)
y = (r0 + r1) * sin(t) - d * sin(((r0 + r1) * t) / r1)

The parameters:

  • r0 is the radius of the fixed circle
  • r1 is the radius of the moving circle
  • t is the increasing angle
  • d is the distance from the center of the moving circle to the drawing point (the same as r1 for an epicycloid).

We can use this in some code to draw an epitrochoid in one shot.

width = 800
height = 800
canvas(width, height)

translate(width/2, height/2)

r0 = 180
r1 = 60
d = 60
res = 0.01

circle(0, 0, r0)
stroke()

for (t = 0.0; t < PI * 2; t += res) {
  x = (r0 + r1) * cos(t) - d * cos(((r0 + r1) * t) / r1)
  y = (r0 + r1) * sin(t) - d * sin(((r0 + r1) * t) / r1)	lineTo(x, y)
}
stroke()

We set r0, r1 and d, and then draw the fixed circle just for reference.

Then we just loop t from 0 to 2 * PI, get an x,y point, and draw a line to that point.

And this gives us:

One thing to note here is the ratio of r0 to r1. 180:60 or 3:1. And we got three nodes. If we change r1 to 45, making the ratio 4:1, we get four nodes. (I changed d to 45 as well to match r1.)

And here, the ratio is 12:1:

If the second number in the ratio is not one, what happens? Let’s set r0 to 150 and r1 to 100. Now the ratio is 3:2.

Not surprisingly, we get one and a half nodes. To finish this curve, we’d have to go around again, having t go from 0 to PI * 4. Then we get:

If you use whole numbers, the cycle will always complete eventually. In the next example, the ratio was 11:7, so I had to go to PI * 14:

And 111:70, meaning I had to go up to 140 * PI:

It’s somewhat of a pain to work out the math manually, but you can create a function to simplify the fraction and use the denominator * PI * 2 as your for loop count. Here’s a couple of functions that will help you:

function gcd(x, y) {
  result min(x, y)
  while (result > 0) {
    if (x % result == 0 && y % result == 0) {
      break
    }
    result--
  }
  return result
}

function simplify(x, y) {
  g = gcd(x, y)
  return x / g, y / g
}

Just put r0 and r1 into the simplify method, to get the simplified ratio. Take the second number times 2 * PI for your for loop limit. Here’s an example…

width = 800
height = 800
canvas(width, height)

translate(width/2, height/2)

r0 = 111
r1 = 70
d = 70
res = 0.01
num, den = simplify(r0, r1)

circle(0, 0, r0)
stroke()

for (t = 0.0; t < PI * 2 * den; t += res) {
  x = (r0 + r1) * cos(t) - d * cos(((r0 + r1) * t) / r1)
  y = (r0 + r1) * sin(t) - d * sin(((r0 + r1) * t) / r1)	lineTo(x, y)
}
stroke()

In line 11, I get the simplified fraction and use the denominator in line 16 to make sure we loop enough times to complete the curve.

So far, all of these have been epicycloids. Let’s change that d parameter to make some other epitrochoids.

Making d larger than r1:

And smaller:

I stopped drawing the inner circle in these cases.

Well, you don’t need me to supply you with more examples. Just change the numbers and see what comes up.

Special Epitrochoids

A Limaçon is an epitrochoid where the two circles have the same radius. If the points is exactly on the circle, it’s a special Limaçon called a Cardioid.

Here’s a Limaçon where both radii are 80 and the value of d is 160.

And a cardioid where all three are 80.

A nephroid is an epitrochoid where the fixed circle’s radius is twice the radius of the moving circle and the point is on the edge of the moving circle. Here’s a nephroid:

Of course you can make the point distance more or less than the radius, but it’s no longer technically a nephroid at that point. Still a nice curve though:

Now let’s move on to hypotrochoids!

Hypotrochoids

As mentioned, this is when the moving circles is rolling along the inside of the fixed circle.

And if we trace the curve that dot draws, we have a hypotrochoid:

In fact, because the point is on the edge of the circle, we also have a hypocycloid.

When the point moves out, we get something like this:

And when we move it further in…

As before, the animation was made in a brute force way with multiple steps, but there is a one-shot formula.

The formula for a hypotrochoid:

x = (r0 - r1) * cos(t) + d * cos(((r0 - r1) * t) / r1)
y = (r0 - r1) * sin(t) - d * sin(((r0 - r1) * t) / r1)

The parameters are the same as for epitrochoids. In fact the formula is almost exactly the same, just a different sign on some of the terms.

Here’s how we can use it.

width = 800
height = 800
canvas(width, height)

translate(width/2, height/2)

r0 = 300
r1 = 50
d = 50
res = 0.01
num, den = simplify(r0, r1)

for (t = 0.0; t < PI * 2 * den; t += res) {
  x = (r0 - r1) * cos(t) + d * cos(((r0 - r1) * t) / r1)
  y = (r0 - r1) * sin(t) - d * sin(((r0 - r1) * t) / r1)
}
stroke()

Running this code gives us:

Note that the ratio of r0 to r1 is 6:1, and we have six points. This works the same way as it did for epitrochoids.

Here’s one with the radii of 312 and 76:

Here’s the same one with the d point moved out:

And moved in…

I’m sure you can build something fun to explore these with.

Special Hypotrochoids

I’ll mention just a couple of special, named hypotrochoids. Again, they have to do with ratios between the two circles.

A deltoid is a hypocycloid with a circle ratio of 3:1.

And an astroid has a ratio of 4:1.

There’s one more interesting named ratio, that of 2:1. This is called a Tusi couple, named after the 13th century Persian astronomer who first described it. Here’s an animation of it.

As you can see, the path formed is a straight line.

If you make a bunch of these, each a bit out of phase with the last one, you get a picture like this:

Each dot is plainly moving back and forth in a straight line. But if you remove all those lines and inner circles, you get quite an illusion of a rotating circle. In fact, it looks like yet another hypocycloid!

Spirograph!

If you were a nerdy kid like me, growing up you had one of these in your house:

This one I bought brand new a few weeks ago in preparation for this chapter. It came in a very cool tin box with some drawing paper and an instruction/inspiration booklet.

It’s got one large gear which is fixed to the paper. In the old days with small pins, and now with some sticky-tack stuff (a great improvement!) And a bunch of smaller gears with holes in them. You put a pen in one of the holes and trace it around in a circle till you are back where you started.

Voila! A hypotrochoid!

Or, you put the smaller gear on the outside and get an epitrochoid. Great fun, though to be honest, I enjoy doing it in code a lot more! Here, the pen slipped and messed up my epitrochoid before I could finish it. Also, the pen was rather skippy.

As I was playing with this, I realized that you can only make “curtate” epi- and hypotrochoids. The drawing point is always less than the radius of the moving gear.

Coding Curves 08: Bézier Curves

misc

Chapter 8 of the Coding Curves Series

I had to hold myself back here. Bézier curves are fun to program, fascinating to explore, and you can go down a deep hole in explaining how they are constructed and what the formulas mean. The thing is, I’ve already done that a few times. In books and in videos. Here’s a couple of my own videos you might want to check out to learn more:

https://www.youtube.com/watch?v=dXECQRlmIaE: Coding Curves 08: Bézier Curves https://www.youtube.com/watch?v=2hL1LGMVnVM: Coding Curves 08: Bézier Curves

And here are a couple of amazing videos by Freya Holmer:

https://www.youtube.com/watch?v=aVwxzDHniEw: Coding Curves 08: Bézier Curves https://www.youtube.com/watch?v=jvPPXbo87ds: Coding Curves 08: Bézier Curves

So I’m going to limit myself here to the bare basics, some functions, and some cool, practical tips and tricks I’ve discovered over the years.

The Basics

A Bézier curve is defined by two end points and one or more control points. It starts at one end point, curves towards (but not through) the control point(s) and ends at the other end point. By moving any of the points, you alter the shape of the curve. These curves are generally visually pleasing and are used in all kinds of design tools, and are a key part of the shapes of things from fonts to cars.

There are two types of Bézier curves that you’ll run into most often.

Quadratic Bézier Curves

These are defined by two end points and a single control point. Here’s an example:

The control point is the one near the bottom of the canvas. If I move that to the right, it changes the curve:

The lighter lines and black points I just threw in for visual context.

Cubic Bézier Curves

Cubic Bézier Curves have two end points and two control points. Example:

It is possible to have higher order Bézier curves with more control points, but the math gets more and more costly. See Freya’s videos above for some explanations about this.

Most drawing apis have methods for both quadratic and cubic curves, but how they are names varies greatly.

I’ve seen quadratic Bézier curve methods named:

  • curveTo
  • quadraticCurveTo

And cubic Bézier curve methods named:

  • curveTo
  • cubicCurveTo
  • bezierCurveTo

So make sure you know what the method names are for your api. A common strategy, as seen in the examples above, is to have the starting point defined by using a moveTo, or the last known position of the drawing cursor, and then have the curve method just define the control points and the final end point. So you’d do something like:

moveTo(100, 100)
cubicCurveTo(200, 100, 200, 500, 100, 300)
stroke()

But some apis may have other methods that allow you to specify all the points at once.

That’s about it for the basics and built-ins, but of course, we’ll now leave the apis behind and code up some curves ourselves.

Coding Bézier Curves

We’ll start with quadratic curves and then move on to cubic. But before we create the methods that draw the paths, we’re going to create another, more basic method. This will give us the point at any interval along the Bézier curve.

Quadratic

Interestingly, the basic formulas for Bézier curves are one-dimensional. To make two-dimensional, three-dimensional or higher Bézier curves, you just apply the formula once for each dimension. We’ll be sticking two 2D here, so we’ll be doing this twice. The single parametric formula is:

x = (1 - t) * (1 - t) * x0 + 2 * (1 - t) * t * x1 + t * t * x2

Here, x0, x1, and x2 are the end and control “points” and t is a value that ranges from 0.0 to 1.0. This returns the x value along the Bézier path corresponding to the value of t. When t is 0, x is equal to x0. When t is 1, x is equal to x2. When t is between 0 and 1, x will be interpolated.

So to make a 2D quadratic Bézier point function we do this:

function quadBezierPoint(x0, y0, x1, y1, x2, y2) {
  x = (1 - t) * (1 - t) * x0 + 2 * (1 - t) * t * x1 + t * t * x2
  y = (1 - t) * (1 - t) * y0 + 2 * (1 - t) * t * y1 + t * t * y2
  return x, y
}

You can do this if your language lets you return multiple values. Otherwise, you’ll have to encode x, y in come kind of point object.

Note that we have a lot of duplication going on there. We can clean it up first factoring out all those 1-ts:

function quadBezierPoint(x0, y0, x1, y1, x2, y2, t) {
  m = (1 - t)
  x = m * m * x0 + 2 * m * t * x1 + t * t * x2
  y = m * m * y0 + 2 * m * t * y1 + t * t * y2
  return x, y
}

And then:

function quadBezierPoint(x0, y0, x1, y1, x2, y2, t) {
  m = (1 - t)
  a = m * m
  b = 2 * m * t
  c = t * t
  x = a * x0 + b * x1 + c * x2
  y = a * y0 + b * y1 + c * y2
  return x, y
}

If nothing else, this makes it much easier to read.

Now that we have this, we can make a function to draw quadratic Bézier curves. To make it perfectly clear, I’ll name this quadCurve and the cubic one will be cubicCurve

function quadCurve(x0, y0, x1, y1, x2, y2, res) {
  moveTo(x0, y0)
  for (t = res; t < 1; t += res) {
    x, y = quadBezierPoint(x0, y0, x1, y1, x2, y2, t)
    lineTo(x, y)
  }
  lineTo(x2, y2)
}

To make sure we explicitly start and end on the start and end points, we’ll start with an explicit moveTo to the first point and end with an explicit lineTo to the last point. The function takes a res parameter that lets you know how many steps to take along the curve. We’ll start t out equaling res because we’ve already moved to the first point, which is what you’d get if t was 0. In the middle of all that, we just get the point corresponding to the current t and draw a line to it.

Of course, you can make a quadCurveTo method by dropping the first two parameters and the moveTo. This will rely on the user using their own moveTo to specify the starting point of the curve (or continuing it from an existing path). Example of this in use:

canvas(800, 800)
quadCurve(100, 100, 200, 700, 700, 300, 0.01)
stroke()

Which gives us:

If we change the res to something larger like 0.1, we see that things get a bit chunky:

So here you’ll have to experiment with a resolution value that works well. Given that the built-in Bézier methods already figure out a good resolution for you, this curve method has questionable value. But, it got us to write the quadBezierPoint function, which has a LOT of value, as we’ll see.

One thing that the point function can do that the built in methods can’t is animation. For this section, like I did in previous chapters, I’m going to assume you have or can make some kind of function that runs repeatedly and can create animations. I’m going to call it loop. What I’m going to do here is rather than drawing the curve from a t of 0 to 1, I’ll have it go from 0 to a value finalT that will change over time.

canvas(400, 400)
x0 = 50
y0 = 50
x1 = 150
y1 = 360
x2 = 360
y2 = 150
finalT = 0
dt = 0.01
res = 0.025

function loop() {
  clearCanvas()
  moveTo(x0, y0)
  for (t = res; t < finalT; t += res) {
    x, y = quadBezierPoint(x0, y0, x1, y1, x2, y2, t)
    lineto(x, y)
  }
  stroke()

  // add to finalT
  finalT += dt

  // if we go past 1, turn it around
  if (finalT > 1) {
    finalT = 1
    dt = -dt
  } else if (finalT < 0) {
    // if we go past 0, turn it back
    finalT = 0
    dt = -dt
  }
}

And this should give you an animation that looks something like this:

Here, the for loop is going from res to finalT so it doesn’t draw the entire curve (unless finalT equals 1). Then we change finalT by adding dt to it. This brings finalT closer and closer to 1, so the curve is drawn more and more fully. Eventually finalT will go beyond 1 so we set it back to 1 and make dt negative, which reverses the whole process until finalT goes below 0, where we bounce it back the other way.

Rather than just drawing a line, we can animate an object along a Bézier path now! Here’s the code for that. It should be pretty clear. I’ll just include the loop function. The rest should be the same.

function loop() {
  clearCanvas()

  x, y = quadBezierPoint(x0, y0, x1, y1, x2, y2, finalT)
  circle(x, y, 10)
  fill()

  // no changes beyond here...
  // add to finalT
  finalT += dt

  // if we go past 1, turn it around
  if (finalT > 1) {
    finalT = 1
    dt = -dt
  } else if (finalT < 0) {
    // if we go past 0, turn it back
    finalT = 0
    dt = -dt
  }
}

Now we’re just getting the x, y point for the current value of finalT and drawing a circle there. This assumes you have a circle drawing function. You can use the one we created in Chapter 3 if you need one.

In this example, I drew a light line for the same quadratic curve using the built-in method of my api, just to show that we’re on track with the standard definitions of these things.

OK, let’s take a break here and jump over to cubic curves

Cubic

Pretty much everything I said above for quadratic curves is going to apply to cubics. It’s just a different formula – a bit more complicated. Here is is for one dimension:

x = (1 - t) * (1 - t) * (1 - t) * x0 + 3 * (1 - t) * (1 - t) * t * x1 + 3 * (1 - t) * t * t * x2 + t * t * t * x3

And the 2D function:

function cubicBezierPoint(x0, y0, x1, y1, x2, y2, x3, y3, t) {
  x = (1 - t) * (1 - t) * (1 - t) * x0 + 3 * (1 - t) * (1 - t) * t * x1 + 3 * (1 - t) * t * t * x2 + t * t * t * x3
  y = (1 - t) * (1 - t) * (1 - t) * y0 + 3 * (1 - t) * (1 - t) * t * y1 + 3 * (1 - t) * t * t * y2 + t * t * t * y3
  return x, y
}

Yikes! That’s a mess. Let’s clean it up again by factoring out the 1 - ts:

function cubicBezierPoint(x0, y0, x1, y1, x2, y2, x3, y3, t) {
  m = 1 - t
  x = m * m * m * x0 + 3 * m * m * t * x1 + 3 * m * t * t * x2 + t * t * t * x3
  y = m * m * m * y0 + 3 * m * m * t * y1 + 3 * m * t * t * y2 + t * t * t * y3
  return x, y
}

That’s a little better. One more step to clean it up:

function cubicBezierPoint(x0, y0, x1, y1, x2, y2, x3, y3, t) {
  m = 1 - t
  a = m * m * m
  b = 3 * m * m * t
  c = 3 * m * t * t
  d = t * t * t
  x = a * x0 + b * x1 + c * x2 + d * x3
  y = a * y0 + b * y1 + c * y2 + d * y3
  return x, y
}

Much better!

Now we can make a cubicCurve function.

function cubicCurve(x0, y0, x1, y1, x2, y2, x3, y3, res) {
  moveTo(x0, y0)
  for (t = res; t < 1; t += res) {
    x, y = cubicBezierPoint(x0, y0, x1, y1, x2, y2, x3, y3, t)
    lineTo(x, y)
  }
  lineTo(x2, y2)
}

That was easy. No explanation or example needed I think.

Now for your assignment: adjust the animations we did above to work for cubic curves. It’s really just a matter of adding a new x3, y3 point and calling the new function.

That’s the basic code for implementing Bézier curves and paths. But there are a few more neat tricks for you.

Drawing Through a Point

At some point, everyone who starts coding Bézier curves is going to say:

This is neat, but I want it to go THROUGH the control point(s).

Me – some time around the year 2000.

Well, we can do that! It’s pretty easy for quadratic curves. What you need to do is create another control point that will pull the curve even further out so that it just goes through the original control point. And that new point is pretty simple to calculate. Where the points are x0, y0, x1, y1, x2, y1, the new control point will be:

x = x1 * 2 - x0 / 2 - x2 / 2
y = y1 * 2 - x0 / 2 - x2 / 2

Now we can make a function, let’s call it quadCurveThrough that implements this behavior. It just gets this new point and uses the built-in function to draw the curve. I’m going to postulate that’s called quadraticCurveTo on your system, but it might be something else.

function quadCurveThrough(x0, y0, x1, y1, x2, y2) {
  xc = x1 * 2 - x0 / 2 - x2 / 2
  yc = y1 * 2 - y0 / 2 - y2 / 2
  moveTo(x0, y0)
  quadraticCurveTo(xc, yc, x2, y2)
}

Here I’ve drawn a regular quadratic curve in red and one using this function in blue. And I drew in the points to prove that it does what I claim it does!

Your next question is how to do the same thing for cubic curves. I don’t have an answer for that one yet. But I will keep digging. I guess there’s a chance that someone will comment the answer here, or tell me it’s not possible. 🙂

Piece-wise Quadratic Bézier Curves

The other question people commonly ask is:

How do I make a Bézier curve with N control points (where N is 3 to infinity)?

Also me around the same time as earlier question

As mentioned earlier, this is mathematically possible, but it gets prohibitively expensive when you move past cubic. That’s why you’ll probably never see quartic or quintic Bézier curve functions. But it’s still a very useful thing to have a smooth curve with an arbitrary number of control points. And of course you’ve seen this kind of thing any time you’ve used a pen tool in a drawing program.

In the splines video above (second one by Freya), she shows making a longer curve by piecing together multiple cubic Bézier curves.

https://www.youtube.com/watch?v=jvPPXbo87ds: Coding Curves 08: Bézier Curves

These are sometimes called piecewise Bézier curves. I’m going to show you something a bit simpler, using only quadratic curves. It’s not too hard to implement and it supports any number of control points you want to throw at it. I’ll even show you a version that creates a closed loop.

This technique is actually covered in depth in the second video I posted above (by me):

https://www.youtube.com/watch?v=2hL1LGMVnVM: Coding Curves 08: Bézier Curves

So I’m not going to do too deep a dive here, but I’ll go over the basics and give you some code and examples.

The basic principle is you’re going to first create a new point that is mid-way between p0 and p1. Call that pA. And another that’s between p1 and p2. Call that pB. Draw a line from p0 to pA, then draw a quadratic curve using pA, p1 and pB.

Then you’ll get midpoint from p2 and p3, called pC and make a curve from pB through p2 and ending at pC.

You’ll continue that until the end where you draw a curve from the second-to-last midpoint, through the second-to-last point, ending at the last midpoint. Then finally a line from the last midpoint to the last point.

And here’s the curve:

The code for this can get a little tricky, but having worked through it a few times, I’m pretty happy with a method that looks like the following. Note, that due to the potentially large number of parameters that will be passed in, it really helps to have some kind of point object. Whether that’s a class, structure, or generic object with x and y properties.. up to you and your language. This method will take an array of these point objects. The code assumes the array has a length property, but there might be something different in your language, like a len method.

function multiCurve(points) {
  // line from the first point to the first midpoint.
  moveTo(points[0].x, points[0].y)
  midX = (points[0].x + points[1].x) / 2
  midY = (points[0].y + points[1].y) / 2
  lineTo(midX, midY)

  // loop through the points array, starting at index 1
  // and ending at the second-to-last point
  for (i = 1; i < points.length - 1; i++) {
    // find the next two points and their midpoint
    p0 = points[i]
    p1 = points[i+1]
    midX = (p0.x + p1.x) / 2
    midY = (p0.y + p1.y) / 2

    // curve through next point to midpoint
    quadraticCurveTo(p0.x, p0.y, midX, midY)
  }

  // we'll be left at the last midpoint
  // draw line to last point
  p = points[points.length - 1]
  lineTo(p.x, p.y)
}

That method seems long, but I added a lot of comments to each part.

For an example, I generated a half dozen random points. I don’t know how you’d do that on your system, so I’m going to say we have a function randomPoint(xmin, ymin, xmax, ymax). (Actually I do have such a method in my own library!) Once you have your points array, you just pass that array to your new function and stroke the resulting path:

context(800, 800)
points = []
for (i = 0; i < 6; i++) {
  points.push(randomPoint(0, 0, 800, 800))
}

multiCurve(points)
stroke()

The glorious result:

Quite nice. I also drew all the points to give a bit more context on why the curve looks like it does.

Closed Curves

The final thing in this section will be to alter the method to allow for a closed loop. Mainly that’s going to be getting rid of the starting and final line segments, and figuring out the curve from the end of the curve back to the beginning.

function multiLoop(points) {
  // find the first midpoint and move to it.
  // we'll keep this around for later
  midX0 = (points[0].x + points[1].x) / 2
  midY0 = (points[0].y + points[1].y) / 2
  moveTo(midX0, midY0)

  // the for loop doesn't change
  for (i = 1; i < points.length - 1; i++) {
    p0 = points[i]
    p1 = points[i+1]
    midX = (p0.x + p1.x) / 2
    midY = (p0.y + p1.y) / 2
    quadraticCurveTo(p0.x, p0.y, midX, midY)
  }

  // we'll be left at the last midpoint
  // find the midpoint between the last and first points
  p = points[points.length - 1]
  midX1 = (p.y + points[0].x) / 2
  midY1 = (p.y + points[0].y) / 2

  // curve through the last point to that new midpoint
  quadraticCurveTo(p.x, p.y, midX1, midY1)

  // then curve through the first point to that first midpoint you saved earlier
  quadraticCurveTo(points[0].x, points[0].y, midX0, midY0)
}

To explain a bit more here’s the first step and the for loop:

We move to the first midpoint, then loop through the rest, finding midpoints and doing quadratic curves. This leaves us at the last midpoint. Then…

We find the midpoint between the first and last points, and execute the last two curves to close up the shape. The following image is made with the same setup as the last one, but calling multiLoop instead of multiCurve (and different random points).

These are two of my favorite functions, and I’m happy to share them with you.

Even Distribution

The final trick I want to share has to do with evenly distributing objects on a quadratic curve. One use case for this is when you want to have text follow a curve. You want to be able to space the letters evenly on the curve. You’ll also want to rotate them to follow the curve, but that’s beyond the scope of this article.

At first this might seem like a trivial problem. You have your t value that you use to divide up the curve. If you want, say, 20 object spaced out on the curve, just have them 0.05 apart. 20 x 0.05 = 1.0. Work done. Well, let’s try it.

canvas(800, 800)

x0 = 100
y0 = 700
x1 = 100
y1 = 100
x2 = 700
y2 = 400

moveTo(x0, y0)
quadraticCurveTo(x1, y1, x2, y2)
stroke()

// 20 evenly spaced t values (21 counting the end one)
for (t = 0; t <= 1; t += 0.05) {
  x, y = quadBezierPoint(x0, y0, x1, y1, x2, y2, t)
  circle(x, y, 6)
  fill()
}

Here’s what that gives us.

Not evenly spaced at all. The ones near the ends are spaced out and in the middle they are closer together. That’s just the way it is with Bézier curves. So, we have to figure out a way to get those points evenly spaced.

Sadly, there’s not a very easy way to do this. I’m going to give you a horribly un-optimized, brute force way to get decent results, and a couple of hints on how to make it better.

So, to get evenly spaced points along a curve, it makes sense that you need to know the length of the curve. If the length is 200 pixels, and you want 20 points, you put one point every 10 pixels along the length of the curve.

Surprisingly, there is no simple formula to get the length of a Bézier curve. But we can get pretty close by sampling a bunch of points along the curve and getting the distance between each pair. This would look something like this:

function quadBezLength(x0, y0, x1, y1, x2, y2, count) {
  length = 0.0
  dt = 1.0 / count
  x, y = x0, y0
  for (t = dt; t < 1; t += dt) {
    xn, yn = quadBezierPoint(x0, y0, x1, y1, x2, y2, t)
    length += distance(x, y, xn, yn)
    x, y = xn, yn
  }
  length == distance(x, y, x2, y2)
  return length
}

Here, count is how many samples we want to take. The more samples, the more accurate we’ll be.

Then, dt is the amount to increase t by as we loop through the curve.

We keep track of the last point, x, y, which will start as x0, y0. Then we loop through the curve getting each new point, xn, yn and finding the distance between the last point and the new point, then making the new point the last point. I’m not going to show you how to find the distance between two points, just assuming that you have a function for that. You add that distance to the accumulating length.

Then you do one last length addition for x2, y2. Then return the length.

Make sure that all makes sense first, because I’m going to throw something else in there.

It’s going to be very useful to keep track of what the length was at each point as we work through the curve. So we’re going to store each successive value in an array. And rather than returning the total length, we’ll return the array.

function quadBezLengths(x0, y0, x1, y1, x2, y2, count) {
  lengths = []
  length = 0.0
  dt = 1.0 / count
  x, y = x0, y0
  for (t = dt; t < 1; t += dt) {
    xn, yn = quadBezierPoint(x0, y0, x1, y1, x2, y2, t)
    length += distance(x, y, xn, yn)
    lengths.push(length)
    x, y = xn, yn
  }
  length == distance(x, y, x2, y2)
  lengths.push(length)
  return lengths
}

Now the full length of the curve is in the last element, but we have a whole bunch of other sub-lengths too. Here’s what we do.

count = 500
lengths = quadBezLengths(x0, y0, x1, y1, x2, y2, count)
length = lengths[count-1]

for (i = 0.0; i <= 1; i += 0.05) {
  // the length of the curve up to the next point
  targetLength = i * length

  // loop through the array until the length is higher than the target length
  for (j = 0; j < count; j++) {
    if (lengths[j] > targetLength) {
      // t is now the percentage of the way we got through the array.
      // this is the t value we need to get the next point
      t = j / count

      // get the point and draw the next circle.
      x, y = quadBezierPoint(x0, y0, x1, y1, x2, y2, t)
      circle(x, y, 6)
      fill()
      break
    }
  }
}

OK, a bit complex, but let’s go through it. We get the lengths using a count of 500, and capture the total length.

The we loop through from 0 to 1 by 0.05, like we were doing before. But rather than using that as the t value of the Bézier curve, we use it to find a fraction of the length of the curve. Say the curve was 500 pixels long and i was 0.5, then the target length we are looking for is 250 for the next point.

Now we loop through the array with j and get the length values until we go above our target length. If we divide j by count, we wind up with the t value that created this particular length. We plug that back into the Bézier point function to get the next point and draw it. We also want to break out of the inner loop at that point, since we are done and can move on to the next point.

This gives us the following:

We missed the last point (because we never went beyond that length), but we could just draw another point at x2, y2. This is pretty close to evenly spaced. The higher the value you use for count, the more accurate it will be. If we move count down to 100, we can see it’s off:

Now, there is a LOT wrong with this code. Mostly in its optimization.

First of all, rather than looping through from 0 on each point, we could use a binary search.

Secondly, we have to add a lot of points to get any kind of accuracy. Because we are just grabbing the last point we find in that inner loop. Instead of taking a predefined point, which we know is a bit too high, we could interpolate between that point and the previous point.

For example, say our target length was 150. And at index 87 we got a length of 160. We look back at index 86 and we find 140. OK, so we want a value half way in between 86 and 87. So instead of calculating t as 87 / count, or even 86 / count, we interpolate 50% and say 86.5 / count. It still won’t be perfect, but you can use a lower count and still get good results.

I’m going to leave all of that as an exercise for you.

If you want more info on this technique, and a fuller explanation, check out this site:

http://www.planetclegg.com/projects/WarpingTextToSplines.html

Summary

So there, you go, a few basics, a few tips, a few tricks. Until next time…

Goodbye Twitter

misc

I’ll keep this brief. I’m leaving Twitter after 16 years. A tough decision. 8000 followers. Lots of history. It was an overwhelmingly positive experience. But the platform is in a very different place now and I just don’t want to have any connection to it going forward. I’ll disable my account some time on December 13, 2022 probably.

You can now find me over at Mastodon at bit101 (@bit101@mstdn.social) – Mastodon 🐘

I’ve flirted with Mastodon for several years now. There is now a critical mass of people I know there. I’m having a great time there and it feels very much like the early days of Twitter. Good conversation, good creativity, low drama.

I met so many people and maintained so many good relationships on Twitter. I hope to continue to stay in touch with a lot of you, either here, on Mastodon, or on other avenues.

Coding Curves 07: Parabolas

coding curves

Chapter 7 of the Coding Curves Series

I’ll admit that when I started brainstorming this chapter, I began to second guess the decision to include it at all. The curve we’ll cover here, the parabola, is fairly simple, basic, even plain compared to the crazy curves we’ve been generating in the last three chapters. But as I got into it I realized that this curve is pretty cool and has lots of interesting properties. In fact, I was going to cover another curve, the hyperbola, in this chapter as well, but I got so deep into parabolas that I decided to save hyperbolas for another time.

Parabolas

Just to make sure we’re on the same page here, this is what we’re talking about when we’re talking about parabolas:

One of the first things you’ll learn about parabolas is that they are a type of conic section. OK. But since we’re just doing two-dimensional curves here, their relation to a three-dimensional form is a bit irrelevant. But I’ll throw this image in here anyway and be done with the subject:

From https://en.wikipedia.org/wiki/Parabola

Since we’ve already covered circles and ellipses, this chapter and the next will fill out all of the conic sections. That was never a goal, but I’ll never pass up an achievement gained.

In the first picture above, note that the parabola is mirrored across the y-axis, and that the vertex (the highest or lowest point, depending on which way it opens) is just touching the x-axis.

As with most of these geometric shapes, there are various formulas used to describe them. Here’s a fairly simple and usable one:

y = a * x * x

In this case, a is a parameter that controls how wide the parabola opens up – and whether it opens on the top or bottom. So that’s pretty simple. Let’s draw one. But first, let’s set things up correctly and even create a couple of utility functions that will help us do so.

Set Up

As you can see in the initial image, we want the origin of the space to be in the center of the canvas so we can clearly see the parabola. So we’ll want to translate the canvas there.

Also, there’s a very good chance that your drawing api will be drawing things upside down, compared to normal Cartesian coordinates. In other words, y values will increase positively towards the bottom of the screen and negatively towards the top. So it will be good to flip it the other way around.

Finally, it will be nice if we can actually see the axes sometimes. We can create a function to draw them easily.

To center the canvas, you’ll use something like the following:

translate(width / 2, height / 2)

This is assuming your drawing api has transformation functions built in. It’s possible to do everything here without them, but it would be cumbersome to try to explain things twice, so I’m just going to assume that your drawing api, like most decent ones, does have transformation built in.

For the above, we can simplify that by making a center function:

function center() {
  translate(width / 2, height / 2)
}

As always, this is making the pseudocode assumption that your canvas or surface or whatever object you draw things on is a global object and you can just call these functions like that. Processing works that way, but for may other systems, these functions will be methods on a canvas object, so it might be more like

canvas.translate(canvas.width / 2, canvas.height / 2)

But I’m sure you can figure all that out.

To flip the y-axis, just do a scale:

scale(1, -1)

This keeps the x-axis as is, and flips the y-axis. I know, I know. Sometimes in this series I flip the canvas, other times I do not. Once again, these posts are meant to be practical drawing tutorials, not necessarily mathematically rigorous. So you should know how to flip the canvas and do it when and if you think you should.

Finally, I like to have a function around to draw the axes. For the purposes here, we can do something like:

function drawAxes() {
  lineWidth = 0.25
  moveTo(-width, 0)
  lineTo(width, 0)
  moveTo(0, -height)
  lineTo(0, height)
  stroke()
  lineWidth = 1
} 

This draws horizontal and vertical lines well beyond the bounds of the canvas, but that’s usually ok. It also sets the line width very low so it’s just a hint of a line, and then resets it to 1 when it’s done. In your actual function, you’ll probably want to do something like pushing and popping or saving and restoring the state of the canvas, so that your line width goes back to what it was before calling the function, which may not have been 1.

OK, now we can set up our canvas like so:

width = 800
height = 800
canvas(width, height)

center()
scale(1, -1)
drawAxes()

That gets you here, which is a good start:

You might want to combine all that stuff into a more comprehensive setup function. Up to you.

Drawing the Parabola

Building on the last code we can now loop x from -width/2 to width/2 to go across the canvas from left to right. We’ll use a low value of a like 0.003 because we are operating on pixel values that are in the several hundreds. I just chose that value because it worked visually.

a = 0.003
for (x = -width / 2; x <= width / 2; x++) {
  y = a * x * x
  lineTo(x, y)
}
stroke()

And here’s our parabola:

If we change a to something larger, like 0.3, we get a much narrower parabola:

And something smaller like 0.0003, we get a much wider one:

Technically, a should not be 0 in the parabola formula. But if you do try it, you’ll just get a line. Then when a goes negative, the parabola opens up on the opposite side. Here’s -0.003:

And that’s all there is to parabolas.

Oh, no, wait. There’s a few more things!

Focus and Directrix

Parabolas have what is known as a focus point. That point is defined as:

x = 0
y = 1 / (4 * a)

Let’s draw that point:

a = 0.003
for (x = -width / 2; x <= width / 2; x++) {
  y = a * x * x
  lineTo(x, y)
}
stroke()

// draw focus
focusX = 0
focusY = 1/(4 * a)
circle(focusX, focusY, 4)
fill()

And it has another property called the directrix. This is a horizontal line where the y value is:

y = -1 / (4 * a)

We can draw that:

a = 0.003
for (x = -width / 2; x <= width / 2; x++) {
  y = a * x * x
  lineTo(x, y)
}
stroke()

// draw focus
focusX = 0
focusY = 1/(4 * a)
circle(focusX, focusY, 4)
fill()

// draw directrix
directrixY = -1/(4 * a)
moveTo(-width / 2, directrixY)
lineTo(width / 2, directrixY)
stroke()

What’s obvious here is that the distance from that the vertex to the focus is the same as the distance from the vertex to the directrix. Obvious because the formula for the y value of both is the same but reversed in sign.

But here’s an interesting fact – that statement of equal distance holds true for any point on the parabola! We can show that. Building on top of the last code sample…

lineWidth = 0.5
for (x = -width / 2; x <= width/2; x += 40) {
    // find a point on the parabola
	y = a * x * x
	circle(x, y, 4)

    // draw a vertical line from the directrix to that point
    // then from the that point to the focus
	moveTo(x, directrixY)
	lineTo(x, y)
	lineTo(0, focusY)
}

We sample the parabola at a number of points. We draw the point and then draw a line up from the directrix, and then over to the focus. Both lines emanating from each point will be the same length. You may not have a practical use for this straight off, but it looks neat anyway!

Tangent Line

Another thing you can do is find the tangent line at any point on the parabola. This will represent the slope of the curve at that point. The formula for the tangent at point x0, y0 is:

y = 2 * a * x0 * x - a * x0 * x0

This might be a bit confusing because we have x and x0. Again, x0 is a fixed point on the parabola. And x is one of the points the defines the tangent line. The formula gives you the y for that x. Let’s code that up for a single point on the curve. We’ll draw the parabola, then choose an x0, y0 point and use the formula to find two other points that make up the line.

// draw the parabola
a = 0.003
for (x = -width / 2; x <= width / 2; x++) {
  y = a * x * x
  lineTo(x, y)
}
stroke()

// find a point on the parabola
x = -80
y = a * x * x
circle(x, y, 4)
fill()

// find a point on the far left of the canvas
x1 = -width / 2
y1 = 2 * a * x0 * x1 - a * x0 * x0

// and one on the far right
x2 = width / 2
y2 = 2 * a * x0 * x2 - a * x0 * x0

// and draw a line
moveTo(x1, y1)
lineTo(x2, y2)
stroke()

And that gives us this point and line:

Now we can sample various points on the parabola (like we did above) and draw the tangent line for each.

lineWidth = 0.5
for (x0 = -width / 2; x0 <= width/2; x0 += 40) {
    // find a point on the parabola
	y0 = a * x0 * x0
	circle(x0, y0, 4)
    fill()

    // find a point on the far left of the canvas
    x1 = -width / 2
    y1 = 2 * a * x0 * x1 - a * x0 * x0

    // and one on the far right
    x2 = width / 2
    y2 = 2 * a * x0 * x2 - a * x0 * x0

    // and draw a line
    moveTo(x1, y1)
    lineTo(x2, y2)
    stroke()
}

This is pretty much all the same code just moved inside the loop. But then we get this:

You could remove the original parabola and tighten up the interval and have some nice pseudo-string-art.

Parabolic Mirrors

One thing you’ll always read about with parabolas is that rays coming straight in to the parabola will all reflect onto a single point. This is used in antennas and radio telescopes, to focus the received signal on the receiver, and in various solar devices to focus the rays of the sun onto a single point (which becomes incredibly hot). Not surprisingly, the point they are focused on is the focus point.

https://en.wikipedia.org/wiki/Parabolic_antenna#/media/File:Erdfunkstelle_Raisting_2.jpg

In fact, I remember when I was a kid, my step-father had a little solar cigarette lighter like this:

While more of a novelty than something you’d use day to day, it actually worked.

If you wanted to do all the math, you could find out the point where an incoming ray hits the parabola, find the tangent line for that point. Then find the normal at that point (the vector perpendicular to the tangent line) and reflect the incoming ray across that normal. If you did all that, you’d see that the reflected ray hits the focus point. I’m not going to do that exercise with you, but let’s draw some of these rays and you should be able to see that it looks correct. We’ll just draw some rays from the top of the canvas, down to where they hit the parabola and then over to the focus point.

a = 0.003

// assume code for drawing the parabola is here...

// draw the focus
focusX = 0
focusY = 1/(4 * a)
circle(focusX, focusY, 4)
fill()

lineWidth = 0.5
for (x = -width / 2; x <= width/2; x += 40) {
    // find a point on the parabola
	y = a * x * x
	
    moveTo(x, -height / 2)
    lineTo(x, y)
    lineTo(focusX, focusY)
    stroke()
}

You can follow any ray down to the parabola and from there to the focus point. And you can see it reflects off the curve at a believable angle. Again, I cheated and just drew it directly, but you could work out the physics and get the same thing.

Note that this only works for rays that are parallel to the y-axis in this configuration. Rays coming in at any other angle will not converge on the focus. This is why you have to carefully aim the cigarette lighter at the sun to get enough heat, or aim the satellite dish at the satellite to get a good signal.

Another Formula

Of course, parabolas are not always centered on the y-axis and just touching the x-axis. They are so far because we’re using a very simplified formula. Here’s one that we can do more with:

y = a * x * x + b * x + c

We’ve added a couple new parameters here. It’s pretty clear that the c parameter just gets added on to the rest, so has a direct influence on the final y position of the vertex and the rest of the curve. The b parameter is a bit more complex, so let’s code it up and see what it does. We’ll start by making a parabola function

function parabola(a, b, c, x0, x1) {
	for (x = x0; x <= x1; x++) {
		y = a*x*x + b*x + c
		lineTo(x, y)
	}
}

You might want to improve on this, but this will work for our purposes here. We just loop through from x0 to x1, find a y for each x, and draw a line to it.

Let’s see what different values give us.

parabola(0.01, 0, 0, -width/2, width/2)

This gives us what we’re already used to:

Here’s two more with a c of -200 and +200.

parabola(0.01, 0, -200, -width/2, width/2)
parabola(0.01, 0, 200, -width/2, width/2)

Nothing surprising. Now let’s set c back to 0 and give b a positive value.

parabola(0.01, 3, 0, -width/2, width/2)

And a negative value…

parabola(0.01, -3, 0, -width/2, width/2)

So for a positive a, changing b moves the parabola down and to the left or right. What about for negative a?

parabola(-0.01, -3, 0, -width/2, width/2)
parabola(-0.01, 3, 0, -width/2, width/2)

Not too surprising that they move up now.

Note that for any of the formulas we’ve used so far, we can swap the x’s and y’s and get a parabola that opens to the left or the right.

In this case, the simple formula is:
x = a * y * y

And the deluxe version would be:

x = a * y * y + b * y + c

Summary

Well, that’s all we’ll cover about parabolas for now. I don’t know if you’ll ever have the need to draw a parabola with code, but if you do, you will now be well equipped!

Coding Curves 06: Pintographs

coding curves

Chapter 6 of the Coding Curves Series

Another physical device which renders complex curves that we can simulate!

This one is called a Pintograph, and I’ve actually built one of these myself.

We’ll start with a video – this is literally the first video that came up when I searched Youtube for “pintograph”, but it does the job.

A pintograph can be considered to be a type of harmonograph, but rather than being based on pendulums, pintographs are usually driven by electric motors (though some are hand-cranked). There are disks attached to the motors and arms attached to the disks and a pen attached to the arms. You can change the size of the disks and where the arms are attached, the length of the arms and where they pivot, and the relative speed and offset of the motors to create a bunch of different types of curves.

For a long time I didn’t know where the word “pintograph” came from. I finally discovered it from the person who coined the term. Actually, their daughter coined the term. It comes from a pantograph, and the idea that the spinning wheels look like a Ford Pinto. Read more here: http://www.fxmtech.com/harmonog.html

In case you are not familiar with a pantograph, it’s a device that is often used for copying drawings. It has a few pivoting arms. You pin one of the pivot points so it doesn’t move, then put a pointer in one of the points and a pen in another. As you move the pointer along the original drawing, it moves the pen along the same shape. You can configure it to copy the drawing at the same size, or scale it up or down.

From https://en.wikipedia.org/wiki/Pantograph

The Simulation

The pintograph we’ll simulate will be pretty simple. It will have two rotating disks with one arm attached to each. The arms will be attached to each other on the opposite end and that’s where the virtual pen will be.

First of all, we’ll need to simulate the two disks. They’ll each have an x, y position, a radius, a speed and a phase. We’ll make these visible to start with so you can get a feeling for what’s going on. Otherwise it would just be a long complicated formula.

We’ll create two circles and rotate them and show the point where the arms will be attached.

Each disk will be represented by this kind of structure:

disk: {
  x,
  y,
  radius,
  speed,
  phase,
}

Like last time, it doesn’t matter if this is a generic object, a struct, a class or what. I’m going to assume that you’ll also have a function that will create one of these disks like this:

d0 = disk(100, 200, 100, 2, 0.5)
d1 = disk(400, 200, 60, 3, 0.0)

Now we can set up an animation like so:

width = 600
height = 400
canvas = (width, height)
t = 0

d0 = disk(100, 200, 100, 2, 0.5)
d1 = disk(400, 200, 60, 3, 0.0)

function loop() {
  clearScreen()

  circle(d0.x, d0.y, d0.radius)
  stroke()
  x0 = circle.d0.x + cos(t * d0.speed + d0.phase) * d0.radius
  y0 = circle.d0.y + cos(t * d0.speed + d0.phase) * d0.radius
  circle(x0, y0, 4)
  fill()

  circle(d0.x, d0.y, d0.radius)
  stroke()
  x1 = circle.d1.x + cos(t * d1.speed + d1.phase) * d1.radius
  y1 = circle.d1.y + cos(t * d1.speed + d1.phase) * d1.radius
  circle(x1, y1, 4)
  fill()

  t += 0.1
}

Again, loop is a theoritical function that will be fun over and over again so that we get an animation. I’ve set mine up to create the frames for an animated gif, but this can be done as a real time animation as well. Here’s what I got:

This doesn’t smoothly loop, but that’s ok. You can see that we have two disks of different sizes moving at different speeds. The code itself shouldn’t be too complex. We clear the screen and for each disk we draw a circle at its position and with its radius. Then we take an offset from that circle to calculate the point where the arm will be attached. This uses basic trig: x = cos(angle) * radius, y = sin(angle) * radius. Here, the angle is t times the disk’s speed, plus it’s phase. When we get that point, we fill a smaller circle there.

Attaching the Arms and Finding the Pen

OK. Now things get a bit more mathy. We’re going to have two arms. each one is going to be attached to one of those spinning points at one and, and they’ll be attached to each other on the other end. At this point, a sketch is in order:

We have our two disks. Through their radii and rotation and position, we know the positions of points p0 and p1. That’s what we just did above. We can also define the lengths of those arms. They’ll be the same for now, but they could be different. We’ll call them a0 and a1. (I know, it looks like 90 and 91 in the sketch. Sorry!) What we need to do is get the position of that top point where the two arms join.

We can easily get the distance between p0 and p1 with the Pythagorean theorem. We can call that d. It’s represented by the dotted line in the diagram.

Now we have a triangle whose sides are a0, a1 and d.

There’s a trigonometric rule called the “law of cosines” that will help us here. With it, if you you know the lengths of all three sides of a triangle, a, b, and c, you can get any angle of that triangle. Usually the way it’s written is as follows though:

c = sqrt(a*a + b*b - 2*a*b*cos(y))

… where y is the angle opposite side c. So if you know the length of two sides and the angle between them, you can find the length of the opposite side. Also useful in some cases, but not what we need here.

But we can rearrange that formula and put the unknown variable by itself on one side, and the other side will be the formula you need to calculate that value.

What I want is to know is one of those angles so I can use it to find the location of the pen. Here’s what I came up with, showing the formula, some hand-wavey algebra and the resulting formula we’ll need.

In our case, the sides we have correspond to the sides in the formula like so:

a = a0
b = d
c = a1

So if we apply this formula, we’ll get the angle between the p1, p0, and the pen.

We can also use atan2 to get the angle from p0 to p1. And if we subtract them, we get the actual angle that goes from p0 to the pen.

Again, a drawing is in order:

The big angle that goes from p1 to p0 to the pen, which we get with the law of cosines, we call p1_p0_pen. The smaller angle we get with atan2, we call p0toP1. Subtract them and we have the angle we have to go towards to get the location of the pen. I’m sure there are different ways to do this, and probably much better ones, but this one works. Once you have all the steps, you can further simplify it to something more concise, but I wanted to show the steps to hopefully have it make some sense.

Anyway, we can now code this up:

width = 600
height = 600
canvas = (width, height)
t = 0

d0 = disk(150, 450, 100, 2, 0.5)
d1 = disk(450, 450, 60, 3, 0.0)

function loop() {
  clearScreen()

  circle(d0.x, d0.y, d0.radius)
  stroke()
  x0 = d0.x + cos(t * d0.speed + d0.phase) * d0.radius
  y0 = d0.y + sin(t * d0.speed + d0.phase) * d0.radius
  circle(x0, y0, 4)
  fill()

  circle(d0.x, d0.y, d0.radius)
  stroke()
  x1 = d1.x + cos(t * d1.speed + d1.phase) * d1.radius
  y1 = d1.y + sin(t * d1.speed + d1.phase) * d1.radius
  circle(x1, y1, 4)
  fill()

  // the length of the arms
  a0 = 350
  a1 = 350

  // get the distance between p0 and p1
  dx = x1 - x0
  dy = y1 - y0
  d = sqrt(dx * dx + dy * dy)

  // find the two key angles and subtract them
  p1_p0_pen = acos((a1 * a1 - a0 * a0 - d * d) / (-2 * a0 * d))
  p0toP1 = atan2(y1 - y0, x1 - x0)
  angle = p0toP1 - p1_p0_pen

  // find the pen point
  pX = x0 + cos(angle) * a0
  pY = y0 + sin(angle) * a0
  
  // draw the arms
  moveTo(x0, y0)
  lineTo(pX, pY)
  lineTo(x1, y1)
  stroke()

  t += 0.1
}

With that, you should be able to render something like what you see below.

Again, this is a non-cleanly looping gif, but shows the idea at work. I’ve commented the code to explain what’s going on in each step. Hopefully that helps. Note that I also changed the size of the canvas and moved the disks towards the bottom to make room for the arms.

One thing to be careful of is making sure that the two arms are long enough so that they can always reach between the two connection points on the disks. In the real world, if they were too short, they’t probably just break or jam up the motors. In code, you’ll probably just get a NaN (not a number) error when you go to do the acos and you’ll be sitting there wondering what’s wrong. I speak from experience.

Drawing the Curve

Finally, let’s see what this thing draws. For this, I’m going to abandon the animation and stop drawing the circles and arms. I’ll just track where the pen is on each iteration and use that to draw a long, looping, Lissajous-ish curve.

width = 800
height = 600
canvas = (width, height)

function render() {
  t = 0
  
  d0 = disk(250, 550, 141, 2.741, 0.5)
  d1 = disk(650, 550, 190 0.793, 0.0)

  // the length of the arms
  a0 = 400
  a1 = 400
  
  for (i = 0; i < 50000; i++) {
    x0 = d0.x + cos(t * d0.speed + d0.phase) * d0.radius
    y0 = d0.y + sin(t * d0.speed + d0.phase) * d0.radius

    x1 = d1.x + cos(t * d1.speed + d1.phase) * d1.radius
    y1 = d1.y + sin(t * d1.speed + d1.phase) * d1.radius

    // get the distance between p0 and p1
    dx = x1 - x0
    dy = y1 - y0
    d = sqrt(dx * dx + dy * dy)

    // find the two key angles and subtract them
    p1_p0_pen = acos((a1 * a1 - a0 * a0 - d * d) / (-2 * a0 * d))
    p0toP1 = atan2(y1 - y0, x1 - x0)
    angle = p0toP1 - p1_p0_pen

    // find the pen point
    pX = x0 + cos(angle) * a0
    pY = y0 + sin(angle) * a0
  
    lineTo(pX, pY)
    t += 0.01
  }
  stroke()
}

I’ve left lots of room for you to optimize this, so go for it! Even in its rough form though, this code should draw something like this:

There’s lots of things to play around and experiment with here. But the truth is that this is a pretty simple pintograph, so mostly the shapes are going roughly look like what you see above. But you can build on these principles and make all kinds of more complex devices. This site has several demos to inspire you:

https://michaldudak.github.io/pintograph/demo/

Disks on disks, and rotary pintographs and 3-wheel pintographs, etc.

And if you’re into this stuff, just search “harmonograph”, “pintograph” and “drawing machines” on Youtube to get an endless supply of inspiration – either for coding or actually building. Some of the more interesting ones (in my mind) are the ones that draw the curves on a piece of paper that itself is slowing rotating on a turntable.

Summary

This wraps up our discussion of Lissajous curves and the simulation of physical drawing machines. At least for a while. Next time we’ll go back a bit more to basic, standard geometric curves.