1
$\begingroup$

This is a follow up question to my previous one regarding objects moving in a circle (which was a simple test case I had put together). I really appreciate the answer and detailed explanation provided by Blunder. I do understand the issue now, but I am still no closer to getting my more complex animation (involving Mobius strip paths) to work correctly.

Mobius animation

Here the object is rotating along all 3 axes and I am getting the Euler angles from a "normal" direction vector (joining the location of the object pnt to the center circle of the Mobius strip cnt). The angles sometimes change abruptly by 180 or 360 degrees which I assume is due to the looping and the nature of the Mobius strip itself. I am just not sure how these angles can be manipulated to avoid the weird unintended rotations that are happening as the objects move along the Mobius strip paths. I tried using quaternions as well but manipulating them is less intuitive than working with Eulers.

While my intent is to animate objects moving along 9 paths defined by mathematical equations, the code I have attached simplifies this by defining just 3 paths. I would greatly appreciate any insights I can get into addressing the weird rotations in this situation. I have been struggling with this for several weeks now. Thanks!

import numpy as np
import mathutils
from mathutils import Vector
import bpy
import sys

bpy.ops.object.select_all(action='SELECT')
bpy.ops.object.delete(use_global=True)

matr = bpy.data.materials.new("Red")
matr.diffuse_color = (1,0,0,1)
matg = bpy.data.materials.new("Green")
matg.diffuse_color = (0,1,0,1)
matb = bpy.data.materials.new("Blue")
matb.diffuse_color = (0,0,1,1)

columns = 42
rows = 42

theta = np.radians(np.linspace(0,360,columns+1)[:-1])
phi = np.radians(np.array([0,120,240]))
#phi = np.radians(np.array([0,30,90,120,150,210,240,270,330]))
#phi = np.radians(np.array([0,11,30,60,90,109,120,131,150,180,210,229,240,251,270,300,330,349]))
gamma = np.radians(np.linspace(0,360,rows+1)[:-1])

phi = phi[:,None]
gamma2 = gamma[:,None]
gamma = gamma[:,None,None]

n = 3
a = 1
t = 2

# Cross section is an equilateral triangle with circumradius r
r = np.cos(np.pi/n) / np.cos(np.mod(phi,2*np.pi/n) - np.pi/n)

# 1D array
x1 = (4/3) * r * np.cos(phi)
y1 = (4/3) * r * np.sin(phi)

# 2D array
x2 = (a + (x1 / 2) * np.cos(t * theta / 2) - (y1 / 2) * np.sin(t * theta / 2)) * np.cos(theta)
y2 = (a + (x1 / 2) * np.cos(t * theta / 2) - (y1 / 2) * np.sin(t * theta / 2)) * np.sin(theta)
z2 = (x1 / 2) * np.sin(t * theta / 2) + (y1 / 2) * np.cos(t * theta / 2)

# 3D array
x = x2 * np.cos(gamma) - y2 * np.sin(gamma)
y = x2 * np.sin(gamma) + y2 * np.cos(gamma)
z = np.tile(z2,(rows,1,1))

cx2 = a * np.cos(theta)
cy2 = a * np.sin(theta)
cz2 = np.tile(0,columns)

cx = cx2 * np.cos(gamma2) - cy2 * np.sin(gamma2)
cy = cx2 * np.sin(gamma2) + cy2 * np.cos(gamma2)
cz = np.tile(cz2,(rows,1))

# Loop over j for different paths
# Loop over k for different points on each path
# Create cones
for j in range(phi.size):
    for k in range(theta.size):
        bpy.ops.mesh.primitive_cone_add(vertices=4)
        cone = bpy.context.active_object
        if j % 3 == 0:
            cone.active_material = matr
        elif j % 3 == 1:
            cone.active_material = matg
        elif j % 3 == 2:
            cone.active_material = matb
        cone.name = 'cone_{}_{}'.format(j,k)
        cone.scale = (0.075, 0.075, 0.03)

frame_cnt = 336

key_frame_cnt = rows

frame_step = frame_cnt / key_frame_cnt

bpy.context.scene.frame_end = frame_cnt

frame_num = 0

for i in range(gamma.size):
    for j in range(phi.size):
        for k in range(theta.size):
            if k - i < 0:
                m = k - i + columns
            else:
                m = k - i

            cone = bpy.data.objects['cone_{}_{}'.format(j,k)]
            cone.location = (x[i][j][m],y[i][j][m],z[i][j][m])

            cnt = Vector((cx[i][m],cy[i][m],cz[i][m]))
            pnt = Vector((x[i][j][m],y[i][j][m],z[i][j][m]))

            rot = pnt - cnt
            rot.normalize()

            rot_quat = rot.to_track_quat('Z','Y')
            rot_euler = rot_quat.to_euler()

            rot_x = rot_euler[0]
            rot_y = rot_euler[1]
            rot_z = rot_euler[2]

            # Manipulate rot_x, rot_y and rot_z to resolve flips but how?

            cone.rotation_euler = (rot_x,rot_y,rot_z)

            cone.keyframe_insert(data_path="location", frame=frame_num)
            cone.keyframe_insert(data_path="rotation_euler", frame=frame_num)
    frame_num += frame_step

