SplineFrames

Generates reference frames along a 3D spline path for DNA construction.

mdna.spline.SplineFrames

Source code in mdna/spline.py

class SplineFrames:

    def __init__(self, control_points, degree=3, closed=False, up_vector=[0, 0, 1],frame_spacing=0.34, twist=True, bp_per_turn=10.5, frame_tolerance=0.5, verbose=False, num_points=1000, initial_frame=None, modified_ranges=[],n_bp=None,dLk=None):
        """Initialize a spline path with Bishop-frame coordinate systems.

        The spline is fitted to the control points, evenly-spaced frames are
        placed along the arc length, continuity is verified, and an optional
        twist is applied.

        Args:
            control_points (numpy.ndarray): 3-D control points, shape ``(n, 3)``.
            degree (int): B-spline degree (1–5).
            closed (bool): If True, close the spline into a loop.
            up_vector (list or numpy.ndarray): Reference up direction for
                initial frame orientation.
            frame_spacing (float): Target arc-length distance between frames
                (nm).  Defaults to 0.34 nm (one base-pair rise).
            twist (bool): Apply a helical twist to the frames.
            bp_per_turn (float): Base pairs per helical turn (default 10.5).
            frame_tolerance (float): Maximum allowed angular deviation when
                testing frame continuity.
            verbose (bool): Print diagnostic messages.
            num_points (int): Number of spline evaluation points for internal
                interpolation.
            initial_frame (tuple or numpy.ndarray, optional): Explicit initial
                frame ``(origin, d1, d2, d3)`` to align the first frame.
            modified_ranges (list): List of ``(start, end, twist_angle)`` tuples
                specifying local twist modifications.
            n_bp (int, optional): Desired number of base pairs.  When set, the
                spline is rescaled to match this count.
            dLk (int, optional): Change in linking number applied via uniform
                twist redistribution.
        """
        self.control_points = control_points
        self.n = num_points
        self.twist = twist
        self.degree = degree
        self.closed = closed
        self.up_vector = np.array(up_vector)
        self.frame_spacing = frame_spacing
        self.bp_per_turn = bp_per_turn
        self.n_bp = n_bp
        self.dLk = dLk
        self.frame_tolerance = frame_tolerance
        self.initial_frame = initial_frame  # Added variable for the initial frame if we need to align the first frame with a given frame
        self.tck = None
        self.curve = None
        self.der1 = None
        self.der2 = None
        self.arc_length = None
        self.point_indices = None
        self.frames = []
        self.verbose = verbose
        self._initialize_spline()
        self._evaluate_spline()
        self.distribute_points()
        self.test_frames()

        if self.n_bp is not None:
            print(f"""\nStart rescaling spline based on requested number of base pairs.\n\tThis requires recomputation of the control points to match the desired number of base pairs.""")
            self._scale_to_nbp()
        else:
            self.n_bp = self.frames.shape[0]

        if self.twist:
            self.twist_frames(modified_ranges=modified_ranges)

    def update_initial_frame(self, initial_frame) -> 'SplineFrames':
        """
        Updates the initial frame.

        Args:
            initial_frame (tuple): Tuple containing the position, right, up, and forward vectors of the frame.

        Returns:
            SplineFrames: The updated instance.
        """
        self.initial_frame = initial_frame
        self._compute_frames()
        return self

    def update_control_points(self, control_points) -> None:
        """
        Updates the control points defining the path.

        Args:
            control_points (numpy.ndarray): Control points defining the path.
        """
        self.control_points = control_points
        self.frames = []
        self._initialize_spline()
        self._evaluate_spline()
        self.distribute_points()
        self.test_frames()
        # return self

    def update_spline_degree(self, degree) -> 'SplineFrames':
        """
        Updates the degree of the spline.

        Args:
            degree (int): Degree of the spline.

        Returns:
            SplineFrames: The updated instance.
        """
        self.degree = degree
        self._initialize_spline()
        self._evaluate_spline()
        return self

    def update_closed(self, closed) -> 'SplineFrames':
        """
        Updates the closed property indicating if the path is closed.

        Args:
            closed (bool): Indicates if the path is closed.

        Returns:
            SplineFrames: The updated instance.
        """
        self.closed = closed
        self._initialize_spline()
        self._evaluate_spline()
        return self

    def update_up_vector(self, up_vector) -> 'SplineFrames':
        """
        Updates the up vector used for frame computation. Not used at the moment.

        Args:
            up_vector (list or numpy.ndarray): Up vector for frame computation.

        Returns:
            SplineFrames: The updated instance.
        """
        self.up_vector = np.array(up_vector)
        self._compute_frames()
        return self

    def _initialize_spline(self) -> 'SplineFrames':
        """
        Initializes the spline with the control points.

        Returns:
            SplineFrames: The initialized instance.
        """
        if self.verbose:
            print("Initializing spline")
        control_points = self.control_points
        # if self.closed:
        #     control_points = np.vstack((self.control_points, self.control_points[0]))

        #     # Calculate tangent vectors at the first and last control points
        #     t_first = self.control_points[1] - self.control_points[0]
        #     t_last = self.control_points[-1] - self.control_points[-2]

        #     # Adjust the control points to ensure consistent tangent directions
        #     if np.dot(t_first, t_last) < 0:
        #         control_points[-1] = control_points[-1] * -1

        tck, _ = splprep(control_points.T, s=0, k=self.degree)#, per=self.closed)
        self.tck = tck
        return self

    def _evaluate_spline(self) -> 'SplineFrames':
        """
        Evaluates the spline curve and its derivatives.

        Returns:
            SplineFrames: The evaluated instance with ``curve``, ``der1``, and ``arc_length`` attributes set.
        """
        if self.verbose:
            print("Evaluating spline")
        dt = 1 / self.n
        self.v = np.arange(0, 1 + dt, dt)
        self.curve = np.array(splev(self.v, self.tck)).T
        self.der1 = np.array(splev(self.v, self.tck, der=1)).T
        if self.verbose:
            print("Calculating arc length")
        # Calculate the arc length of the spline
        arc_diffs = np.diff(self.curve, axis=0)

        arc = np.sqrt((arc_diffs**2).sum(axis=1))
        self.arc_length = np.hstack([0, cumulative_trapezoid(arc)])
        return self

    def distribute_points(self) -> 'SplineFrames':
        """
        Distributes equidistant points along the spline and evaluates the derivatives
        at these points. Adjusts the spacing to match the first and last points of the spline.

        Returns:
            SplineFrames: The instance with ``positions`` and ``derivatives`` attributes set.

        Raises:
            ValueError: If a suitable segment length cannot be found within the tolerance.
        """

        # Calculate the new segment length within the specified tolerance
        adjusted_arc_length = self.arc_length[-1] - self.frame_spacing  # Subtract the first and last segments
        num_segments = np.round(adjusted_arc_length / self.frame_spacing)
        new_segment_length = adjusted_arc_length / num_segments

        # If the new segment length is within the tolerance, proceed with the distribution
        if abs(new_segment_length - self.frame_spacing) <= self.frame_tolerance:
            self.segment_lengths = np.linspace(new_segment_length, adjusted_arc_length, int(num_segments))

            # Including the first and last points (and not the last segment if closed topology)
            self.t_values = [0] + list(np.interp(self.segment_lengths, self.arc_length, self.v[:len(self.arc_length)]))
            if not self.closed:
                self.t_values += [1] 

            if self.verbose:
                print(f"Evenly distributing {len(self.t_values)} points along the spline and computing derivatives.")

            # Evaluate the spline and its derivatives at the reparametrized t values
            equidistant_points = np.array(splev(self.t_values, self.tck)).T
            # Use difference vectors of equidistant points to calculate the derivatives
            derivatives = np.diff(equidistant_points, axis=0)
            #print(derivatives.shape,equidistant_points.shape)
            # Add the last derivative to the end
            derivatives = np.vstack([derivatives, derivatives[-1]])
            # Last point some numerical error that inverts the direction of the last point? Need to do more testing...
            #derivatives = np.array(splev(self.t_values, self.tck, der=1)).T
        else:
            raise ValueError(f"Cannot find a suitable segment length within the tolerance of {self.frame_tolerance} nm rise.")

        # Store or return the results
        self.positions = equidistant_points
        self.derivatives = derivatives
        # Assuming self._compute_frames() updates self with the new frames
        self._compute_frames()
        return self

    def rotation_matrix_from_vectors(self,vec1, vec2):
        """ Find the rotation matrix that aligns vec1 to vec2 """
        a, b = (vec1 / np.linalg.norm(vec1)).reshape(3), (vec2 / np.linalg.norm(vec2)).reshape(3)
        v = np.cross(a, b)
        c = np.dot(a, b)
        s = np.linalg.norm(v)
        kmat = np.array([[0, -v[2], v[1]], [v[2], 0, -v[0]], [-v[1], v[0], 0]])
        rotation_matrix = np.eye(3) + kmat + kmat.dot(kmat) * ((1 - c) / (s ** 2))
        return rotation_matrix

    def _compute_initial_frame(self):
        # Compute the initial frame based on the first derivative and up_vector
        T = self.derivatives[0] / np.linalg.norm(self.derivatives[0])
        N = np.cross(T, self.up_vector)
        if np.linalg.norm(N) < 1e-6: # Fallback if T and up_vector are parallel
            N = np.cross(T, np.array([1, 0, 0]) if abs(T[2]) < 0.9 else np.array([0, 1, 0]))
        N = N / np.linalg.norm(N)
        B = np.cross(T,N)
        #print('det initial frame',np.linalg.det(np.array([T,N,B])))
        self.frames.append((self.positions[0], N, B, T))
        #self.frames.append((self.positions[0], T, N, B))

    def _slide_frames(self):
        """
        Compute the coordinate frames along the curve by sliding the frame from the previous position to the current position.

        This method iterates through each position along the curve and computes the coordinate frame (position, T', N', B') at each position.
        The tangent vector (T') is computed as the normalized first derivative of the curve at the current position.
        The normal vector (N') is computed by rotating the previous normal vector (N) using the computed axis and angle.
        The binormal vector (B') is computed as the cross product of T' and N'.
        The computed coordinate frames are stored in the `frames` list for later use.

        Returns:
            None
        """

        # Iterate through each position along the curve, starting from the second position
        for i in range(1, len(self.positions)):
            # Get the last tangent vector (T) from the previously computed frames
            T_prev = self.frames[-1][3]  
            # Compute the new tangent vector (T') analytically as the normalized first derivative of the curve at the current position
            T = self.derivatives[i] / np.linalg.norm(self.derivatives[i])

            # Compute the axis of rotation as the normalized cross product of T and T', indicating the direction to rotate N to N'
            axis = np.cross(T_prev, T)
            # Check if the axis is significant to avoid division by zero and unnecessary rotation calculations
            if np.linalg.norm(axis) > 1e-6:
                axis = axis / np.linalg.norm(axis)  # Normalize the axis to ensure it has unit length
                # Compute the angle of rotation by the arccos of the dot product of T and T', clipped to [-1, 1] to avoid numerical issues
                angle = np.arccos(np.clip(np.dot(T_prev, T), -1.0, 1.0))
                # Rotate the normal vector (N) using the computed axis and angle to find the new normal vector (N')
                N = RigidBody.rotate_vector(self.frames[-1][1], axis, angle)
            else:
                # If the axis of rotation is negligible, use the previous normal vector (N) without rotation
                N = self.frames[-1][1]

            # Compute the binormal vector (B') as the cross product of T' and N', completing the coordinate frame
            B = np.cross(T, N)
            # Append the new coordinate frame (position, T', N', B') to the frames list for later use
            # print('det slide:',i,np.linalg.det(np.array([T,N,B])))
            self.frames.append((self.positions[i], N, B, T)) # original


    def _compute_frames(self):
        """
        Computes the coordinate frames on a spline with minimal torsion.
        Source: https://jbrd.github.io/2011/06/19/computing-coordinate-frames-on-a-spline-with-minimal-torsion.html

        """
        if self.initial_frame is not None:
            # Directly use the provided custom initial frame
            self.frames.append((self.positions[0], *self.initial_frame[1:]))
        else:
            # Compute the initial frame based on the first derivative and up_vector
            self._compute_initial_frame()

        # Compute frames for the rest of the points along the path
        self._slide_frames()

        # Convert the list of frames to a NumPy array for efficient storage and manipulation
        self.frames = np.array(self.frames)
        # # Swap the T and B vectors to match the expected order for the DNA generation
        # self.frames[:, [1, 3]] = self.frames[:, [3, 1]] 


    def twist_frames(self, modified_ranges=[], plot=False):
        self.twister = Twister(frames=self.frames, bp_per_turn=self.bp_per_turn, modified_ranges=modified_ranges, plot=plot, circular=self.closed,dLk=self.dLk)
        # This can be a separate call if you don't want to twist immediately upon calling twist_frames
        #self.twister._compute_and_plot_twists()

    def _scale_to_nbp(self):
        """Recompute control points to scale the spline to match the number of base pairs"""

        # Set the target number of base pairs, current number of base pairs, and total length of the spline
        if self.closed:
            target_bp = self.n_bp
        else:
            target_bp = self.n_bp - 1
        current_n_bp = self.frames.shape[0]
        total_length = self.arc_length[-1] 

        # Compute the ratio to scale the spline to match the target number of base pairs
        ratio = target_bp*self.frame_spacing/total_length
        new_control_points = self.control_points * ratio

        # Update the control points and reinitialize the spline
        self.update_control_points(new_control_points)  

        # Check if the number of base pairs matches the target
        if self.frames.shape[0] == self.n_bp:
            print(f"\tSpline scaled to match the target number of base pairs: {self.n_bp}\n")
        else:
            print(f"\tError: Spline could not be scaled to match the target number of base pairs: {target_bp}")
            print("\tNew number of base pairs:", self.frames.shape[0])


    def plot_frames(self, fig=False, equal_bounds=False, equal=True, spline=False, control_points=False, triads=True, transparent=False, legend=False) -> None:
        """
        Plots the frames along the spline.

        Note:
            This method needs to be called after the frames are computed.

        Args:
            fig (bool): If True, return the figure and axes objects instead of None.
            equal_bounds (bool): If True, set equal axis bounds based on data extent.
            equal (bool): If True, set equal aspect ratio on all axes.
            spline (bool): If True, overlay the continuous spline curve.
            control_points (bool): If True, plot the control points.
            triads (bool): If True, draw the right/up/forward triads at each frame.
            transparent (bool): If True, make the figure background transparent.
            legend (bool): If True, suppress the automatic legend.
        """
        fig = plt.figure()
        ax = fig.add_subplot(111, projection='3d')
        if triads:
            for frame in self.frames:
                position, right, up, forward = frame
                ax.quiver(*position, *right, length=0.2, color='g')
                ax.quiver(*position, *up, length=0.2, color='b')
                ax.quiver(*position, *forward, length=0.2, color='r')

        _ = self.tck[0]  # Spline parameters for plotting
        u = np.linspace(0, 1, 100)
        spline_points = np.array(splev(u, self.tck)).T
        if spline:
            ax.plot(*spline_points.T, color='gray', label='Spline')
        if equal_bounds:
            # Compute bounds
            all_points = np.vstack([spline_points, *[frame[0] for frame in self.frames]])
            max_bound = np.max(np.abs(all_points))

            ax.set_xlim([-max_bound, max_bound])
            ax.set_ylim([-max_bound, max_bound])
            ax.set_zlim([-max_bound, max_bound])
        if control_points:
            ax.scatter(*self.control_points.T, color='black', label='Control Points')
        if equal:
            ax.axis('equal')

        ax.set_xlabel('X')
        ax.set_ylabel('Y')
        ax.set_zlabel('Z')

        if spline or control_points:
            if not legend:
                ax.legend()



        if transparent:
            fig.patch.set_alpha(0)
            ax.patch.set_alpha(0)
            ax.grid(False)

        if fig:
            return fig, ax 
        else:
            return None 

    def test_frames(self) -> None:
        """
        Tests the computed frames for correctness.

        Checks consecutive frames for abrupt flips in the right and up vectors
        and prints warnings when angle deviations exceed 90 degrees.
        """
        for i in range(len(self.frames)-1):
            frame1 = self.frames[i]
            frame2 = self.frames[i+1]

            right1, up1 = frame1[1:3]
            right2, up2 = frame2[1:3]

            # Check if right vectors have flipped
            right_dot_product = np.dot(right1, right2)
            if right_dot_product < 0:
                angle_deviation = np.degrees(np.arccos(np.clip(np.dot(right1, right2), -1.0, 1.0)))
                print(f"Warning: Right Vectors may have flipped. Frame {i+1} to Frame {i+2}. Angle Deviation: {angle_deviation:.2f} degrees")

            # Check if up vectors have flipped
            up_dot_product = np.dot(up1, up2)
            if up_dot_product < 0:
                angle_deviation = np.degrees(np.arccos(np.clip(np.dot(up1, up2), -1.0, 1.0)))
                print(f"Warning: Up Vectors may have flipped. Frame {i+1} to Frame {i+2}. Angle Deviation: {angle_deviation:.2f} degrees")

        return None

