On February 26th I was working on the draft of a future post, likely some time in May. I was writing a script that would allow for the plotting of most everything I had looked at/covered in this series of spirograph posts. That included the images covered in the earlier post this week. Which I had just begun working on in that future post.

That night I woke up with those kinds of plots on my mind. And, what I thought might be a better approach for the orthogonal line version popped into my head. The following morning I got it to work. So I thought, I’d post a quick update to this past Monday’s post. The future code is/will be different than what will be shown here. That’s because in the future post I can use six different shapes for the rotating wheels of the spirograph.

Original Method/Algorithm

My first attempt went something like this. I estimated a reasonable cycle length from the number of plot points. Then I randomly selected a maximum line length (lots of playing around with those numbers) and assigned to a variable. Then I created a multiplier value based on the maximum line length and the cycle size. And a variable indicating when I was halfway through the cycle so I would start reducing line widths

  # length of a cycle, need an integer, min 20 so doesn't get too small
  lc_frq = max(20, t_pts // 20)
  # sort line lengths and multipiler for some of the variations
  max_ln_len = np.random.random()*(.2 - .175) + .175
  lc_mult = 2 * max_ln_len / lc_frq
  # want to 1 index cycle loops
  lc_diff = lc_frq + 1
  # cycle halfway mark, go other way (size)
  lc_dir = int((lc_frq / 2) + 1)

Then I generated complex points for the curve. And, we are only interested in the final curve points. Since I am getting the slope for the current curve point based on the ones either side of it, we have a couple of variables pointing at those points. Then we loop over the all the data points. Getting the cartesian coordinates for the current point.

  r_pts = splt.f(t)
  c_fi = r_pts[-1]
  prv_pt = (np.real(c_fi[-1]), np.imag(c_fi[-1]))
  nxt_pt = (np.real(c_fi[1]), np.imag(c_fi[1]))
  
  for i in range(t_pts):
    c_x = np.real(c_fi[i])
    c_y = np.imag(c_fi[i])

Then for the case where we are plotting the lines simulating the width orthogonal to the curve, we get the estimated tangent using the two points either side of the current point. Assuming, of course, we are not going to divide by zero. We get the slope orthogonal to that tangent, and calculate the intercept. We use those values to create a function for the line in question. And update the next and previous point variables.

    elif use_axis == 's':
      is_0_div = (nxt_pt[0] - prv_pt[0]) == 0
      if not is_0_div:
        est_tngt = (nxt_pt[1] - prv_pt[1]) / (nxt_pt[0] - prv_pt[0])
        m1 = -1 / est_tngt
        b1 = c_y - (m1 * c_x)
        pln = lambda x: m1*x + b1
      prv_pt = (c_x, c_y)
      nxt_ndx = (i + 2) % t_pts
      nxt_pt = (np.real(c_fi[nxt_ndx]), np.imag(c_fi[nxt_ndx]))

Then I use that multiplier to create to new x values slightly to either side of the current point. Use the line function to get new y values for the two new x values. But those values could be rather large depending on the slope. So there was some messing around to reduce the size if necessary. Whatever necessary meant in my mind.

      if not is_0_div:
        x_1, x_2 = (c_x - lc_mult*10, c_x + lc_mult*10)
        y_1, y_2 = (pln(x_1), pln(x_2))
        ln_len = len_ln(x_1, y_1, x_2, y_2)
        if ln_len > 2:
          cnt_rsz += 1
          max_rsz = max(ln_len, max_rsz)
          min_rsz = min(ln_len, min_rsz)
          while (not ln_len < 2) and (ln_len >= 2.25): 
            y_1 *= .95
            y_2 *= .95
            ln_len = len_ln(x_1, y_1, x_2, y_2)

Then we plotted the line and moved on until we covered all the points. Got some interesting images, but not exactly what I was aiming for.

Revised Method/Algorithm

Well, numerous nights followed by a return to the subject, a period of wakefulness that night and bingo. Well, at least I think so.

Again, we start by estimating and generating some values. But a little differently this time. The cycle length is estimated the same way. But this time, I generate the maximum line length as a fraction of the maximum x or y values in the plot. The divisor, 6, took a bit of playing around until it suited my tastes.

Then, I generated an array of widths, going from near zero to the maximum (~ halfway) and back down to near zero.

      r_pts = splt.f(t)
      t_xs = np.real(r_pts[-1])
      t_ys = np.imag(r_pts[-1])

      # length of a cycle, need an integer, min 20 so doesn't get too small
      lc_frq = max(20, t_pts // 20)
      hf_frq = lc_frq // 2

      # sort line lengths
      max_ln_len = max(max(t_xs), max(t_ys)) / 6
      w_cycle = [(max_ln_len * i) / hf_frq for i in range(1, hf_frq+1)] + [(max_ln_len * i) / hf_frq for i in range(hf_frq+1, 1, -1)]

      # track the points to use for estimating the tangent
      prv_pt = (t_xs[-1], t_ys[-1])
      nxt_pt = (t_xs[1], t_ys[1])

The thing I realized during that period of wakefulness is that the width I wanted to generate was essentially the diameter of a suitably sized circle centered on the point of interest. I just had to get the angle right. And since we were calculating the slope of a line othogonal to the curve’s tangent at that point, getting the angle was all to simple. From there I could find the value of x that generated the point on the circumference at that angle. Add/subtract that x value to my current x position, and use the line equation to get the two y values I was after. Just had to deal with a couple of edge cases. Well I thought a couple of edge cases. But, might be able to cut that down a little as well.

Oh, and one other possible issue: limitations of computer arithmetic for very large or very small values. So, I add a little buffer when checking for possible zero division, and zero or ininfinite slopes.

      for i in range(t_pts):
        c_x = t_xs[i]
        c_y = t_ys[i]

        elif use_axis == 'o':
          # let's give ourselves some room sorting vertical or horizontal lines
          t_rng = 0.001

          # let's avoid div by zero
          is_0_div = (abs(nxt_pt[0] - prv_pt[0]) >= 0 - t_rng) and (abs(nxt_pt[0] - prv_pt[0]) <= 0 + t_rng)
          if not is_0_div:
            est_tngt = (nxt_pt[1] - prv_pt[1]) / (nxt_pt[0] - prv_pt[0])
            m1 = -1 / est_tngt
            b1 = c_y - (m1 * c_x)
            pln = lambda x: m1*x + b1

          prv_pt = (c_x, c_y)
          nxt_ndx = (i + 2) % t_pts
          nxt_pt = (t_xs[nxt_ndx], t_ys[nxt_ndx])

          if not is_0_div:
            t_tht = math.atan(m1)
            a_tht = abs(t_tht)
            if (a_tht >= 0 - t_rng and a_tht <= 0 + t_rng):
              x_1, x_2 = c_x, c_x
              y_1, y_2 = (c_y-(w_cycle[i % lc_frq] / 2), c_y+(w_cycle[i % lc_frq] / 2))
            elif (a_tht >= (np.pi / 2) - t_rng and a_tht <= (np.pi / 2) + t_rng):
              x_1, x_2 = (c_x-(w_cycle[i % lc_frq] / 2), c_x+(w_cycle[i % lc_frq] / 2))
              y_1, y_2 = c_y, c_y
            else:
              x_adj = math.cos(t_tht) * (w_cycle[i % lc_frq] / 2)
              x_1, x_2 = c_x - x_adj, c_x + x_adj
              y_1, y_2 = pln(x_1), pln(x_2)
          else:
            x_1, x_2 = c_x, c_x
            y_1, y_2 = (c_y-(w_cycle[i % lc_frq] / 2), c_y+(w_cycle[i % lc_frq] / 2))

And there you have it.

Examples Old and New

I am plotting the actual curve over the simulations to show that they are the same.

(ani-3.10) PS R:\learn\py_play\spirograph> python spiro_ortho.test.py
wheels: 5, k_f: 4, cgv: 2, points: 800
widths: [1, 0.7430092754274557, 0.652010706113512j, 0.4571170138062064j, 0.3446664837085567j],
heights: [1, 0.6153207444987907j, 0.4626712020626752j, 0.4286033164020989, 0.26958187196577493],
freqs: [-2, 10, -2, 18, -6]
ln wd: 4, ln alpha: 0.61
colour map: plt.cm.GnBu_r (16), max line height: 0.19353894113427933 (0.009676947056713967), line width: 4
freq: 40, max line height: 0.33179543750367846, line width: 4, w_cycle: 40
w_cycle: [0.01658977187518392, 0.03317954375036784, 0.04976931562555177, 0.06635908750073569, 0.08294885937591961, 0.09953863125110354, 0.11612840312628747, 0.13271817500147137, 0.1493079468766553, 0.16589771875183923, 0.18248749062702316, 0.1990772625022071, 0.215667034377391, 0.23225680625257494, 0.24884657812775884, 0.26543635000294274, 0.2820261218781267, 0.2986158937533106, 0.31520566562849456, 0.33179543750367846, 0.34838520937886236, 0.33179543750367846, 0.31520566562849456, 0.2986158937533106, 0.2820261218781267, 0.26543635000294274, 0.24884657812775884, 0.23225680625257494, 0.215667034377391, 0.1990772625022071, 0.18248749062702316, 0.16589771875183923, 0.1493079468766553, 0.13271817500147137, 0.11612840312628747, 0.09953863125110354, 0.08294885937591961, 0.06635908750073569, 0.04976931562555177, 0.03317954375036784]

The old way:

old style attempt at simulating cycling line with using lines orthogonal to the curve

The new way:

new style attempt at simulating cycling line with using lines orthogonal to the curve

Sorry about the differing plot sizes and/or aspect ratios between the two plot types.

You know, I kind of like the look of the old way as much as that of the new way. Though the new algorithm does what I was intending to do all along. Guess I will be modifying the code for that future post too include both approaches.

And another example.

(ani-3.10) PS R:\learn\py_play\spirograph> python spiro_ortho.test.py
wheels: 7, k_f: 4, cgv: 3, points: 800
widths: [1, 0.636277876552616, 0.5187067743890382j, 0.5151675479498962, 0.4649507042465718, 0.2761619129613484j, 0.14350939344935063j],
heights: [1, 0.6812722314912222, 0.4997507716378505, 0.33985517227592205j, 0.22171889908638162j, 0.16896385568419967, 0.12744155941437194],
freqs: [3, 7, -5, 11, -5, -5, -5]
ln wd: 7, ln alpha: 0.68
colour map: plt.cm.GnBu_r (16), max line height: 0.18593507524266892 (0.009296753762133446), line width: 7
freq: 40, max line height: 0.45206650022706163, line width: 7, w_cycle: 40
w_cycle: [0.022603325011353082, 0.045206650022706164, 0.06780997503405925, 0.09041330004541233, 0.11301662505676542, 0.1356199500681185, 0.15822327507947156, 0.18082660009082466, 0.20342992510217775, 0.22603325011353084, 0.24863657512488388, 0.271239900136237, 0.2938432251475901, 0.31644655015894313, 0.3390498751702962, 0.3616532001816493, 0.3842565251930024, 0.4068598502043555, 0.4294631752157086, 0.4520665002270617, 0.47466982523841467, 0.4520665002270617, 0.4294631752157086, 0.4068598502043555, 0.3842565251930024, 0.3616532001816493, 0.3390498751702962, 0.31644655015894313, 0.2938432251475901, 0.271239900136237, 0.24863657512488388, 0.22603325011353084, 0.20342992510217775, 0.18082660009082466, 0.15822327507947156, 0.1356199500681185, 0.11301662505676542, 0.09041330004541233, 0.06780997503405925, 0.045206650022706164]
old style attempt at simulating cycling line with using lines orthogonal to the curve new style attempt at simulating cycling line with using lines orthogonal to the curve

Done

Nice short post. Hope you found it of interest.

Should you be interested the code to generate the above images is available at Code Refactor Cycling Line Widths.

Well! Not Done

While working on another spirograph module (and perhaps a draft post) it occured to me that I was using a cycle of 16 colours for a line length cycle of 40 changes. Though it would look a lot better if I used that n_vals= parameter in set_colour_map() to get the same number of colours as lines in each cycle. Simple enough.

# rcm, cycle = splt.set_colour_map(ax)
rcm, cycle = splt.set_colour_map(ax, n_vals=lc_frq)

I have modified my code so that I can specify the colour map and have hard coded the curve parameters to match the ones for the curves above. Here’s the result of that simple change above.

old style attempt at simulating cycling line with using lines orthogonal to the curve new style attempt at simulating cycling line with using lines orthogonal to the curve

You can make out fairly clearly on the last image that the colours cycle in common with the line widths. But, can’t figure out why I am getting that line break in the bottom right of the curve. Thought I had dealt with the edge cases?

Tighter Gradients

And, looking at that I wondered how things would look if the colour gradient was just a little less extreme. I have a solution, but it ain’t pretty. And, took a bit more work than I was expecting.

I am going to create a colour cycle larger than needed. Then I am going to select a portion of that larger cycle to get the colour cycle to match the length of the changing line cycle. That should keep the colours closer together. I am selecting a larger set of colours in multiples of 50 based on the specified number desired. Figured that would be more than enough. Still making sure the second half of the colour cycle is the reverse of the first. Oh yes, I also decided not to have multiple lines with hard coded colour map instances. And, have made the selection of the set of colours somewhat random.

So, the new set_colour_map looks like the following. Bit long!

def set_colour_map(ax, n_vals=None, cc=None):
  # version 0.3.1
  clr_opt = ['default', 'plt.cm.GnBu_r', 'plt.cm.PuBu_r', 'plt.cm.viridis', 'plt.cm.hot', 'plt.cm.twilight_shifted']
  clr_cmap = ['default', plt.cm.GnBu_r, plt.cm.PuBu_r, plt.cm.viridis, plt.cm.hot, plt.cm.twilight_shifted]
  c_cycle = []
  u_cmap = None

  if cc and cc in clr_opt:
    u_cmap = clr_opt.index(cc)
  elif n_vals and n_vals < 20:  
    u_cmap = np.random.randint(0,6)
  elif n_vals and n_vals >= 20:
    u_cmap = np.random.randint(1,5)
  else:
    u_cmap = np.random.randint(0,6)

  if not n_vals or n_vals < 20:
    c_rng = 8
    if u_cmap >= 1:
      c_cycle = [clr_cmap[u_cmap](i) for i in np.linspace(0, 1, c_rng)] + [clr_cmap[u_cmap](i) for i in np.linspace(0, 1, c_rng+2)][1:-1]
  else:
    c_rng = n_vals // 2
    t_r = 50 * ((n_vals // 50) + 1)
    h_t_r = t_r // 2
    t_cycle = [clr_cmap[u_cmap](i) for i in np.linspace(0, 1, t_r)]
    s1, e1 = h_t_r - (c_rng // 2), h_t_r + (c_rng // 2)
    # s1, e1 = t_r - c_rng, t_r
    t_shft = np.random.randint(10, h_t_r + 1)
    s1, e1 = t_shft - (c_rng // 2), t_shft + (c_rng // 2)
    s2, e2 = e1 - 1, s1 - 1
    c_cycle = t_cycle[s1:e1]
    c_cycle = c_cycle + c_cycle[::-1]
    # print(f"colour cycle ({t_r}): {s1}:{e1}, {s2}:{e2}:-1 ({len(c_cycle)}, {c_rng})")

  # nothing to do if default colour map
  if u_cmap > 0:
    ax.set_prop_cycle('color', c_cycle)

  return clr_opt[u_cmap], c_cycle

And here’s how that works for a single random colour cycle. Won’t look quite like the ones above, but will still be comparable.

old style attempt at simulating cycling line with using lines orthogonal to the curve new style attempt at simulating cycling line with using lines orthogonal to the curve

I think there might be something to be said for the tighter colour gradient in the colour cycles for this kind of plot. Though may have to play with the random starting point being used to make the colour selection. The smaller it is, the darker the gradient. The bigger it is, the lighter. But given the use of various colour maps, there may not be a happy compromise.

Done, for sure!

Have an unlisted post with the updated code, should you be interested, Code Refactor Cycling Line Widths II.

I will have to look at perhaps making the gradient reduction optional or random. And, I expect I will need to update the function in spiro_plotlib. Or perhaps offer an alternate function and let the user decide. That is likely the simpler solution.

Refactor at Later Date

Have added the tight gradient function to the package spiro_plotlib. Reworked the test module for this post to use the package functions instead of the test function I had originally added to the module. Also added a slightly darker plot background for some of the colour maps.

Have added an unlisted post with the updated code, should you be interested, for the module and spiro_plotlib at Code Refactor Cycling Line Widths III.