3D image engine implementation from scratch: (9) Four Element Library

1. Mathematical Analysis

1) What is a quaternary element?

This abacus, matrix, and plural is a type of thing, that is, a mathematical tool. mathematicians have created this thing to solve some mathematical problems. In fact, the Quaternary element is a super-plural number. It is not only a plural number of virtual numbers, but a plural number of three virtual numbers. Let's review the plural.

2) Sources of virtual data

There is no square root of-1 in a real number set. Because no real number's Square is equal to-1, mathematicians create it-virtual I and define I * I =-1.

So we can calculate SQRT (-4), SQRT (-4) = 2 * I

3) plural

Definition: The sum of a real number and a virtual number. Z = a + bi.

A is called the real part and B is called the virtual part.

Where (A, B) is the point on the complex plane. This is also why 3D graphic operations can be linked to the Quaternary element.

Algorithm:

Z1 = a + B * I

Z2 = C + D * I

Multiply by a scalar K * z1 = K * A + K * B * I

Plural addition Z1 + Z2 = A + C + (B + d) * I

Multiply the numbers by z1 * Z2 = (a + B * I) * (C + D * I) = A * C + A * D * I + C * B * I-B * D = A * C-B * D + (A * D + B * C) * I

Plural division Z1/Z2 = (a + B * I)/(C + D * I)

Z x z * is a real number.

So the plural division Z1/Z2 = (a + B * I) * (c-d * I)/(C2 + D2) =... = (A * C + B * D)/(A2 + b2) + (B * C-A * D) * I)/(A2 + b2)

The last 1/z = 1/(a + B * I) of z1 can also be converted to: A/(A2 + b2) + B * I/(A2 + b2)

Z1 norm | z | = SQRT (A2 + b2) = SQRT (Z × z *)

4) superplural

Q = q0 + Q1 * I + q2 * j + Q3 * k

The above is the representation of the Quaternary element, where q0 is the real part, while the virtual part <Q1, q2, Q3> exactly forms the vector in the 3D Complex Plane.

The abbreviation is Q = q0 + QV, where QV is a vector <Q1, q2, Q3>.

The operation rule of the supercomplex number is the same as that of the plural number, so we will not repeat it here, but it must be particularly noted that the reciprocal q-1 of the Quaternary element.

Q × Q-1 = 1

Multiply both sides by Q *:

Q × Q-1 × Q * = Q *

Because: Q x Q * = | q | 2

Q-1 = Q */| q | 2

It can be noted that if Q is a unit of four elements, then the reciprocal of Q is equal to the concatenation of Q.

This feature is very important because it should be used in the rotation of the Quaternary element.

5) rotate the element

For a 3D vector <x, y, z>, we put it in the form of a quaternary number, which is: VQ = <0, x, y, z>.

Given a Q unit that represents the rotation axis and rotation angle, the result of vector VQ rotating around the specified axis and the specified angle is as follows:

Right coordinate system:

Clockwise rotation: vq' = Q ** VQ × Q

Clockwise rotation: vq' = Q × VQ × Q *

Left-hand coordinate system:

Clockwise rotation: vq' = Q × VQ × Q *

Clockwise rotation: vq' = Q ** VQ × Q

For the unit Q, which is used to store the rotation axis and rotation angle, the relationship between it and the rotation axis v = <x0, y0, z0> is as follows:

Q = cos (theta/2) + sin (theta/2) × V

2. Code Implementation

1) struct Definition

View plaincopy to clipboardprint?
Typedef struct quat_type // Quaternary
{
Union
{
Double M [4];
Struct
{
Double w, x, y, z;
};
Struct
{
Double q0;
Vector3d QV;
};
};
} Quat, * quat_ptr;
Typedef struct quat_type // Quaternary
{
Union
{
Double M [4];
Struct
{
Double w, x, y, z;
};
Struct
{
Double q0;
Vector3d QV;
};
};
} Quat, * quat_ptr;

2) implementation of common operation functions of the Quaternary Element

View plaincopy to clipboardprint?
Void _ cppyin_math: quatcreate (quat_ptr Q, vector3d_ptr V, double theta) // create a quaternary Q for rotation, and V must be a unit vector
{
Double theta_div_2 = (0.5) * Theta;
Double sin_theta = sin (theta_div_2 );

Q-> X = sin_theta * V-> X;
Q-> Y = sin_theta * V-> Y;
Q-> Z = sin_theta * V-> Z;
Q-> W = cos (theta_div_2 );
}

Void _ cppyin_math: quatgetvectorandtheta (quat_ptr Q, vector3d_ptr V, double * theta)
{
* Theta = ACOs (Q-> W );
Double sin_theta_inv = 1.0/sin (* theta );

V-> X = Q-> X * sin_theta_inv;
V-> Y = Q-> y * sin_theta_inv;
V-> Z = Q-> Z * sin_theta_inv;

* Theta * = 2;
}

