Rotation: using quaternions in Igor Pro

Quaternions are essentially magic. With just four numbers, they can describe the position and orientation of an object in space. For example, the rotation of an object about one axis. Consecutive rotations can also be described by a quaternion.

There are many explainer articles which help to understand quaternions and their uses (they are used extensively in computer graphics, particularly in gaming). Quaternions are tough to understand but when you do, you can see why Hamilton, the Irish Mathematician who discovered them, engraved his eureka moment on the side of a bridge in Dublin.

In Igor, 3D objects are displayed in a Gizmo window. There are three keywords that allow rotation: setQuaternion, setRotationMatrix and euler. The focus here is on setQuaternion. Note that it is possible to use appendRotation to add successive rotations (specified using a quaternion). Here I am going to show some code to do three distinct rotations in a single quaternion.

First let’s make an object:

The code to make this object is below:

Function MakeATestObject()
	String gizName = "testGiz"
	KillWindow/Z $gizName
	NewGizmo/N=$gizName
	AppendToGizmo/N=$gizName/D pieWedge=pieWedge0
	ModifyGizmo/N=$gizName ModifyObject=pieWedge0,objectType=pieWedge,property={ calcNormals,1}
	ModifyGizmo/N=$gizName ModifyObject=pieWedge0,objectType=pieWedge,property={ rMin,0.3}
	ModifyGizmo/N=$gizName ModifyObject=pieWedge0,objectType=pieWedge,property={ colorType,2}
	ModifyGizmo/N=$gizName ModifyObject=pieWedge0,objectType=pieWedge,property={ topRGBA,1,0.499947,0.250019,1}
	ModifyGizmo/N=$gizName ModifyObject=pieWedge0,objectType=pieWedge,property={ sidesRGBA,0.533333,0.533333,0.533333,1}
	AppendToGizmo/N=$gizName/D Axes=BoxAxes,name=axes0
	ModifyGizmo/N=$gizName setOuterBox={-1,1,-1,1,-1,1}
	ModifyGizmo/N=$gizName scalingOption=0
	ModifyGizmo/N=$gizName ModifyObject=pieWedge0,objectType=pieWedge,property={ zMin,0}
	ModifyGizmo/N=$gizName ModifyObject=pieWedge0,objectType=pieWedge,property={ zMax,-1}
//	AppendToGizmo/N=$gizName/D freeAxesCue={0,0,0,1},name=freeAxesCue0
	ModifyGizmo/N=$gizName showAxisCue=1
End

Now that we have an object we can rotate it using a command like this:

ModifyGizmo/N=$gizName setQuaternion={0.5,0.5,0.5,0.5}

This is where x,y,z and w terms of the quaternion are 0.5. This will rotate the object such that the Z-axis points up, the Y-axis points right and the X-axis points forward. But how do we specify those terms?

A function to set the four terms is given below. First, a 4 point wave needs to be made and fed to the function specifying the vector of the axis for rotation, and the angle. So for example, MakeQuaternion(1,0,0,pi/2,wave0) will put into wave0, the rotation of 90 degrees CCW about the X-axis.

Function MakeQuaternion(Ax,Ay,Az,theta,w)
	Variable Ax,Ay,Az,theta
	Wave w
	Variable N = sqrt(Ax^2 + Ay^2 + Az^2)
	w[0] = Ax * sin(theta/2)/N
	w[1] = Ay * sin(theta/2)/N
	w[2] = Az * sin(theta/2)/N
	w[3] = cos(theta/2)
End

As mentioned above quaternions can be combined in a product quaternion. This code is modified from the Igor Help File.

Function MultiplyQuaternions(q1,q2,qr)
	Wave q1,q2,qr
	
	Variable w2 = q2[3]
	Variable w1 = q1[3]
	qr[3] = w2 * w1 - (q2[0] * q1[0] + q2[1] * q1[1] + q2[2] * q1[2])
	Make/N=4/FREE vcross = 0
 	vcross[0] = (q1[1] * q2[2]) - (q1[2] * q2[1])
 	vcross[1] = (q1[2] * q2[0]) - (q1[0] * q2[2])
 	vcross[2] = (q1[0] * q2[1]) - (q1[1] * q2[0])
 	MatrixOP/FREE aa = w2 * q1 + w1 * q2 + vcross
	qr[0] = aa[0]
	qr[1] = aa[1]
	qr[2] = aa[2]
	Variable NN = norm(qr)
	qr /= NN
End

In that function, three waves are specified. q1 and q2 contain the first and second rotations as quaternions and the result goes into qr.

What is really neat is that we can combine two quaternions, but then also use the result with a further quaternion to do three rotations. This gives us a single quaternion that specifies the desired position of our object.

Function TestIt(GizName)
	String GizName
	Make/O/N=4/FREE p1,p2,pr
	MakeQuaternion(1,0,0,pi/2,p1) // 90 CCW about X
	MakeQuaternion(0,1,0,0.75 * pi,p2) // 135 CCW about Y
	MultiplyQuaternions(p1,p2,pr)
	p1[] = pr[p]
	MakeQuaternion(1,0,0,1.75 * pi,p2) // 315 CCW about X
	MultiplyQuaternions(p1,p2,pr)
	ModifyGizmo/N=$gizName setQuaternion={pr[0],pr[1],pr[2],pr[3]}
End

Notice how p1 gets reset as the result of the combination of p1 and p2, then combined with a third rotation (in p2) to give the final result, pr.

So TestIt() applies a quaternion that is the product of three rotations. The first is 90 degrees CCW about X, the next is 135 degrees CCW about Y (after the first rotation) and finally 315 degrees CCW about X (the new X axis after the previous two rotations). The results is the piewedge now being in the front lower middle of the cube.

This code is part of a larger project to visualise point sets in space. It needs programmatic handling of the rotations to capture images of many point sets. Quaternions are tough to understand but give excellent control over these objects. In the Igor Help File under ModifyGizmo General Settings, it offers the following advice:

You may find it easier to use the euler keyword.

Hopefully this post will help other programmers out there.

—

The post title comes from “Rotation” by Tom Verlaine off his “Cover” LP.