Geometry Classes

class Obj(MeshData):
    """Load and save Wavefront OBJ geometry as CPU-side mesh data."""

    def _parse_vertex(self, tokens: list[str]) -> None:
        try:
            if len(tokens) not in (4, 7):
                raise ValueError
            self.vertex.append(
                Vec3(float(tokens[1]), float(tokens[2]), float(tokens[3]))
            )
            if len(tokens) == 7:
                self.colour.append(
                    Vec3(float(tokens[4]), float(tokens[5]), float(tokens[6]))
                )
        except (IndexError, ValueError) as error:
            raise ObjParseVertexError from error

    def _parse_normal(self, tokens: list[str]) -> None:
        try:
            self.normals.append(
                Vec3(float(tokens[1]), float(tokens[2]), float(tokens[3]))
            )
        except (IndexError, ValueError) as error:
            raise ObjParseNormalError from error

    def _parse_uv(self, tokens: list[str]) -> None:
        try:
            if len(tokens) not in (3, 4):
                raise ValueError
            self.uv.append(
                Vec3(
                    float(tokens[1]),
                    float(tokens[2]),
                    float(tokens[3]) if len(tokens) == 4 else 0.0,
                )
            )
        except (IndexError, ValueError) as error:
            raise ObjParseUVError from error

    @staticmethod
    def _index(value: str, length: int) -> int:
        index = int(value)
        if index == 0:
            raise ObjParseFaceError
        return length + index if index < 0 else index - 1

    def _parse_face(self, tokens: list[str]) -> None:
        if len(tokens) < 2:
            raise ObjParseFaceError
        parts = [token.split("/") for token in tokens[1:]]
        shapes = {len(part) for part in parts}
        if len(shapes) != 1 or shapes.pop() not in (1, 2, 3):
            raise ObjParseFaceError
        kind = len(parts[0])
        if any(
            not part[0] or (kind == 2 and not part[1]) or (kind == 3 and not part[2])
            for part in parts
        ):
            raise ObjParseFaceError
        if (
            kind == 3
            and any(part[1] == "" for part in parts)
            and not all(part[1] == "" for part in parts)
        ):
            raise ObjParseFaceError
        try:
            face = Face(
                vertex=[self._index(part[0], len(self.vertex)) for part in parts]
            )
            if kind == 2 or (kind == 3 and parts[0][1]):
                face.uv = [self._index(part[1], len(self.uv)) for part in parts]
            if kind == 3:
                face.normal = [
                    self._index(part[2], len(self.normals)) for part in parts
                ]
            self.faces.append(face)
        except (ValueError, IndexError, ObjParseFaceError) as error:
            raise ObjParseFaceError from error

    def load(self, file: str) -> bool:
        """Append geometry from an OBJ file and update its dimensions."""
        with open(file) as obj_file:
            for line in obj_file:
                tokens = line.strip().split()
                if not tokens:
                    continue
                if tokens[0] == "v":
                    self._parse_vertex(tokens)
                elif tokens[0] == "vn":
                    self._parse_normal(tokens)
                elif tokens[0] == "vt":
                    self._parse_uv(tokens)
                elif tokens[0] == "f":
                    self._parse_face(tokens)
        self.calc_dimensions()
        return True

    @classmethod
    def from_file(cls, fname: str) -> "Obj":
        """Create an Obj instance from a Wavefront file."""
        obj = cls()
        obj.load(fname)
        return obj

    def add_vertex(self, vertex: Vec3) -> None:
        """Append a vertex position."""
        self.vertex.append(vertex)

    def add_vertex_colour(self, vertex: Vec3, colour: Vec3) -> None:
        """Append a vertex position and its colour."""
        self.vertex.append(vertex)
        self.colour.append(colour)

    def add_normal(self, normal: Vec3) -> None:
        """Append a vertex normal."""
        self.normals.append(normal)

    def add_uv(self, uv: Vec3) -> None:
        """Append a texture coordinate."""
        self.uv.append(uv)

    def add_face(self, face: Face) -> None:
        """Append a face."""
        self.faces.append(face)

    def save(self, filename: str) -> None:
        """Save valid CPU mesh data as a Wavefront OBJ file."""
        self.validate()
        with open(filename, "w") as obj_file:
            obj_file.write("# This file was created by nccapy/Geo/Obj.py exporter\n")
            self._write_vertices(obj_file)
            self._write_uvs(obj_file)
            self._write_normals(obj_file)
            self._write_faces(obj_file)

    def _write_vertices(self, obj_file: TextIO) -> None:
        for index, vertex in enumerate(self.vertex):
            values = f"v {vertex.x} {vertex.y} {vertex.z}"
            if self.colour:
                colour = self.colour[index]
                values += f" {colour.x} {colour.y} {colour.z}"
            obj_file.write(f"{values}\n")

    def _write_uvs(self, obj_file: TextIO) -> None:
        for uv in self.uv:
            obj_file.write(f"vt {uv.x} {uv.y}\n")

    def _write_normals(self, obj_file: TextIO) -> None:
        for normal in self.normals:
            obj_file.write(f"vn {normal.x} {normal.y} {normal.z}\n")

    def _write_faces(self, obj_file: TextIO) -> None:
        for face in self.faces:
            corners = []
            for index, vertex in enumerate(face.vertex):
                corner = str(vertex + 1)
                if face.uv:
                    corner += f"/{face.uv[index] + 1}"
                if face.normal:
                    corner += (
                        f"/{face.normal[index] + 1}"
                        if face.uv
                        else f"//{face.normal[index] + 1}"
                    )
                corners.append(corner)
            obj_file.write(f"f {' '.join(corners)}\n")

class MeshData:
    """CPU-side mesh geometry shared by rendering back ends."""

    def __init__(self) -> None:
        """Create an empty CPU-side mesh."""
        self.vertex: list[Vec3] = []
        self.normals: list[Vec3] = []
        self.uv: list[Vec2 | Vec3] = []
        self.faces: list[Face] = []
        self.colour: list[Vec3] = []
        self.bbox: BBox | None = None
        self.min_x = self.max_x = 0.0
        self.min_y = self.max_y = 0.0
        self.min_z = self.max_z = 0.0

    def is_triangular(self) -> bool:
        """Return True when every face has exactly three vertices."""
        return all(len(face.vertex) == 3 for face in self.faces)

    def validate(self) -> None:
        """Validate indices and optional vertex attributes."""
        if self.colour and len(self.colour) != len(self.vertex):
            raise MeshValidationError("colour data must be empty or match vertex data")
        for face in self.faces:
            for name, values, limit in (
                ("vertex", face.vertex, len(self.vertex)),
                ("normal", face.normal, len(self.normals)),
                ("uv", face.uv, len(self.uv)),
            ):
                if name != "vertex" and not values:
                    continue
                if name == "vertex" and len(values) != 3:
                    raise MeshValidationError(
                        "faces must contain exactly three vertices"
                    )
                if name != "vertex" and len(values) != len(face.vertex):
                    raise MeshValidationError(f"face {name} indices must be complete")
                if any(index < 0 or index >= limit for index in values):
                    raise MeshValidationError(f"face {name} index is out of range")

    def calc_dimensions(self) -> None:
        """Calculate extents and bounding box from the vertex positions."""
        if not self.vertex:
            self.min_x = self.max_x = self.min_y = self.max_y = 0.0
            self.min_z = self.max_z = 0.0
            self.bbox = None
            return
        self.min_x = self.max_x = self.vertex[0].x
        self.min_y = self.max_y = self.vertex[0].y
        self.min_z = self.max_z = self.vertex[0].z
        for vertex in self.vertex[1:]:
            self.min_x, self.max_x = (
                min(self.min_x, vertex.x),
                max(self.max_x, vertex.x),
            )
            self.min_y, self.max_y = (
                min(self.min_y, vertex.y),
                max(self.max_y, vertex.y),
            )
            self.min_z, self.max_z = (
                min(self.min_z, vertex.z),
                max(self.max_z, vertex.z),
            )
        self.bbox = BBox.from_extents(
            self.min_x, self.max_x, self.min_y, self.max_y, self.min_z, self.max_z
        )

    def triangle_vertex_data(self, *, flip_v: bool = False) -> np.ndarray:
        """Pack expanded triangle data as position, normal and UV float32 values."""
        self.validate()
        if not self.faces:
            return np.empty(0, dtype=np.float32)
        data: list[float] = []
        for face in self.faces:
            for corner, vertex_index in enumerate(face.vertex):
                vertex = self.vertex[vertex_index]
                normal = self.normals[face.normal[corner]] if face.normal else Vec3()
                texcoord = self.uv[face.uv[corner]] if face.uv else Vec2()
                data.extend(
                    (
                        vertex.x,
                        vertex.y,
                        vertex.z,
                        normal.x,
                        normal.y,
                        normal.z,
                        texcoord.x,
                        1.0 - texcoord.y if flip_v else texcoord.y,
                    )
                )
        return np.ascontiguousarray(data, dtype=np.float32)