__init__(control_points, degree=3, closed=False, up_vector=[0, 0, 1], frame_spacing=0.34, twist=True, bp_per_turn=10.5, frame_tolerance=0.5, verbose=False, num_points=1000, initial_frame=None, modified_ranges=[], n_bp=None, dLk=None)

Initialize a spline path with Bishop-frame coordinate systems.

The spline is fitted to the control points, evenly-spaced frames are placed along the arc length, continuity is verified, and an optional twist is applied.

Parameters:

Name	Type	Description	Default
control_points	ndarray	3-D control points, shape (n, 3).	required
degree	int	B-spline degree (1–5).	3
closed	bool	If True, close the spline into a loop.	False
up_vector	list or ndarray	Reference up direction for initial frame orientation.	[0, 0, 1]
frame_spacing	float	Target arc-length distance between frames (nm). Defaults to 0.34 nm (one base-pair rise).	0.34
twist	bool	Apply a helical twist to the frames.	True
bp_per_turn	float	Base pairs per helical turn (default 10.5).	10.5
frame_tolerance	float	Maximum allowed angular deviation when testing frame continuity.	0.5
verbose	bool	Print diagnostic messages.	False
num_points	int	Number of spline evaluation points for internal interpolation.	1000
initial_frame	tuple or ndarray	Explicit initial frame (origin, d1, d2, d3) to align the first frame.	None
modified_ranges	list	List of (start, end, twist_angle) tuples specifying local twist modifications.	[]
n_bp	int	Desired number of base pairs. When set, the spline is rescaled to match this count.	None
dLk	int	Change in linking number applied via uniform twist redistribution.	None