# Last key frame
for j in range(phi.size):
    for k in range(theta.size):
        cone = bpy.data.objects['cone_{}_{}'.format(j,k)]
        cone.location = (x[0][j][k],y[0][j][k],z[0][j][k])

        cnt = Vector((cx[0][k],cy[0][k],cz[0][k]))
        pnt = Vector((x[0][j][k],y[0][j][k],z[0][j][k]))

        rot = pnt - cnt

        rot.normalize()

        rot_quat = rot.to_track_quat('Z','Y')
        rot_euler = rot_quat.to_euler()

        rot_x = rot_euler[0]
        rot_y = rot_euler[1]
        rot_z = rot_euler[2]

        cone.rotation_euler = (rot_x,rot_y,rot_z)

        cone.keyframe_insert(data_path="location", frame=frame_num)
        cone.keyframe_insert(data_path="rotation_euler", frame=frame_num)

$\endgroup$
4
  • $\begingroup$ It looks like there are two problems that cause these "random" rotations: 1) the groups of cones that have strange interpolated rotations are caused by gimbal locks. You can prevent this if you use Quaternion instead of Euler. 2) At some point, the individual cones rotate 180° one after the other. This seems to be similar to blender.stackexchange.com/a/218657/107598 - you make the cones tracking a target postion with the to_track_quat function. The calculation works in the range of -180° to 180°. It's ok for poses but for interpolated animations the rotation direction is messed up ... $\endgroup$
    – Blunder
    Commented Mar 12 at 11:07
  • $\begingroup$ ... because the direction of rotation is ignored and you can achieve a rotary position on one axis by making several turns in different directions. 270° is the same as -90° (=270°-360°) or 630° (=270°+360°). For example, your rotation animation uses the poses 0°, 90°, 180°, 270°, 360° (rotating clockwise), but to_track_quat() calculates -90° (counterclockwise) instead of 270° (clockwise). The resulting key frames are then: 0°, 90°, 180°, -90, 0° instead of 0°, 90°, 180°, 270°, 360°. The individual poses/key frames are the same but the interpolated movement is not. $\endgroup$
    – Blunder
    Commented Mar 12 at 11:32
  • $\begingroup$ To fix this you will have to consider the rotation of the previous frame (=direction of rotation) and correct the calculated rotation if it is in the range of -180° to 0° instead of 180° to 360°. I.e. at -90° add the correction value of 360° (-90+360=270). The same for the other direction. -Or- do without interpolation and calculate a pose/key frame for each frame. $\endgroup$
    – Blunder
    Commented Mar 12 at 11:34
  • $\begingroup$ Thank you Blunder for once again coming to my rescue! To avoid gimbal locks I am assuming I should define the key frames using rotation_quaternion instead of rotation_euler. I had tried that and had the same weird behavior. I am assuming that is because of the second issue. But if I am using quaternions I can't manipulate individual angles. I am assuming what I need to do instead is negate the quaternion under certain conditions ? Can you please confirm ? Thanks again ! $\endgroup$
    – Deep Shen
    Commented Mar 12 at 13:14

1 Answer 1

1
$\begingroup$

The first problem is the XY axes of the cones rotate 90 degrees when Z crosses the global vertical axis. This is a real thing that happens at the keyframes, not an interpolation artifact, so it won't be changed whether you use eulers, quaternions, or whatever.

This happens because you used to_track_quat. Conceptually you cannot ever make a rotation by just knowing where the up (Z) axis is, because the XY plane is still free to rotate around that up axis. You also need a "forward" direction.

To put it another way, you need a consistent moving frame on the curve. You have your consistent normal vector, now you need a consistent tangent vector. (The bitangent is then fixed as their cross product, up to sign.)

I don't know how to compute the tangent for your curve, so I'll just approximate it with the vector pointing to the next cone on the path.

Once we have the local XYZ axes, ie the rotation matrix, we can easily compute the rotation that will give those axes.

