An object-oriented SVG torus

from shape import Shape
import numpy as np

class Torus(Shape):
    def __init__(self, c, a, ntheta=30, nphi=60):
        self.c, self.a = c, a
        # Number of quads per 2.pi radians around each angular coordinate.
        self.ntheta, self.nphi = ntheta, nphi

        # Create a mesh of points on the angular coordinates, theta and phi.
        theta = np.linspace(0, 2.*np.pi, self.ntheta)
        phi = np.linspace(0, 2.*np.pi, self.nphi)
        self.theta, self.phi = np.meshgrid(theta, phi)

    def get_xyz(self):
        """Get the cartesian coordinates of a point on the torus."""

        x = (self.c + self.a*np.cos(self.theta)) * np.cos(self.phi)
        y = (self.c + self.a*np.cos(self.theta)) * np.sin(self.phi)
        z = self.a * np.sin(self.theta)
        self.xyz = np.stack((x, y, z), axis=-1)
        self.y = self.xyz[:,:,1]

    def setup_torus(self):
        self.get_xyz()

    def get_quads(self):
        # Calculate the coordinates of the quads, and keep track of each quad's
        # mean position in the y-axis (depth-axis) in ypos.

        # To look inside the torus, we may not want to draw it all. Specify the
        # indexes of the phi coordinates of the quads to draw as j1, j2 here:
        j1, j2 = self.nphi//2, self.nphi
        j1, j2 = 0, self.nphi
        ni, nj = self.ntheta, j2-j1
        self.quads = np.empty((ni*nj,4,2))
        self.ypos = np.empty(ni*nj)
        for i in range(self.ntheta):
            for j in range(j1, j2):
                ip, jp = (i+1) % self.ntheta, (j+1) % self.nphi
                idx = i*nj+j-j1
                self.quads[idx] = ((self.bx[j,i], self.by[j,i]),
                             (self.bx[j,ip], self.by[j,ip]),
                             (self.bx[jp,ip], self.by[jp,ip]),
                             (self.bx[jp,i], self.by[jp,i]))
                ym = np.mean( (self.y[j,i], self.y[j,ip],
                               self.y[jp,ip], self.y[jp,i]) )
                self.ypos[idx] = ym
        # Get the indexes of the quads sorted by their distance from the camera
        self.idy = np.argsort(self.ypos)

import numpy as np

class Shape:
    def __init__(self):
        pass


    def distort(self, noise=0.1):
        self.xyz += noise * (np.random.random(self.xyz.shape)-0.5)
        self.y = self.xyz[:,:,1]

    def rotate_xyz(self, alpha, beta, gamma):
        """Rotate the point(s) P=xyz by Rz(alpha)Ry(beta)Rx(gamma).P"""

        alpha, beta, gamma = (np.radians(alpha), np.radians(beta),
                              np.radians(gamma))
        ca, sa = np.cos(alpha), np.sin(alpha)
        Rz = np.array(((1,0,0),(0,ca,-sa),(0,sa,ca)))
        cb, sb = np.cos(beta), np.sin(beta)
        Ry = np.array(((cb,0,sb),(0,1,0),(-sb,0,cb)))
        cg, sg = np.cos(gamma), np.sin(gamma)
        Rx = np.array(((cg,-sg,0),(sg,cg,0),(0,0,1)))
        R = Rz @ Ry @ Rx
        self.xyz = self.xyz @ R.T
        self.y = self.xyz[:,:,1]

    def get_perspective_view(self, C, E):
        """
        Simple perspective projection of point A onto display plane at E
        as seen by camera at C. Returns bx,by the image plane coordinates
        of the projection (which are the x,z coordinates of B).

        """

        D = self.xyz - C
        r = E[1] / D[:,:,1]
        B = r[:,:,None] * D + E
        self.bx, self.by = B[:,:,0], B[:,:,2]

    # Some different colour schemes for our torus.
import random

def rgb_to_html(rgb):
    return '#{:02x}{:02x}{:02x}'.format(*rgb)

def html_to_frac(html):
    r, g, b = int(html[1:3], 16), int(html[3:5], 16), int(html[5:7], 16)
    return r/255, g/255, b/255

pastel_colours = ['#e0bbe4', '#957dad', '#d291bc', '#fec8d8', '#ffdfd3']

elmer_colours = [[68,90,233], [112,0,128], [103,204,33], [217,17,29],
         [228,151,200], [239,240,233], [23,16,16], [242,97,9], [251,236,21]]
for i, rgb in enumerate(elmer_colours):
    elmer_colours[i] = rgb_to_html(rgb)

def depth_shaded(i, idy):
    """Map the ramp i=0->b onto a triangle i=0->1->0."""
    b = len(idy)
    r = b - abs(i - b/2)
    rgb = (int(r/b*255),)*3
    return rgb_to_html(rgb)

def depth_gradient(i, idy, c1, c2):
    c1, c2 = html_to_frac(c1), html_to_frac(c2)
    d = c2[0]-c1[0], c2[1]-c1[1], c2[2]-c1[2]
    b = len(idy)
    r = 1 - abs(i - b/2) / b
    rgb = (c1[0] + r*d[0], c1[1] + r*d[1], c1[2] + r*d[2])
    rgb = (int(s*255) for s in rgb)
    return rgb_to_html(rgb)

colour_funcs = {
    'grey': lambda i, _: '#dddddd',
    'white': lambda i, _: '#ffffff',
    'pastels': lambda i, _: random.choice(pastel_colours),
    'elmer': lambda i, _: random.choice(elmer_colours),
    'depth shaded': depth_shaded,
    'depth gradient': depth_gradient,
}

import numpy as np
from torus import Torus
from palettes import colour_funcs

# The major and minor radius of the torus.
c, a = 2.5, 0.5
# Pick a colour scheme from the get_colours dictionary
palette = 'depth gradient'
palette_args = '#0000ff', '#00ff00'
# Image dimensions and scaling factors from torus units to image units.
width, height = 800, 600
scalex, scaley = 140, 140
# Tait-Bryan angles for intrinsic rotation of the figure.
alpha, beta, gamma = 90, 25, 10
# Camera position, C, and projection plane position, E (relative to C).
cx, cy, cz = 0, 6, 0
ex, ey, ez = 0, 3, 0
C, E = np.array((cx,cy,cz)), np.array((ex,ey,ez))
# Distortion noise factor
noise = 0.1


def preamble(fo):
    """The SVG preamble and styles."""

    print('<?xml version="1.0" encoding="utf-8"?>\n'

    '<svg xmlns="http://www.w3.org/2000/svg"\n' + ' '*5 +
       'xmlns:xlink="http://www.w3.org/1999/xlink" width="{}" height="{}" >'
            .format(width, height), file=fo)

    print("""
        <defs>
        <style type="text/css"><![CDATA[""", file=fo)

    print('path {stroke-width: 0.5px; stroke: #000;}', file=fo)

    print("""]]></style>
    </defs>""", file=fo)

torus = Torus(c, a)
torus.setup_torus()
torus.rotate_xyz(alpha, beta, gamma)
torus.distort(noise)

torus.get_perspective_view(C, E)
torus.get_quads()


# Draw the torus as a SVG image.
get_colour = colour_funcs[palette]
with open('torus.svg', 'w') as fo:
    preamble(fo)
    for i in torus.idy:
        quad = torus.quads[i] * (scalex, scaley) + (width/2, height/2)
        colour = get_colour(i, torus.idy, *palette_args)
        print('<path d="M{},{} L{},{} L{},{} L{},{} Z" fill="{}"/>'.format(
            *quad.ravel(), colour),
            file=fo)
    print('</svg>', file=fo)