Source code in mdna/spline.py

def __init__(self, control_points, degree=3, closed=False, up_vector=[0, 0, 1],frame_spacing=0.34, twist=True, bp_per_turn=10.5, frame_tolerance=0.5, verbose=False, num_points=1000, initial_frame=None, modified_ranges=[],n_bp=None,dLk=None):
    """Initialize a spline path with Bishop-frame coordinate systems.

    The spline is fitted to the control points, evenly-spaced frames are
    placed along the arc length, continuity is verified, and an optional
    twist is applied.

    Args:
        control_points (numpy.ndarray): 3-D control points, shape ``(n, 3)``.
        degree (int): B-spline degree (1–5).
        closed (bool): If True, close the spline into a loop.
        up_vector (list or numpy.ndarray): Reference up direction for
            initial frame orientation.
        frame_spacing (float): Target arc-length distance between frames
            (nm).  Defaults to 0.34 nm (one base-pair rise).
        twist (bool): Apply a helical twist to the frames.
        bp_per_turn (float): Base pairs per helical turn (default 10.5).
        frame_tolerance (float): Maximum allowed angular deviation when
            testing frame continuity.
        verbose (bool): Print diagnostic messages.
        num_points (int): Number of spline evaluation points for internal
            interpolation.
        initial_frame (tuple or numpy.ndarray, optional): Explicit initial
            frame ``(origin, d1, d2, d3)`` to align the first frame.
        modified_ranges (list): List of ``(start, end, twist_angle)`` tuples
            specifying local twist modifications.
        n_bp (int, optional): Desired number of base pairs.  When set, the
            spline is rescaled to match this count.
        dLk (int, optional): Change in linking number applied via uniform
            twist redistribution.
    """
    self.control_points = control_points
    self.n = num_points
    self.twist = twist
    self.degree = degree
    self.closed = closed
    self.up_vector = np.array(up_vector)
    self.frame_spacing = frame_spacing
    self.bp_per_turn = bp_per_turn
    self.n_bp = n_bp
    self.dLk = dLk
    self.frame_tolerance = frame_tolerance
    self.initial_frame = initial_frame  # Added variable for the initial frame if we need to align the first frame with a given frame
    self.tck = None
    self.curve = None
    self.der1 = None
    self.der2 = None
    self.arc_length = None
    self.point_indices = None
    self.frames = []
    self.verbose = verbose
    self._initialize_spline()
    self._evaluate_spline()
    self.distribute_points()
    self.test_frames()

    if self.n_bp is not None:
        print(f"""\nStart rescaling spline based on requested number of base pairs.\n\tThis requires recomputation of the control points to match the desired number of base pairs.""")
        self._scale_to_nbp()
    else:
        self.n_bp = self.frames.shape[0]

    if self.twist:
        self.twist_frames(modified_ranges=modified_ranges)

distribute_points()

Distributes equidistant points along the spline and evaluates the derivatives at these points. Adjusts the spacing to match the first and last points of the spline.

Returns:

Name	Type	Description
SplineFrames	SplineFrames	The instance with positions and derivatives attributes set.

Raises:

Type	Description
ValueError	If a suitable segment length cannot be found within the tolerance.

Source code in mdna/spline.py

def distribute_points(self) -> 'SplineFrames':
    """
    Distributes equidistant points along the spline and evaluates the derivatives
    at these points. Adjusts the spacing to match the first and last points of the spline.

    Returns:
        SplineFrames: The instance with ``positions`` and ``derivatives`` attributes set.

    Raises:
        ValueError: If a suitable segment length cannot be found within the tolerance.
    """

    # Calculate the new segment length within the specified tolerance
    adjusted_arc_length = self.arc_length[-1] - self.frame_spacing  # Subtract the first and last segments
    num_segments = np.round(adjusted_arc_length / self.frame_spacing)
    new_segment_length = adjusted_arc_length / num_segments

    # If the new segment length is within the tolerance, proceed with the distribution
    if abs(new_segment_length - self.frame_spacing) <= self.frame_tolerance:
        self.segment_lengths = np.linspace(new_segment_length, adjusted_arc_length, int(num_segments))

        # Including the first and last points (and not the last segment if closed topology)
        self.t_values = [0] + list(np.interp(self.segment_lengths, self.arc_length, self.v[:len(self.arc_length)]))
        if not self.closed:
            self.t_values += [1] 

        if self.verbose:
            print(f"Evenly distributing {len(self.t_values)} points along the spline and computing derivatives.")

        # Evaluate the spline and its derivatives at the reparametrized t values
        equidistant_points = np.array(splev(self.t_values, self.tck)).T
        # Use difference vectors of equidistant points to calculate the derivatives
        derivatives = np.diff(equidistant_points, axis=0)
        #print(derivatives.shape,equidistant_points.shape)
        # Add the last derivative to the end
        derivatives = np.vstack([derivatives, derivatives[-1]])
        # Last point some numerical error that inverts the direction of the last point? Need to do more testing...
        #derivatives = np.array(splev(self.t_values, self.tck, der=1)).T
    else:
        raise ValueError(f"Cannot find a suitable segment length within the tolerance of {self.frame_tolerance} nm rise.")

    # Store or return the results
    self.positions = equidistant_points
    self.derivatives = derivatives
    # Assuming self._compute_frames() updates self with the new frames
    self._compute_frames()
    return self

plot_frames(fig=False, equal_bounds=False, equal=True, spline=False, control_points=False, triads=True, transparent=False, legend=False)

Plots the frames along the spline.

Note

This method needs to be called after the frames are computed.

Parameters:

Name	Type	Description	Default
fig	bool	If True, return the figure and axes objects instead of None.	False
equal_bounds	bool	If True, set equal axis bounds based on data extent.	False
equal	bool	If True, set equal aspect ratio on all axes.	True
spline	bool	If True, overlay the continuous spline curve.	False
control_points	bool	If True, plot the control points.	False
triads	bool	If True, draw the right/up/forward triads at each frame.	True
transparent	bool	If True, make the figure background transparent.	False
legend	bool	If True, suppress the automatic legend.	False

Source code in mdna/spline.py

def plot_frames(self, fig=False, equal_bounds=False, equal=True, spline=False, control_points=False, triads=True, transparent=False, legend=False) -> None:
    """
    Plots the frames along the spline.

    Note:
        This method needs to be called after the frames are computed.

    Args:
        fig (bool): If True, return the figure and axes objects instead of None.
        equal_bounds (bool): If True, set equal axis bounds based on data extent.
        equal (bool): If True, set equal aspect ratio on all axes.
        spline (bool): If True, overlay the continuous spline curve.
        control_points (bool): If True, plot the control points.
        triads (bool): If True, draw the right/up/forward triads at each frame.
        transparent (bool): If True, make the figure background transparent.
        legend (bool): If True, suppress the automatic legend.
    """
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    if triads:
        for frame in self.frames:
            position, right, up, forward = frame
            ax.quiver(*position, *right, length=0.2, color='g')
            ax.quiver(*position, *up, length=0.2, color='b')
            ax.quiver(*position, *forward, length=0.2, color='r')

    _ = self.tck[0]  # Spline parameters for plotting
    u = np.linspace(0, 1, 100)
    spline_points = np.array(splev(u, self.tck)).T
    if spline:
        ax.plot(*spline_points.T, color='gray', label='Spline')
    if equal_bounds:
        # Compute bounds
        all_points = np.vstack([spline_points, *[frame[0] for frame in self.frames]])
        max_bound = np.max(np.abs(all_points))

        ax.set_xlim([-max_bound, max_bound])
        ax.set_ylim([-max_bound, max_bound])
        ax.set_zlim([-max_bound, max_bound])
    if control_points:
        ax.scatter(*self.control_points.T, color='black', label='Control Points')
    if equal:
        ax.axis('equal')

    ax.set_xlabel('X')
    ax.set_ylabel('Y')
    ax.set_zlabel('Z')

    if spline or control_points:
        if not legend:
            ax.legend()



    if transparent:
        fig.patch.set_alpha(0)
        ax.patch.set_alpha(0)
        ax.grid(False)

    if fig:
        return fig, ax 
    else:
        return None 

rotation_matrix_from_vectors(vec1, vec2)

Find the rotation matrix that aligns vec1 to vec2

Source code in mdna/spline.py

def rotation_matrix_from_vectors(self,vec1, vec2):
    """ Find the rotation matrix that aligns vec1 to vec2 """
    a, b = (vec1 / np.linalg.norm(vec1)).reshape(3), (vec2 / np.linalg.norm(vec2)).reshape(3)
    v = np.cross(a, b)
    c = np.dot(a, b)
    s = np.linalg.norm(v)
    kmat = np.array([[0, -v[2], v[1]], [v[2], 0, -v[0]], [-v[1], v[0], 0]])
    rotation_matrix = np.eye(3) + kmat + kmat.dot(kmat) * ((1 - c) / (s ** 2))
    return rotation_matrix

test_frames()

Tests the computed frames for correctness.

Checks consecutive frames for abrupt flips in the right and up vectors and prints warnings when angle deviations exceed 90 degrees.

Source code in mdna/spline.py

def test_frames(self) -> None:
    """
    Tests the computed frames for correctness.

    Checks consecutive frames for abrupt flips in the right and up vectors
    and prints warnings when angle deviations exceed 90 degrees.
    """
    for i in range(len(self.frames)-1):
        frame1 = self.frames[i]
        frame2 = self.frames[i+1]

        right1, up1 = frame1[1:3]
        right2, up2 = frame2[1:3]

        # Check if right vectors have flipped
        right_dot_product = np.dot(right1, right2)
        if right_dot_product < 0:
            angle_deviation = np.degrees(np.arccos(np.clip(np.dot(right1, right2), -1.0, 1.0)))
            print(f"Warning: Right Vectors may have flipped. Frame {i+1} to Frame {i+2}. Angle Deviation: {angle_deviation:.2f} degrees")

        # Check if up vectors have flipped
        up_dot_product = np.dot(up1, up2)
        if up_dot_product < 0:
            angle_deviation = np.degrees(np.arccos(np.clip(np.dot(up1, up2), -1.0, 1.0)))
            print(f"Warning: Up Vectors may have flipped. Frame {i+1} to Frame {i+2}. Angle Deviation: {angle_deviation:.2f} degrees")

    return None

update_closed(closed)

Updates the closed property indicating if the path is closed.

Parameters:

Name	Type	Description	Default
closed	bool	Indicates if the path is closed.	required

Returns:

Name	Type	Description
SplineFrames	SplineFrames	The updated instance.

Source code in mdna/spline.py

def update_closed(self, closed) -> 'SplineFrames':
    """
    Updates the closed property indicating if the path is closed.

    Args:
        closed (bool): Indicates if the path is closed.

    Returns:
        SplineFrames: The updated instance.
    """
    self.closed = closed
    self._initialize_spline()
    self._evaluate_spline()
    return self

update_control_points(control_points)

Updates the control points defining the path.

Parameters:

Name	Type	Description	Default
control_points	ndarray	Control points defining the path.	required

Source code in mdna/spline.py

def update_control_points(self, control_points) -> None:
    """
    Updates the control points defining the path.

    Args:
        control_points (numpy.ndarray): Control points defining the path.
    """
    self.control_points = control_points
    self.frames = []
    self._initialize_spline()
    self._evaluate_spline()
    self.distribute_points()
    self.test_frames()

update_initial_frame(initial_frame)

Updates the initial frame.

Parameters:

Name	Type	Description	Default
initial_frame	tuple	Tuple containing the position, right, up, and forward vectors of the frame.	required