for i in range(gamma.size):
    for j in range(phi.size):
        for k in range(theta.size):
            if k - i < 0:
                m = k - i + columns
            else:
                m = k - i

            next_k = 0 if k == theta.size - 1 else k+1
            next_m = next_k - i + columns if next_k - i < 0 else next_k - i

            cone = bpy.data.objects['cone_{}_{}'.format(j,k)]
            cone.location = (x[i][j][m],y[i][j][m],z[i][j][m])

            cnt = Vector((cx[i][m],cy[i][m],cz[i][m]))
            pnt = Vector((x[i][j][m],y[i][j][m],z[i][j][m]))
            next_pnt = Vector((x[i][j][next_m],y[i][j][next_m],z[i][j][next_m]))

            up = pnt - cnt
            up.normalize()

            forward = next_pnt - pnt
            forward.normalize()
            
            # Compute rotation matrix from local axes
            rot = mathutils.Matrix.Identity(3)
            rot[0] = forward.cross(up)
            rot[1] = forward
            rot[2] = up
            rot.transpose()  # set columns, not rows

            rot_euler = rot.to_euler()
            
            cone.rotation_euler = rot_euler

            cone.keyframe_insert(data_path="location", frame=frame_num)
            cone.keyframe_insert(data_path="rotation_euler", frame=frame_num)
    frame_num += frame_step

It is still not fixed. The second problem should be familiar though, it was what you asked about in your last question. Sometimes the cone will make unnecessary full rotations, or go backwards instead of forwards.

To fix this, you can use e2.make_compatible(e1) to make each euler, e2, consistent with the euler from the previous frame, e1. "Compatible" means it takes the shortest path to get to e2.

I'll store the euler from the previous frame in a dict indexed by (j,k) so we can access it on the next iteration.

prevs = {}  # NEW

for i in range(gamma.size):
    for j in range(phi.size):
        for k in range(theta.size):
            if k - i < 0:
                m = k - i + columns
            else:
                m = k - i

            next_k = 0 if k == theta.size - 1 else k+1
            next_m = next_k - i + columns if next_k - i < 0 else next_k - i

            cone = bpy.data.objects['cone_{}_{}'.format(j,k)]
            cone.location = (x[i][j][m],y[i][j][m],z[i][j][m])

            cnt = Vector((cx[i][m],cy[i][m],cz[i][m]))
            pnt = Vector((x[i][j][m],y[i][j][m],z[i][j][m]))
            next_pnt = Vector((x[i][j][next_m],y[i][j][next_m],z[i][j][next_m]))

            up = pnt - cnt
            up.normalize()

            forward = next_pnt - pnt
            forward.normalize()
            
            # Compute rotation matrix from local axes
            rot = mathutils.Matrix.Identity(3)
            rot[0] = forward.cross(up)
            rot[1] = forward
            rot[2] = up
            rot.transpose()  # set columns, not rows

            rot_euler = rot.to_euler()
            
            if i > 0:                                    # NEW
                rot_euler.make_compatible(prevs[(j,k)])  # NEW
            prevs[(j,k)] = rot_euler.copy()              # NEW
            
            cone.rotation_euler = rot_euler

            cone.keyframe_insert(data_path="location", frame=frame_num)
            cone.keyframe_insert(data_path="rotation_euler", frame=frame_num)
    frame_num += frame_step

Is it fixed now? No! You will notice the cones still "wobble" sometimes, not going backwards, but "twisting" on the way to their target instead of moving "straight".

This is characteristic gimbal lock. The solution is to use quaternions instead of Euler angles. That's easy.

prevs = {}

for i in range(gamma.size):
    for j in range(phi.size):
        for k in range(theta.size):
            if k - i < 0:
                m = k - i + columns
            else:
                m = k - i

            next_k = 0 if k == theta.size - 1 else k+1
            next_m = next_k - i + columns if next_k - i < 0 else next_k - i

            cone = bpy.data.objects['cone_{}_{}'.format(j,k)]
            cone.rotation_mode = 'QUATERNION'  # NEW, you actually only need to do this once of course
            cone.location = (x[i][j][m],y[i][j][m],z[i][j][m])

            cnt = Vector((cx[i][m],cy[i][m],cz[i][m]))
            pnt = Vector((x[i][j][m],y[i][j][m],z[i][j][m]))
            next_pnt = Vector((x[i][j][next_m],y[i][j][next_m],z[i][j][next_m]))

            up = pnt - cnt
            up.normalize()

            forward = next_pnt - pnt
            forward.normalize()
            
            # Compute rotation matrix from local axes
            rot = mathutils.Matrix.Identity(3)
            rot[0] = forward.cross(up)
            rot[1] = forward
            rot[2] = up
            rot.transpose()  # set columns, not rows

            rot_quat = rot.to_quaternion()              # NEW
            
            if i > 0:
                rot_quat.make_compatible(prevs[(j,k)])  # NEW
            prevs[(j,k)] = rot_quat.copy()              # NEW
            
            cone.rotation_quaternion = rot_quat         # NEW

            cone.keyframe_insert(data_path="location", frame=frame_num)
            cone.keyframe_insert(data_path="rotation_quaternion", frame=frame_num)  # NEW
    frame_num += frame_step

Is it fixed now?! Yes, I think so.

Note that I didn't do the last frame though.

$\endgroup$

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .