SOLVED Gimbal Lock Safe Target Expression

Hi,
I have seen and as expected this was discussed already in the past, but I could not find any working example in this forum.

I'd simply like to have a (gimbal lock safe) python replacement for the Target (Expression) Tag.

I'm attaching what I have come up with until now...
two versions #1 without , #2 with quaternion --> sadly #1 has gimbal issues and #2 fails completely
both with a pitch lock like in the Target Tag --> sure this could be done more elegant .-)

It was stated in one of the (below) mentioned threads, that support can't solve such problems for the users,
but this would be such a valuable tool in the quiver that it might be an exception (?)

best, index

def lookat(op, target, pitchlock=0):
    ''' gimbal lock UNSAFE target expression
        this basically works, but results in axis flips sometimes '''
    r1 = op.GetRelRot() # backup

    p1 = op.GetMg().off # from
    p2 = target.GetMg().off # to
    tmp = c4d.Vector(0, 1, 0) # temp-up
    fo = (p2 - p1).GetNormalized() # forward
    ri = tmp.Cross(fo).GetNormalized() # right
    up = fo.Cross(ri).GetNormalized() # up

    # https://plugincafe.maxon.net/topic/11152/object-target-expression-in-python/2
    # try to mitigate gimbal-lock --> UNSAFE, RESULTS IN AXIS FLIPS !
    if (fo + tmp).z < 0: up.x = -up.x
    if abs(ri.x - 1.0) < 0.001: ri = -ri

    m = c4d.Matrix()
    m.v1 = ri
    m.v2 = up
    m.v3 = fo
    m.off = p1
    op.SetMg(m)

    if pitchlock:
        r = op.GetRelRot()
        op.SetRelRot(c4d.Vector(r.x, r1.y, r.z))
def lookatQ(op, target, pitchlock=0):
    ''' quaternion based (gimbal lock safe ???) target expression
        sadly doesn't work and results in a wrong rotation '''

    r1 = op.GetRelRot() # backup

    p1 = op.GetMg().off # from
    p2 = target.GetMg().off # to
    tmp = c4d.Vector(0, 1, 0) # temp-up
    fo = (p2 - p1).GetNormalized() # forward
    ri = tmp.Cross(fo).GetNormalized() # right
    up = fo.Cross(ri).GetNormalized() # up

    # copied from
    # https://www.gamedev.net/forums/topic/613595-quaternion-lookrotationlookat-up/
    w = math.sqrt(1.0 + ri.x + up.y + fo.z) * 0.5
    w4_recip = 1.0 / (4.0 * w)
    x = (up.z - fo.y) * w4_recip
    y = (fo.x - ri.z) * w4_recip
    z = (ri.y - up.x) * w4_recip

    # given the above is correct, is this correctly applied here ???
    q = c4d.Quaternion()
    v = c4d.Vector(x, y, z)
    q.SetAxis(v, w)
    m = q.GetMatrix()
    m.off = p1
    op.SetMg(m)

    if pitchlock:
        r = op.GetRelRot()
        op.SetRelRot(c4d.Vector(r.x, r1.y, r.z))

Hey @indexofrefraction,

thank you for reaching out to us. We will answer this question more thoroughly tomorrow. In the mean time I would ask you to:

  1. Provide executable code: You provide only two functions. I can make up the context, a main function of a Python tag which calls one of these functions, but that is not the goal of us asking for executable examples. Please give us the whole tag reduced to a level where your problem still occurs. There is also a problem in the line abs(ri - 1.0) < 0.001 as abs(ri - 1.0) will evaluate to a Vector, and you cannot compare a vector against a float in our API (in fact vectors have no comparison operators at all besides __eq__ and __ne__).
  2. Explain what you mean by axis flip or gimbal lock. Unless I am overlooking something fundamental, a gimbal lock cannot occur in your scenario. What can however happen, is that your ri collapses as you do not test for tmp not being (anti)parallel to fo. It is mandatory to do this, as otherwise ri can become the null vector. The Matrix Manual has an example for that.

edit: Eh, now I see it. Riccardo started the gimbal lock talk in that thread. Unless I am overlooking some fundamental context here, both the solution and explanation provided there does not seem to be correct.

Cheers,
Ferdinand

edit 2: Here is a solution 🙂

"""Demonstrates constructing a look-at transform inside a Python tag.

Add a Python tag to an object, give the tag a single user data link field, link the target object
in the link field, and paste this code into the tag to run the example.
"""
import c4d

op: c4d.BaseTag

EPSILON = 1E-5 # The floating point precision we are going to assume, i.e., 0.00001

def GetLookAtTransform(host: c4d.Matrix, target: c4d.Matrix) -> c4d.Matrix:
    """Returns a transform which orients its z/k/v3 component from #host to #target.
    """
    # Get the position of both transforms.
    p: c4d.Vector = host.off
    q: c4d.Vector = target.off

    # The normalized delta vector will become the z-axis of the frame.
    z: c4d.Vector = ~(q - p)

    # We compute an up-vector which is not (anti)parallel to #z.
    up: c4d.Vector = c4d.Vector(0, 1, 0) \
                     if 1. - abs(z * c4d.Vector(0, 1, 0)) > EPSILON else \
                     c4d.Vector(EPSILON, 1, 0)

    # We compute the rest of the frame.
    x: c4d.Vector = ~(up % z)
    y: c4d.Vector = ~(z % x)

    # Return the frame (x, y, z) plus the offset of #host as the look-at transform.
    return c4d.Matrix(off=p, v1=x, v2=y, v3=z)


def main() -> None:
    """
    """
    # The object to orient (#host) and the object to orient #host towards (#target)
    host: c4d.BaseObject = op.GetObject()
    target: c4d.BaseObject = op[c4d.ID_USERDATA, 1]

    # Make #host look at #target
    host.SetMg(GetLookAtTransform(host.GetMg(), target.GetMg()))

Hi Ferdinand,

thanks a lot, seems I'm not the only one missing the good weather outside .-)

The code in the first lookat() function should be abs(ri.x - 1.0),
this comes from the post in the mentioned link. I corrected that in the original post.
(it was just missing in the post, my code had the ri.x)

what I'm on to is a 2 axis ceiling spot.
its a null with a python tag having a target link userdata, then inside a hierarchy with two axes.
atm they are driven by 2 target tags (initialized by the pyTag) one needs pitch on, one needs pitch off.
this works fine, but as i have 900 of these in a file,
i am trying to get rid of the target tags and do all in the pyTag.
I can post the file tomorrow.

what happens with my lookat() is that if you move the target,
the spot sometimes turns upside down 180 degrees and sticks in the ceiling.

with lookatQ() the spot simply doesn't point in the correct direction.
but it would be very interesting to learn how to do it with quaternions.

edit: your solution seems to work nicely, btw!
p.s. i guess ~ is normalise and % is the cross product

Hey @indexofrefraction,

so, the problem is solved for you now? Sorry for asking so bluntly, but despite Riccardo introducing the term gimbal lock, it is simply not the correct one, and I am struggling a bit understanding what your problem is. The most obvious answer was that you were constructing a collapsed frame and the Cinema 4D more or less randomly correcting that frame (resulting in the "flips"). The other answer would be the parallel transport problem, i.e., the question of what the direction "up" means.

I think from what was written in the other thread and what you write here, the first item is the obstacle you are encountering. What happens here is that when you have for example the following inputs:

up: Vector = (0, 1, 0)
a: Vector (200, 0, 0)
b: Vector (200, -200, 0) or 

Then (b-a) is antiparallel to up, i.e., b-a and up span exactly a 180° angle. But you want to take the cross product of b-a and up and Wikipedia says:

If two vectors have the same direction or have the exact opposite direction from each other (that is, they are not linearly independent), or if either one has zero length, then their cross product is zero.

So, if b is either exactly 'above' or 'below' a, the cross product is sort of not defined in this case (at least in the capacity we are using it here, to construct a normal to plane defined by two vectors). But this depends of course on what up is. If up would be (1, 0, 0) the same thing would happen if b is either exactly 'left' or 'right' to b. When you then construct a frame (a coordinate system) in sich case, your frame will look something like this:

{ z: ~(b-a),
  y: (0, 0, 0),
  x: (0, 0, 0) }

Which is of course not a valid frame (y is also the null vector because it depends on x and the cross product of any vector and a null vector is also the null vector) So, when you then set this frame to an object, Cinema will try to correct it (and make up the missing components), resulting in the "jump".

i guess ~ is normalise and % is the cross product

Yes, that is true.

with lookatQ() the spot simply doesn't point in the correct direction.
but it would be very interesting to learn how to do it with quaternions.

I do not want to be rude, but I do not see much value in doing it with quaternions here. If you are interested in these subjects, how to convert back and forth from matrices to quaternions (and other rotation representations), you should look at an informative book on that subject. I would recommend

Schneider, P. J. & Eberly, D. H. (2003). Appendix A Numerical Methods: Representations of 3D Rotations (p. 857-869). In: Geometric Tools for Computer Graphics. Morgan Kaufmann Publishers, Elsevier Science. doi: 10.5555/863093.

as it is a delightful book in general. The quoted chapter goes over conversion formulas between different forms of orientation/rotation representation, and the computational advantages and disadvantages of them. Quaternions are not as amazing as people often think them to be. They have a better memory footprint than rotation matrices and allow for unambiguous interpolation between two states (the gimbal lock thing) but are more expensive to compute. But matrices are usually the way to go, and they are not inherently bad.

Cheers,
Ferdinand

Hi Ferdinand, yes all solved
and sorry if the gimbal lock was a false interpretation of what was happening.
From reading the old threads I thought the quaternion solution was mandatory to solve this.

for clarity and as an example I'm including a comparison of
your solution (working) and riccardos solution (showing the flip)

spotLight.c4d

but for me the case is closed...
Thanks a lot for helping with this !

I have found another way to do this, which is maybe better (?)
because it seems to work around tackling with parallel vectors....

LookAtCamera::Execute() translated to python ... 🙂
https://github.com/PluginCafe/cinema4d_cpp_sdk_extended/blob/master/plugins/cinema4dsdk/source/tag/lookatcamera.cpp

def lookat(op, target, pitch=True):
  ''' orient op to look at target '''
  local = (~(op.GetUpMg() * op.GetFrozenMln())) * target.GetMg().off - op.GetRelPos()
  hpb = c4d.utils.VectorToHPB(local)
  if not pitch: hpb.y = op.GetRelRot().y
  hpb.z = op.GetRelRot().z
  op.SetRelRot(hpb)
  c4d.EventAdd()

It seems to work fine... but for clarity... after checking the Matrix Manual (see end of the page) ...
https://developers.maxon.net/docs/Cinema4DPythonSDK/html/manuals/data_algorithms/classic_api/matrix.html?highlight=matrix

Is this correct / covers all cases to get the targets position in local space?
(~(op.GetUpMg() * op.GetFrozenMln())) * target.GetMg().off
specifically.. what is
m = utils.MatrixMove(op.GetAbsPos()) * utils.HPBToMatrix(op.GetAbsRot()) * utils.MatrixScale(op.GetAbsScale())
in this context? (ps. GetAbsScale misses the brackets in the docs)

and... is this the real code of the Target Tag ?

Hello @indexofrefraction,

Is this correct / covers all cases to get the targets position in local space?

There is no "correct" way. The way I showed you is a/the common approach to construct a frame from a vector. A vector in R³ (opposed to a tuple of rotation values which might be stored in a type called Vector) only defines two out of the three degrees of freedom of a rotation in R³. When you construct a frame manually, you can define that third degree of freedom with the up vector (or just right away construct the frame out of two meaningful vectors which are already orthogonal to each other). VectorToHPB() makes this decision for you, as it, quote, '[will always set] the bank [...] to 0.0' because a vector does not store three angles but two. I personally would never go over the Euler angles, as this seems unnecessarily complicated, but this is a question of personal preference and not right or wrong.

Is this correct / covers all cases to get the targets position in local space?
(~(op.GetUpMg() * op.GetFrozenMln())) * target.GetMg().off [- op->GetRelPos();]

What @s_bach did there was specific to that example which operates on Euler angles. I cannot explain all the math here, but in general, the position of an object B relative to an object A, i.e, the local offset of B within the coordinate system of A is just the delta of both offsets transformed by the coordinate system of A.

A: c4d.BaseObject
B: c4d.BaseObject

mgA: c4d.Matrix = A.GetMg()
mgB: c4d.Matrix = B.GetMg()

# This is the location of B in relation to A in world coordinates.
delta: c4d.Vector = mgB.off - mgA.off

# The same point but in the coordinate system of A. If B were parented to A, it would have 
# this offset.
offBLocalToA: c4d.Vector = mgA * delta

I again would have personally done some things differently in the loot-at example, but the solution is just as good as others.

what is m = utils.MatrixMove(op.GetAbsPos()) * utils.HPBToMatrix(op.GetAbsRot()) * utils.MatrixScale(op.GetAbsScale()) in this context? (ps. GetAbsScale misses the brackets in the docs)

Please read the docs, I cannot explain here transforms from scratch again. Thanks for the hint for the missing parenthesis, I will fix that in an upcoming release. As a tip: Do not fixate too much on frozen transforms in our API, you can ignore them more or less, the relevant methods are GetMg and GetMl.

and... is this the real code of the Target Tag ?

I am not quite sure how you mean that? You mean if the LookAtCamera::Execute example is identical to what the ::Execute of Ttargetexpression does internally? No really, the expression ~(op.GetUpMg() * op.GetFrozenMln()) pops up there, but Ttargetexpression also constructs the frame manually with the cross product, but there are more steps involved than in my simple example. But again, there is no real right or wrong here, you can do things differently if you want to.

And before the question comes up, no, we unfortunately cannot share the code of the target expression. I would also point out again our Forum Rules, especially regarding the scope of support and the rules of conduct. We cannot teach you the math and you must implement your things yourself.

Cheers,
Ferdinand

PS: If you want to explore what the frozen transform means (it is mostly smoke and mirrors for the API perspective, unless you want to specifically deconstruct a transform), I would recommend adding a Python tag to objects, adding this code, and then play around with the freeze transform feature.

import c4d

op: c4d.BaseTag

def main() -> None:
    """
    """
    def PrettyMatrixString(m: c4d.Matrix):
        """_summary_
        """
        s: str = f"{m.__class__.__name__}("
        l: int = len(s)
        for i, v in enumerate((m.off, m.v1, m.v2, m.v3)):
            s += f"{' ' * (l if i != 0 else 0)}({round(v.x, 3):>8}, {round(v.y, 3):>8}, {round(v.z, 3):>8}\n"
        return s[:-1] + ")"

    node: c4d.BaseObject = op.GetMain()
    if not isinstance(node, c4d.BaseObject): # Should technically never happen
        return 

    for line in (
        f"GetMg():\n{PrettyMatrixString(node.GetMg())}",
        f"GetMl():\n{PrettyMatrixString(node.GetMl())}",
        f"GetFrozenMln():\n{PrettyMatrixString(node.GetFrozenMln())}"):
        print (line)

09c7c8b7-63d3-4020-8ef3-2881ff262003-image.png