[Python] Setting joint hierarchy rotations to a new bind pose with quaterions

Developing a set of opensource BVH (Biovision Heuristics) animation tools and I've run into a snag I can't seem to solve. I'd like to convert a 'flipped' skeleton to the more common format.

The BVH whitepaper is not too specific and has resulted in 2 different implementations.

I can render both but I'd really like to translate one to the other and do away with case-based rendering.

The following is a small recursive function that builds the world positions for the joints with the more-common and flipped implementations. The only difference is whether we use the cumulative rotation or its conjugate.

(Note: I for clarity I omitted a bunch of if joint.parent is None: logic where the position | rotation is defaulted to 0,0,0)

self and parent are type Joint. w_position is world position of the joint (parent relative position + world position of parent) position is the parent-relative position of the joint. joint.rotation is a parent-relative rotation rotation (the parameter) is the cumulative rotation of the joint's ancestors (world relative rotation)

def compute_frame_wpos( self, frame_no : int, flipped=False, rotation : glm.quat = glm.quat(glm.vec3(0,0,0)) ):
'''This function sets the world position of the joints for a given
frame. Recursively calls its children, passing them the cumulative
rotations of their ancestors'''

if flipped:
self.w_position = self.parent.w_position + \
self.position * rotation
else:
self.w_position = self.parent.w_position + \
self.position * glm.conjugate(rotation)

#Recurse
child_rotation = self.frames[frame_no].rotation * rotation
for child in self.children:
child.compute_frame_wpos( frame_no, flipped, child_rotation )


The initial (parent relative) positions of the two skeletons are different, the flipped one rendering as folded up like a pretzel for the bind pose. What I would like to do is use the 0th frame (the resting pose) in place of the bind pose and correct all the rotations.

I then created an unfolded pose for the flipped skeleton by calling compute_frame_wpos for the 0th (resting) frame, then applying the following recursion to build new local positions from world positions.

def _make_relative_from_world( self ):
'''A recursive call to set the joint's parent relative positions from the world positions.'''

self.position = self.w_position - self.parent.w_position

for child in joint.children:
child._make_relative_from_world( )


Finally, I got stuck, trying to fix the rotations. In my calling code, I loop through the frames like this:

Skeleton is an encapsulating class containing the joint hierarchy and an alias list. skeleton contains the original skeleton with its imported parent-relative joint positions. ref_skeleton contains the new unfolded skeleton with the joint positions and world_positions described above.

for frame_no in range(0, self.skeleton.num_frames):
self.skeleton.get_root().compute_frame_wpos( frame_no, True )
self._reset_resting_rotation( ref_skel, frame_no )


My naive implementation employs a shorest-arc from the frame's world position for a joint to the world position for the new pose.

def _reset_frame_rotation( self, reference_skeleton : Skeleton, frame_no : int ):
'''Recursively change frame's rotations to be relative the reference skeleton'''

#Get the joint matching this one in the reference_skeleton
ref_joint = reference_skeleton.joints[ self.name ]

frame_pos = self.w_position - self.parent.w_position

self.frames[frame_no].rotation = shortest_arc( frame_pos, ref_joint.position )

for child in self.children:
child._reset_frame_rotation( reference_skeleton, frame_no )


Pretty sure the shortest_arc implementation is correct but I'll include it, just in case:

def shortest_arc(point_a : glm.vec3, point_b : glm.vec3) -> glm.quat:
'''Given the vectors a and b return a quaternion representing the shortest arc between'''

result = glm.quat(1.0,0.0,0.0,0.0)

ab_dot = glm.dot(point_a,point_b) #dot product
cross = glm.cross(point_a,point_b) #cross product
cross_sqr = glm.dot(cross,cross) #Squared length of the crossproduct

#Test if args have sufficient magnitude
if ab_dot*ab_dot+cross_sqr != 0.0:
#Test if the arguments are (anti)parllel
if cross_sqr > 0.0:
sqr = math.sqrt(ab_dot * ab_dot + cross_sqr) + ab_dot

#Inverted magnitude of the quaternion.
mag = 1.0 / math.sqrt(cross_sqr + sqr * sqr)
result = glm.quat(sqr * mag, cross[0] * mag, cross[1] * mag, cross[2] * mag)
elif ab_dot < 0.0: #Test if the angle is > than PI/2 (anti parallel)
#The arguments are anti-parallel, we have to choose an axis
point_c = point_a - point_b

#The length is projected on the XY plane
mag = math.sqrt( point_c[0] * point_c[0] + point_c[1] * point_c[1] )

#F32 magnitude threshhold. Probably not needed.
f32_mag_thresh = 0.0000001
if mag > f32_mag_thresh:
result = glm.quat( 0.0, -point_c[1]/mag, point_c[0]/mag, 0.0 )
else:
result = glm.quat( 0.0, 1.0, 0.0, 0.0 )
return result


In writing this, I realize I probably need to get the world rotations of the hierarchy from both the reference and final skeleton and subtract their difference from my final frame rotation and will try that and update if it works.

Continue reading...