catmul interpolation class

catmul_interpolation_3d.py

class CatmullInterpolation():

    def __init__(self, coord_list, num_segments, cyclic):
        self.cyclic = cyclic
        self.coord_list = coord_list
        self.num_segments = num_segments

        # coordlist must be decomposed into lists of just 
        # [x,x,x,x...], [y,y,y,y...] , [z,z,z,z...]

        self.coords_x = [v[0] for v in self.coord_list]
        self.coords_y = [v[1] for v in self.coord_list]
        self.coords_z = [v[2] for v in self.coord_list]


    def eval(self):
        vx = self.interpolate_catmul_1d(self.coords_x)
        vy = self.interpolate_catmul_1d(self.coords_y)
        vz = self.interpolate_catmul_1d(self.coords_z)
        return list(zip(vx, vy, vz))


    def get_sections(self, knots, cyclic):
        sections = []
        if not cyclic:
            if len(knots) == 2:
                sections.append((0, 0, 1, 1))
            else:
                for i in range(len(knots)-1):
                    if i == 0:
                        indices = i, i, i+1, i+2
                    elif i == len(knots)-2:
                        indices = i-1, i, i+1, i+1
                    else:
                        indices = i-1, i, i+1, i+2
                    sections.append(indices)
        else:
            if len(knots) == 2:
                sections.append((1, 0, 1, 0))
            else:
                last_idx = len(knots)-1
                for i in range(len(knots)):
                    if i == 0:
                        indices = last_idx, 0, 1, 2
                    elif i == last_idx:
                        indices = last_idx-1, last_idx, 0, 1
                    else:
                        indices = i-1, i, i+1, (i+2) % len(knots)
                    sections.append(indices)

        return sections


    def interpolate_catmul_1d(self, knots, cyclic, num_segments):
        """
        from https://www.mvps.org/directx/articles/catmull/
        q(t) = 0.5 *( (2 * P1) +
                    (-P0 + P2) * t +
                    (2*P0 - 5*P1 + 4*P2 - P3) * t2 +
                    (-P0 + 3*P1- 3*P2 + P3) * t3)
        """
        cyclic = self.cyclic
        num_segments = self.num_segments

        radii = []
        radii_add = radii.append

        sections = self.get_sections(knots, cyclic)
        theta = 1 / num_segments

        for idx, (p0, p1, p2, p3) in enumerate(sections):

            appendix = -1 if (cyclic or (idx < (len(sections) - 1))) else 0

            P0, P1, P2, P3 = knots[p0], knots[p1], knots[p2], knots[p3]
            for xt in range(num_segments + 1 + appendix):
                t = theta * xt
                t2 = t*t
                t3 = t2*t
                radt = 0.5*((2*P1)+(-P0+P2)*t+(2*P0-5*P1+4*P2-P3)*t2+(-P0+3*P1-3*P2+P3)*t3)
                radii_add(radt)

        return radii