moebius.py

   1
   2 import numpy as np
   3
   4 tau = np.pi * 2
   5
   6 origin =  (0,0,0)
   7 unit_x = (1,0,0)
   8 unit_y = (0,1,0)
   9 unit_z = (0,0,1)
  10
  11 def ParametricCircle:
  12   def __init__(pc, c, r0, r1):
  13     ''' circle centred on c
  14         with theta=0 point at c+r0
  15         and with theta=tau/4 point at c+r1 '''
  16     pc._c  = c
  17     pc._r0 = r0
  18     pc._r1 = r1
  19   def radius(pc, theta):
  20     return pc._r0 * cos(theta) + pc._r1 * sin(theta)
  21   def point(pc, theta):
  22     return pc._c + pc.radius(theta)
  23
  24 def Twirler(ParametricCircle):
  25   def __init__(tw, c, r0, r1, cycles, begin_zeta):
  26     ''' circle centred on c, etc.
  27         but with an orientation at each point, orthogonal to
  28           the circle
  29         the orientation circles round cycles times during the
  30           whole cycle
  31         begin_zeta is the angle from outwards at theta==0
  32           positive meaning in the direction of r0 x r1
  33     '''
  34     ParametricCircle.__init__(tw, c, r0, r1)
  35     tw._cycles = cycles
  36     tw._begin_zeta = begin_zeta
  37     tw._axis = np.cross(r0, r1)
  38   def dirn(tw, theta):
  39     zeta = tw._begin_zeta + theta * tw._cycles
  40     r = radius(tw, theta)
  41     return cos(zeta) * r + sin(zeta) * pc._axis
  42
  43 class Moebius:
  44   def __init__(m, n_u):
  45     m._edge    = Twirler(origin,  unit_z, unit_x, 2, 0)
  46     m._midline = Twirler(-unit_z, unit_z, unit_y, 1, 0)
  47     m._beziers = [ self._bezier(u) for u in np.linspace(0, 1, n_u) ]
  48   def _bezier(u):
  49     theta = u * tau
  50     cp = [None] * 4
  51     cp[0] =               m._edge   .point(theta)
  52     cp[1] = cp[0] + 0.5 * m._edge   .dirn (theta)
  53     cp[3] =               m._midline.point(theta*2)
  54     cp[2] = cp[3] + 0.5 * m._midline.dirn (theta*2)
  55     return BezierSegmentcp)
  56