Definition of terms:

It you want to consider the rotation only then the code for this is shown
here.

Alternatively if we want to consider the possibility that there is also a translation
from the centre, and that the rotation may not be about the centre but may be
about some arbitrary point, then we need to extend the notation as follows:

- A 3D vector has to be added to the quaternion so that together they can
specify the position and orientation.
- The matrix has to be increased from 3x3 to 4x4 so that it can specify the
position and orientation (as explained
here).

## Rotation about a point other than origin

For the derivation of how to rotate about a point see
this page.

Quaternions and 3x3 matrices alone can only represent rotations about the origin.
But if we include a 3D vector with the quaternion we can use this to represent
the point about which we are rotating. Also if we use a 4x4 matrix then this
can hold a translation and therefore can specify a rotation about a point.

The following code generates a 4x4 matrix from a quaternion and a vector. The
derivation is given here.

void setRotate(sfquat q,sfvec3f centre) {
double sqw = q.w*q.w;
double sqx = q.x*q.x;
double sqy = q.y*q.y;
double sqz = q.z*q.z;
m00 = sqx - sqy - sqz + sqw; // since sqw + sqx + sqy + sqz =1
m11 = -sqx + sqy - sqz + sqw;
m22 = -sqx - sqy + sqz + sqw;
double tmp1 = q.x*q.y;
double tmp2 = q.z*q.w;
m01 = 2.0 * (tmp1 + tmp2);
m10 = 2.0 * (tmp1 - tmp2);
tmp1 = q.x*q.z;
tmp2 = q.y*q.w;
m02 = 2.0 * (tmp1 - tmp2);
m20 = 2.0 * (tmp1 + tmp2);
tmp1 = q.y*q.z;
tmp2 = q.x*q.w;
m12 = 2.0 * (tmp1 + tmp2);
m21 = 2.0 * (tmp1 - tmp2);
double a1,a2,a3;
if (centre == null) {
a1=a2=a3=0;
} else {
a1 = centre.x;
a2 = centre.y;
a3 = centre.z;
}
m03 = a1 - a1 * m00 - a2 * m01 - a3 * m02;
m13 = a2 - a1 * m10 - a2 * m11 - a3 * m12;
m23 = a3 - a1 * m20 - a2 * m21 - a3 * m22;
m30 = m31 = m32 = 0.0;
m33 = 1.0;
}

## Rotation about a point - using Quaternion

void setRotate(Quat4d q,sfvec3f centre) {
double sqw = q.w*q.w;
double sqx = q.x*q.x;
double sqy = q.y*q.y;
double sqz = q.z*q.z;
m00 = sqx - sqy - sqz + sqw; // since sqw + sqx + sqy + sqz =1
m11 = -sqx + sqy - sqz + sqw;
m22 = -sqx - sqy + sqz + sqw;
double tmp1 = q.x*q.y;
double tmp2 = q.z*q.w;
m01 = 2.0 * (tmp1 + tmp2);
m10 = 2.0 * (tmp1 - tmp2);
tmp1 = q.x*q.z;
tmp2 = q.y*q.w;
m02 = 2.0 * (tmp1 - tmp2);
m20 = 2.0 * (tmp1 + tmp2);
tmp1 = q.y*q.z;
tmp2 = q.x*q.w;
m12 = 2.0 * (tmp1 + tmp2);
m21 = 2.0 * (tmp1 - tmp2);
double a1,a2,a3;
if (centre == null) {
a1=a2=a3=0;
} else {
a1 = centre.x; a2 = centre.y; a3 = centre.z;
}
m03 = a1 - a1 * m00 - a2 * m01 - a3 * m02;
m13 = a2 - a1 * m10 - a2 * m11 - a3 * m12;
m23 = a3 - a1 * m20 - a2 * m21 - a3 * m22;
m30 = m31 = m32 = 0.0;
m33 = 1.0;
}

This site may have errors. Don't use for critical systems.