Returns:

Name	Type	Description
SplineFrames	SplineFrames	The updated instance.

Source code in mdna/spline.py

def update_initial_frame(self, initial_frame) -> 'SplineFrames':
    """
    Updates the initial frame.

    Args:
        initial_frame (tuple): Tuple containing the position, right, up, and forward vectors of the frame.

    Returns:
        SplineFrames: The updated instance.
    """
    self.initial_frame = initial_frame
    self._compute_frames()
    return self

update_spline_degree(degree)

Updates the degree of the spline.

Parameters:

Name	Type	Description	Default
degree	int	Degree of the spline.	required

Returns:

Name	Type	Description
SplineFrames	SplineFrames	The updated instance.

Source code in mdna/spline.py

def update_spline_degree(self, degree) -> 'SplineFrames':
    """
    Updates the degree of the spline.

    Args:
        degree (int): Degree of the spline.

    Returns:
        SplineFrames: The updated instance.
    """
    self.degree = degree
    self._initialize_spline()
    self._evaluate_spline()
    return self

update_up_vector(up_vector)

Updates the up vector used for frame computation. Not used at the moment.

Parameters:

Name	Type	Description	Default
up_vector	list or ndarray	Up vector for frame computation.	required

Returns:

Name	Type	Description
SplineFrames	SplineFrames	The updated instance.

Source code in mdna/spline.py

def update_up_vector(self, up_vector) -> 'SplineFrames':
    """
    Updates the up vector used for frame computation. Not used at the moment.

    Args:
        up_vector (list or numpy.ndarray): Up vector for frame computation.

    Returns:
        SplineFrames: The updated instance.
    """
    self.up_vector = np.array(up_vector)
    self._compute_frames()
    return self