@dataclass(slots=True)
class Face:
    """Indices for one polygon face."""

    vertex: list[int] = field(default_factory=list)
    uv: list[int] = field(default_factory=list)
    normal: list[int] = field(default_factory=list)

class Primitives:
    """A static class for creating and drawing primitives."""

    # this is effectively a static class so we can use it to store data
    # and generate pipelines for drawing
    _primitives: Dict[str, _primitive] = {}
    _loaded: bool = False

    @classmethod
    def create(cls, type: str, name: str, *args: object, **kwargs: object) -> None:
        """Creates and stores a primitive object of the specified type.

        Prims.SPHERE : (radius: float, precision: int).
        Prims.TORUS : (radius: float, tube_radius: float, precision: int).
        Prims.LINE_GRID : (width: float, depth: float, steps: int).
        Prims.TRIANGLE_PLANE : ( width: float, depth: float, w_p: int, d_p: int, v_n: Vec3).
        Prims.CYLINDER : (radius: float, height: float, slices: int, stacks: int).
        Prims.CAPSULE : (radius: float, height: float, slices: int, stacks: int).
        Prims.CONE : (radius: float, height: float, slices: int, stacks: int).

        Args:
            type (str): The primitive type, typically from the Prims enum (e.g., Prims.SPHERE).
            name (str): The name to associate with the created primitive.
            *args: Positional arguments to pass to the primitive creation function (e.g., radius, precision).
            **kwargs: Keyword arguments to pass to the primitive creation function.

        Raises:
            ValueError: If the primitive type is not recognized.

        Example:
            Primitives.create(Prims.SPHERE, "sphere", 0.3, 32)
            Primitives.create(Prims.SPHERE, "sphere", radius=0.3, precision=32)
        """
        prim_methods = {
            Prims.SPHERE: PrimData.sphere,
            Prims.TORUS: PrimData.torus,
            Prims.LINE_GRID: PrimData.line_grid,
            Prims.TRIANGLE_PLANE: PrimData.triangle_plane,
            Prims.CYLINDER: PrimData.cylinder,
            Prims.DISK: PrimData.disk,
            Prims.CAPSULE: PrimData.capsule,
            Prims.CONE: PrimData.cone,
        }
        # line primitives are position-only GL_LINES data; everything else is
        # 8-float (position, normal, uv) triangle data
        prim_layouts = {
            Prims.LINE_GRID: (gl.GL_LINES, 3),
        }
        try:
            method = prim_methods[type]
        except KeyError:
            raise ValueError(f"Unknown primitive: {name}")

        draw_mode, floats_per_vertex = prim_layouts.get(type, (gl.GL_TRIANGLES, 8))
        cls._primitives[name] = _primitive(
            method(*args, **kwargs), draw_mode, floats_per_vertex
        )

    @classmethod
    def load_default_primitives(cls) -> None:
        """Loads the default primitives from the PrimData directory."""
        logger.info("Loading default primitives...")
        if not cls._loaded:
            for p in Prims:
                try:
                    prim_data = PrimData.primitive(p.value)
                    prim = _primitive(prim_data)
                    cls._primitives[p.value] = prim
                except Exception:
                    pass
            cls._loaded = True

    @classmethod
    def draw(cls, name: str | Prims) -> None:
        """Draws the specified primitive.

        Args:
            name: The name of the primitive to draw, either as a string or a Prims enum.
        """
        key = name.value if isinstance(name, Prims) else name
        try:
            prim = cls._primitives[key]
            with prim.vao:
                prim.vao.draw()
        except KeyError:
            logger.error(f"Failed to draw primitive {key}")
            return