Void _ cppyin_math: quatadd (quat_ptr Q1, quat_ptr Q2, quat_ptr qsum)
{
Qsum-> X = Q1-> X + q2-> X;
Qsum-> Y = Q1-> Y + q2-> Y;
Qsum-> Z = Q1-> Z + q2-> Z;
Qsum-> W = Q1-> W + q2-> W;
}

Void _ cppyin_math: quatsub (quat_ptr Q1, quat_ptr Q2, quat_ptr qdiff)
{
Qdiff-> X = Q1-> X-Q2-> X;
Qdiff-> Y = Q1-> Y-Q2-> Y;
Qdiff-> Z = Q1-> Z-Q2-> Z;
Qdiff-> W = Q1-> W-Q2-> W;
}

Void _ cppyin_math: quatconjugate (quat_ptr Q, quat_ptr qconj) // obtain the concatenation
{
Qconj-> X =-Q-> X;
Qconj-> Y =-Q-> Y;
Qconj-> Z =-Q-> Z;
Qconj-> W = Q-> W;
}

Void _ cppyin_math: quatscale (quat_ptr Q, double scale, quat_ptr QS)
{
QS-> X = scale * q-> X;
QS-> Y = scale * q-> Y;
QS-> Z = scale * q-> Z;
QS-> W = scale * q-> W;
}

Double _ cppyin_math: quatnorm (quat_ptr q)
{
Return SQRT (Q-> W * q-> W + q-> X * q-> X + q-> y * q-> Y + q-> Z * q-> Z );
}

Double _ cppyin_math: quatnorm2 (quat_ptr q)
{
Return Q-> W * q-> W + q-> X * q-> X + q-> y * q-> Y + q-> Z * q-> Z;
}

Void _ cppyin_math: quatnormalize (quat_ptr Q, quat_ptr qn)
{
Double qlength_inv = 1.0/(SQRT (Q-> W * q-> W + q-> X * q-> X + q-> y * q-> Y + q-> Z * q-> Z ));

Qn-> W = Q-> W * qlength_inv;
Qn-> X = Q-> X * qlength_inv;
Qn-> Y = Q-> y * qlength_inv;
Qn-> Z = Q-> Z * qlength_inv;
}

Void _ cppyin_math: quatunitinverse
{
Qi-> W = Q-> W;
Qi-> X =-Q-> X;
Qi-> Y =-Q-> Y;
Qi-> Z =-Q-> Z;
}

Void _ cppyin_math: quatinverse (quat_ptr Q, quat_ptr Qi) // the inverse of a non-unit ry.
{
Double norm2_inv = 1.0/(Q-> W * q-> W + q-> X * q-> X + q-> y * q-> Y + q-> Z * q-> Z );

Qi-> W = Q-> W * norm2_inv;
Qi-> X =-Q-> X * norm2_inv;
Qi-> Y =-Q-> y * norm2_inv;
Qi-> Z =-Q-> Z * norm2_inv;
}

Void _ cppyin_math: quatmul (quat_ptr Q1, quat_ptr Q2, quat_ptr qprod)
{
Double prd_0 = (Q1-> Z-Q1-> Y) * (Q2-> Y-Q2-> Z );
Double prd_1 = (Q1-> W + Q1-> X) * (Q2-> W + q2-> X );
Double prd_2 = (Q1-> W-Q1-> X) * (Q2-> Y + q2-> Z );
Double prd_3 = (Q1-> Y + Q1-> Z) * (Q2-> W-Q2-> X );
Double prd_4 = (Q1-> Z-Q1-> X) * (Q2-> X-Q2-> Y );
Double prd_5 = (Q1-> Z + Q1-> X) * (Q2-> X + q2-> Y );
Double prd_6 = (Q1-> W + Q1-> Y) * (Q2-> W-Q2-> Z );
Double prd_7 = (Q1-> W-Q1-> Y) * (Q2-> W + q2-> Z );
Double prd_8 = prd_5 + prd_6 + prd_7;
Double prd_9 = 0.5*(prd_4 + prd_8 );

Qprod-> W = prd_0 + prd_9-prd_5;
Qprod-> X = prd_1 + prd_9-prd_8;
Qprod-> Y = prd_2 + prd_9-prd_7;
Qprod-> Z = prd_3 + prd_9-prd_6;
}

Void _ cppyin_math: quatmul3 (quat_ptr Q1, quat_ptr Q2, quat_ptr Q3, quat_ptr qprod) // multiply the three elements for Rotation Transformation
{
Quat qtmp;
Quatmul (Q1, q2, & qtmp );
Quatmul (& qtmp, Q3, qprod );
}

 

This article from the csdn blog, reproduced please indicate the source: http://blog.csdn.net/cppyin/archive/2011/02/10/6177742.aspx