class BaseMesh(MeshData):
    """MeshData with the legacy OpenGL drawing API."""

    def __init__(self) -> None:
        """Create an empty mesh with a legacy OpenGL renderer."""
        super().__init__()
        self.texture_id = 0
        self.texture = False
        self._renderer = OpenGLMesh(self)

    @property
    def vao(self) -> AbstractVAO | None:
        """Return the VAO owned by the compatibility renderer."""
        return self._renderer.vao

    def _should_skip_vao_creation(self, reset_vao: bool) -> bool:
        if self.vao is None:
            return False
        if reset_vao:
            logger.warning("VAO exists so returning")
            return True
        logger.warning("Creating new VAO")
        return False

    def _validate_triangular_mesh(self) -> None:
        if not self.is_triangular():
            raise RuntimeError("Can only create VBO from all Triangle data at present")

    def create_vao(self, reset_vao: bool = False) -> None:
        """Create the legacy VAO using the OpenGLMesh adapter."""
        if self._should_skip_vao_creation(reset_vao):
            return
        self._validate_triangular_mesh()
        self._renderer.texture_id = self.texture_id
        self._renderer.upload(force=self.vao is not None)

    def draw(self) -> None:
        """Draw if the VAO has been created, matching the old no-op behaviour."""
        if self.vao is not None:
            self._renderer.texture_id = self.texture_id
            self._renderer.draw()

class OpenGLMesh:
    """Upload and draw a MeshData instance with an OpenGL VAO."""

    def __init__(self, mesh: MeshData, *, texture_id: int = 0) -> None:
        """Store the CPU mesh and an optional borrowed texture ID."""
        self.mesh = mesh
        self.texture_id = texture_id
        self._vao: AbstractVAO | None = None

    @property
    def vao(self) -> AbstractVAO | None:
        """Return the VAO after upload."""
        return self._vao

    def upload(self, *, force: bool = False) -> None:
        """Create the VAO, replacing an existing one only when requested."""
        if self._vao is not None and not force:
            return
        if self._vao is not None:
            self.cleanup()
        data = self.mesh.triangle_vertex_data(flip_v=True)
        if not data.size:
            raise RuntimeError("cannot upload an empty mesh")
        self._vao = VAOFactory.create_vao(VAOType.SIMPLE, gl.GL_TRIANGLES)
        with self._vao as vao:
            count = data.size // 8
            vao.set_data(VertexData(data, count))
            vao.set_vertex_attribute_pointer(0, 3, gl.GL_FLOAT, 32, 0)
            vao.set_vertex_attribute_pointer(1, 3, gl.GL_FLOAT, 32, 12)
            vao.set_vertex_attribute_pointer(2, 2, gl.GL_FLOAT, 32, 24)
            vao.set_num_indices(count)
        self.mesh.calc_dimensions()

    def draw(self) -> None:
        """Draw the uploaded mesh."""
        if self._vao is None:
            raise RuntimeError("mesh must be uploaded before drawing")
        if self.texture_id:
            gl.glBindTexture(gl.GL_TEXTURE_2D, self.texture_id)
        with self._vao as vao:
            vao.draw()

    def cleanup(self) -> None:
        """Release the VAO whilst an OpenGL context is current."""
        if self._vao is not None:
            self._vao.remove_vao()
            self._vao = None

class PrimData:
    """Static methods generating packed vertex data for primitive shapes."""

    @staticmethod
    def line_grid(width: float, depth: float, steps: int) -> np.ndarray:
        """Creates a line grid primitive.

        Args:
            width: The width of the grid.
            depth: The depth of the grid.
            steps: The number of steps in the grid.
        """
        # Calculate the step size for each grid value
        wstep = width / steps
        ws2 = width / 2.0
        v1 = -ws2

        dstep = depth / steps
        ds2 = depth / 2.0
        v2 = -ds2

        # Create a list to store the vertex data
        data = []

        for _ in range(steps + 1):
            # Vertex 1 x, y, z
            data.append([-ws2, 0.0, v1])
            # Vertex 2 x, y, z
            data.append([ws2, 0.0, v1])

            # Vertex 1 x, y, z
            data.append([v2, 0.0, ds2])
            # Vertex 2 x, y, z
            data.append([v2, 0.0, -ds2])

            # Now change our step value
            v1 += wstep
            v2 += dstep

        # Convert the list to a NumPy array
        return np.array(data, dtype=np.float32)

    @staticmethod
    def triangle_plane(
        width: float, depth: float, w_p: int, d_p: int, v_n: Vec3
    ) -> np.ndarray:
        """Creates a triangle plane primitive.

        Args:
            width: The width of the plane.
            depth: The depth of the plane.
            w_p: The number of width partitions.
            d_p: The number of depth partitions.
            v_n: The normal vector for the plane.
        """
        w2 = width / 2.0
        d2 = depth / 2.0
        w_step = width / w_p
        d_step = depth / d_p

        du = 0.9 / w_p
        dv = 0.9 / d_p

        data = []
        v = 0.0
        d = -d2
        for _ in range(d_p):
            u = 0.0
            w = -w2
            for _ in range(w_p):
                # tri 1
                # vert 1
                data.extend([w, 0.0, d + d_step, v_n.x, v_n.y, v_n.z, u, v + dv])
                # vert 2
                data.extend(
                    [w + w_step, 0.0, d + d_step, v_n.x, v_n.y, v_n.z, u + du, v + dv]
                )
                # vert 3
                data.extend([w, 0.0, d, v_n.x, v_n.y, v_n.z, u, v])

                # tri 2
                # vert 1
                data.extend(
                    [w + w_step, 0.0, d + d_step, v_n.x, v_n.y, v_n.z, u + du, v + dv]
                )
                # vert 2
                data.extend([w + w_step, 0.0, d, v_n.x, v_n.y, v_n.z, u + du, v])
                # vert 3
                data.extend([w, 0.0, d, v_n.x, v_n.y, v_n.z, u, v])
                u += du
                w += w_step
            v += dv
            d += d_step

        return np.array(data, dtype=np.float32)

    @staticmethod
    def sphere(radius: float, precision: int) -> np.ndarray:
        """Creates a sphere primitive.

        Args:
            radius: The radius of the sphere.
            precision: The precision of the sphere (number of slices).
        """
        # Sphere code based on a function Written by Paul Bourke.
        # http://astronomy.swin.edu.au/~pbourke/opengl/sphere/
        # the next part of the code calculates the P,N,UV of the sphere for triangles

        # Disallow a negative number for radius.
        if radius < 0.0:
            radius = -radius

        # Disallow a negative number for precision.
        if precision < 4:
            precision = 4

        # Create a numpy array to store our verts
        data = []

        for i in range(precision // 2):
            theta1 = i * 2.0 * np.pi / precision - np.pi / 2.0
            theta2 = (i + 1) * 2.0 * np.pi / precision - np.pi / 2.0

            for j in range(precision):
                theta3 = j * 2.0 * np.pi / precision
                theta4 = (j + 1) * 2.0 * np.pi / precision

                # First triangle
                nx1 = np.cos(theta2) * np.cos(theta3)
                ny1 = np.sin(theta2)
                nz1 = np.cos(theta2) * np.sin(theta3)
                x1 = radius * nx1
                y1 = radius * ny1
                z1 = radius * nz1
                u1 = j / precision
                v1 = 2.0 * (i + 1) / precision
                data.append([x1, y1, z1, nx1, ny1, nz1, u1, v1])

                nx2 = np.cos(theta1) * np.cos(theta3)
                ny2 = np.sin(theta1)
                nz2 = np.cos(theta1) * np.sin(theta3)
                x2 = radius * nx2
                y2 = radius * ny2
                z2 = radius * nz2
                u2 = j / precision
                v2 = 2.0 * i / precision
                data.append([x2, y2, z2, nx2, ny2, nz2, u2, v2])

                nx3 = np.cos(theta1) * np.cos(theta4)
                ny3 = np.sin(theta1)
                nz3 = np.cos(theta1) * np.sin(theta4)
                x3 = radius * nx3
                y3 = radius * ny3
                z3 = radius * nz3
                u3 = (j + 1) / precision
                v3 = 2.0 * i / precision
                data.append([x3, y3, z3, nx3, ny3, nz3, u3, v3])

                # Second triangle
                nx4 = np.cos(theta2) * np.cos(theta4)
                ny4 = np.sin(theta2)
                nz4 = np.cos(theta2) * np.sin(theta4)
                x4 = radius * nx4
                y4 = radius * ny4
                z4 = radius * nz4
                u4 = (j + 1) / precision
                v4 = 2.0 * (i + 1) / precision
                data.append([x4, y4, z4, nx4, ny4, nz4, u4, v4])

                data.append([x1, y1, z1, nx1, ny1, nz1, u1, v1])
                data.append([x3, y3, z3, nx3, ny3, nz3, u3, v3])

        return np.array(data, dtype=np.float32)

    @staticmethod
    def cone(base: float, height: float, slices: int, stacks: int) -> np.ndarray:
        """Creates a cone primitive.

        Args:
            base: The radius of the cone's base.
            height: The height of the cone.
            slices: The number of divisions around the cone.
            stacks: The number of divisions along the cone's height.
        """
        z_step = height / (stacks if stacks > 0 else 1)
        r_step = base / (stacks if stacks > 0 else 1)

        cosn = height / np.sqrt(height * height + base * base)
        sinn = base / np.sqrt(height * height + base * base)

        cs = _circle_table(slices)

        z0 = 0.0
        z1 = z_step

        r0 = base
        r1 = r0 - r_step

        du = 1.0 / stacks
        dv = 1.0 / slices

        u = 1.0
        v = 1.0

        data = []

        for _ in range(stacks):
            for j in range(slices):
                # First triangle
                d1 = [0] * 8
                d1[6] = u
                d1[7] = v
                d1[3] = cs[j, 0] * cosn  # nx
                d1[4] = cs[j, 1] * sinn  # ny
                d1[5] = sinn  # nz
                d1[0] = cs[j, 0] * r0  # x
                d1[1] = cs[j, 1] * r0  # y
                d1[2] = z0  # z
                data.append(d1)

                d2 = [0] * 8
                d2[6] = u
                d2[7] = v - dv
                d2[3] = cs[j, 0] * cosn  # nx
                d2[4] = cs[j, 1] * sinn  # ny
                d2[5] = sinn  # nz
                d2[0] = cs[j, 0] * r1  # x
                d2[1] = cs[j, 1] * r1  # y
                d2[2] = z1  # z
                data.append(d2)

                d3 = [0] * 8
                d3[6] = u - du
                d3[7] = v - dv
                d3[3] = cs[j + 1, 0] * cosn  # nx
                d3[4] = cs[j + 1, 1] * sinn  # ny
                d3[5] = sinn  # nz
                d3[0] = cs[j + 1, 0] * r1  # x
                d3[1] = cs[j + 1, 1] * r1  # y
                d3[2] = z1  # z
                data.append(d3)

                # Second triangle
                d4 = [0] * 8
                d4[6] = u
                d4[7] = v
                d4[3] = cs[j, 0] * cosn  # nx
                d4[4] = cs[j, 1] * sinn  # ny
                d4[5] = sinn  # nz
                d4[0] = cs[j, 0] * r0  # x
                d4[1] = cs[j, 1] * r0  # y
                d4[2] = z0  # z
                data.append(d4)

                d5 = [0] * 8
                d5[6] = u - du
                d5[7] = v - dv
                d5[3] = cs[j + 1, 0] * cosn  # nx
                d5[4] = cs[j + 1, 1] * sinn  # ny
                d5[5] = sinn  # nz
                d5[0] = cs[j + 1, 0] * r1  # x
                d5[1] = cs[j + 1, 1] * r1  # y
                d5[2] = z1  # z
                data.append(d5)

                d6 = [0] * 8
                d6[6] = u - du
                d6[7] = v
                d6[3] = cs[j + 1, 0] * cosn  # nx
                d6[4] = cs[j + 1, 1] * sinn  # ny
                d6[5] = sinn  # nz
                d6[0] = cs[j + 1, 0] * r0  # x
                d6[1] = cs[j + 1, 1] * r0  # y
                d6[2] = z0  # z
                data.append(d6)

                u -= du

            v -= dv
            u = 1.0
            z0 = z1
            z1 += z_step
            r0 = r1
            r1 -= r_step

        return np.array(data, dtype=np.float32)

    @staticmethod
    def _add_cylinder_sides(
        data: list, radius: float, h: float, ang: float, precision: int
    ) -> None:
        """Generates cylinder side geometry."""
        for i in range(2 * precision):
            c = radius * np.cos(ang * i)
            c1 = radius * np.cos(ang * (i + 1))
            s = radius * np.sin(ang * i)
            s1 = radius * np.sin(ang * (i + 1))
            # normals for cylinder sides
            nc = np.cos(ang * i)
            ns = np.sin(ang * i)
            nc1 = np.cos(ang * (i + 1))
            ns1 = np.sin(ang * (i + 1))
            # side top
            data.extend([c1, h, s1, nc1, 0.0, ns1, 0.0, 0.0])
            data.extend([c, h, s, nc, 0.0, ns, 0.0, 0.0])
            data.extend([c, -h, s, nc, 0.0, ns, 0.0, 0.0])
            # side bot
            data.extend([c, -h, s, nc, 0.0, ns, 0.0, 0.0])
            data.extend([c1, -h, s1, nc1, 0.0, ns1, 0.0, 0.0])
            data.extend([c1, h, s1, nc1, 0.0, ns1, 0.0, 0.0])

    @staticmethod
    def _add_hemispherical_caps(
        data: list, radius: float, h: float, ang: float, precision: int
    ) -> None:
        """Generates hemispherical cap geometry."""
        for i in range(2 * precision):
            # longitude
            s = -np.sin(ang * i)
            s1 = -np.sin(ang * (i + 1))
            c = np.cos(ang * i)
            c1 = np.cos(ang * (i + 1))
            for j in range(precision + 1):
                o = h if j < precision / 2 else -h
                # latitude
                sb = radius * np.sin(ang * j)
                sb1 = radius * np.sin(ang * (j + 1))
                cb = radius * np.cos(ang * j)
                cb1 = radius * np.cos(ang * (j + 1))
                if j != precision - 1:
                    nx, ny, nz = sb * c, cb, sb * s
                    data.extend([nx, ny + o, nz, nx, ny, nz, 0.0, 0.0])
                    nx, ny, nz = sb1 * c, cb1, sb1 * s
                    data.extend([nx, ny + o, nz, nx, ny, nz, 0.0, 0.0])
                    nx, ny, nz = sb1 * c1, cb1, sb1 * s1
                    data.extend([nx, ny + o, nz, nx, ny, nz, 0.0, 0.0])
                if j != 0:
                    nx, ny, nz = sb * c, cb, sb * s
                    data.extend([nx, ny + o, nz, nx, ny, nz, 0.0, 0.0])
                    nx, ny, nz = sb1 * c1, cb1, sb1 * s1
                    data.extend([nx, ny + o, nz, nx, ny, nz, 0.0, 0.0])
                    nx, ny, nz = sb * c1, cb, sb * s1
                    data.extend([nx, ny + o, nz, nx, ny, nz, 0.0, 0.0])

    @staticmethod
    def capsule(radius: float, height: float, precision: int) -> np.ndarray:
        """Creates a capsule primitive.

        The capsule is aligned along the y-axis and is composed of a cylinder
        and two hemispherical caps. Based on code from
        https://code.google.com/p/rgine/source/browse/trunk/RGine/opengl/src/RGLShapes.cpp
        and adapted.
        """
        if radius <= 0.0:
            raise ValueError(RAD_POS)
        if height < 0.0:
            raise ValueError(NON_NEG)
        if precision < 4:
            precision = 4

        data = []
        h = height / 2.0
        ang = np.pi / precision

        # Cylinder sides
        PrimData._add_cylinder_sides(data, radius, h, ang, precision)

        # Hemispherical caps
        PrimData._add_hemispherical_caps(data, radius, h, ang, precision)

        return np.array(data, dtype=np.float32)

    @staticmethod
    def cylinder(radius: float, height: float, slices: int, stacks: int) -> np.ndarray:
        """Creates a cylinder primitive.

        The cylinder is aligned along the y-axis.
        This method generates the cylinder walls, but not the top and bottom caps.
        """
        if radius <= 0.0:
            raise ValueError(RAD_POS)
        if height < 0.0:
            raise ValueError(NON_NEG)
        if slices < 3:
            slices = 3
        if stacks < 1:
            stacks = 1

        data = []
        h2 = height / 2.0
        y_step = height / stacks

        cs = _circle_table(slices)

        du = 1.0 / slices
        dv = 1.0 / stacks

        for i in range(stacks):
            y0 = -h2 + i * y_step
            y1 = -h2 + (i + 1) * y_step
            v = i * dv
            for j in range(slices):
                u = j * du

                nx1, nz1 = cs[j, 0], cs[j, 1]
                x1, z1 = radius * nx1, radius * nz1

                nx2, nz2 = cs[j + 1, 0], cs[j + 1, 1]
                x2, z2 = radius * nx2, radius * nz2

                p_bl = [x1, y0, z1, nx1, 0, nz1, u, v]
                p_br = [x2, y0, z2, nx2, 0, nz2, u + du, v]
                p_tl = [x1, y1, z1, nx1, 0, nz1, u, v + dv]
                p_tr = [x2, y1, z2, nx2, 0, nz2, u + du, v + dv]

                # Triangle 1
                data.extend(p_bl)
                data.extend(p_tl)
                data.extend(p_br)
                # Triangle 2
                data.extend(p_br)
                data.extend(p_tl)
                data.extend(p_tr)

        return np.array(data, dtype=np.float32)

    @staticmethod
    def disk(radius: float, slices: int) -> np.ndarray:
        """Creates a disk primitive.

        Args:
            radius: The radius of the disk.
            slices: The number of slices to divide the disk into.
        """
        if radius <= 0.0:
            raise ValueError(RAD_POS)
        if slices < 3:
            slices = 3

        data = []
        cs = _circle_table(slices)

        center = [0, 0, 0, 0, 1, 0, 0.5, 0.5]

        for i in range(slices):
            p1 = [
                radius * cs[i, 0],
                0,
                radius * cs[i, 1],
                0,
                1,
                0,
                cs[i, 0] * 0.5 + 0.5,
                cs[i, 1] * 0.5 + 0.5,
            ]
            p2 = [
                radius * cs[i + 1, 0],
                0,
                radius * cs[i + 1, 1],
                0,
                1,
                0,
                cs[i + 1, 0] * 0.5 + 0.5,
                cs[i + 1, 1] * 0.5 + 0.5,
            ]

            data.extend(center)
            data.extend(p2)
            data.extend(p1)

        return np.array(data, dtype=np.float32)

    @staticmethod
    def torus(
        minor_radius: float,
        major_radius: float,
        sides: int,
        rings: int,
    ) -> np.ndarray:
        """Creates a torus primitive.

        Args:
            minor_radius: The minor radius of the torus.
            major_radius: The major radius of the torus.
            sides: The number of sides for each ring.
            rings: The number of rings for the torus.
        """
        if minor_radius <= 0 or major_radius <= 0:
            raise ValueError(RAD_POS)
        if sides < 3 or rings < 3:
            raise ValueError("Sides and rings must be at least 3")

        d_psi = 2.0 * np.pi / rings
        d_phi = -2.0 * np.pi / sides

        psi = 0.0

        vertices = []
        normals = []
        uvs = []

        for j in range(rings + 1):
            c_psi = np.cos(psi)
            s_psi = np.sin(psi)
            phi = 0.0
            for i in range(sides + 1):
                c_phi = np.cos(phi)
                s_phi = np.sin(phi)

                x = c_psi * (major_radius + c_phi * minor_radius)
                z = s_psi * (major_radius + c_phi * minor_radius)
                y = s_phi * minor_radius
                vertices.append([x, y, z])

                nx = c_psi * c_phi
                nz = s_psi * c_phi
                ny = s_phi
                normals.append([nx, ny, nz])

                u = i / sides
                v = j / rings
                uvs.append([u, v])

                phi += d_phi
            psi += d_psi

        data = []
        for j in range(rings):
            for i in range(sides):
                idx1 = j * (sides + 1) + i
                idx2 = j * (sides + 1) + (i + 1)
                idx3 = (j + 1) * (sides + 1) + i
                idx4 = (j + 1) * (sides + 1) + (i + 1)

                p1 = vertices[idx1] + normals[idx1] + uvs[idx1]
                p2 = vertices[idx2] + normals[idx2] + uvs[idx2]
                p3 = vertices[idx3] + normals[idx3] + uvs[idx3]
                p4 = vertices[idx4] + normals[idx4] + uvs[idx4]

                data.extend(p1)
                data.extend(p3)
                data.extend(p2)

                data.extend(p2)
                data.extend(p3)
                data.extend(p4)

        return np.array(data, dtype=np.float32)

    @staticmethod
    def primitive(name: str | Enum) -> np.ndarray:
        """Load pre-generated vertex data for the named primitive."""
        prim_folder = Path(__file__).parent / "PrimData"
        prims = np.load(prim_folder / "Primitives.npz")
        if isinstance(name, Prims):
            name = name.value

        try:
            return prims[name]
        except KeyError:
            raise ValueError(f"Primitive '{name}' not found")