Is this really worth the trouble?
I don't know. I asked the same question myself, because the light emitted by the phosphorescent surface is not a linear function of the beam intensity or dwell time at all. That is, some variance in the beam intensity or velocity on the surface likely causes much smaller variance in the light emitted.
Mr. Margolin has other vector graphics and vector math related papers on his site.
Thanks; I found them very interesting – except for the Euler angles for 3D rotation, which I consider a loaded footgun with a hair trigger loaded with bullets constructed by a gun nut using very unstable explosives.
What I like to use, is either unit quaternions or bivectors; the actual math turns out to be the same in both cases.
The idea is rotation is described by an unit axis vector \$(x, y, z)\$, \$x^2 + y^2 + z^2 = 1\$, and a rotation around that axis by angle \$2\varphi\$, using four variables:
$$w = \cos \varphi, \quad i = x \sin \varphi, \quad j = y \sin \varphi, \quad k = z \sin \varphi$$
This can be described as a versor (unit quaternion) or as a three-dimensional bivector, with \$w^2 + i^2 + j^2 + k^2 = 1\$. You can always divide the four variables by \$\sqrt{w^2 + i^2 + j^2 + k^2}\$ to normalize it to unit length. Unlike with e.g. rotation matrices, this causes no distortion.
No rotation can be described by \$(w=1, i=0, j=0, k=0)\$ (which is recommended to be used when the sum of the squared components is too close to zero when trying to normalize it to unit length).
The inverse of rotation \$(w, i, j ,k)\$ is \$(w, -i, -j, -k)\$ or \$(-w, i, j, k)\$. When \$w\$ is positive, the rotation is less than 180° ("short way around"); when \$w\$ is negative, the rotation is greater than 180° ("long way around"). \$(w, i, j, k)\$ and \$(-w, -i, -j, -k)\$ represent the same orientation; as rotations, they are the two opposite ways of achieving the same orientation by rotating around the same axis.
To rotate an orientation, you apply
Hamilton product, with the initial orientation rightmost, and applied rotations progressing to the left. (That is, last rotation is leftmost, initial orientation rightmost.) It is 16 real multiplications and 12 additions or subtractions. If we use \$(w_0, i_0, j_0, k_0)\$ for the orientation, and \$(w_1, i_1, j_1, k_1)\$ for the rotation to be applied to it, the result is
$$\begin{aligned}
w &= w_1 w_0 - i_1 i_0 - j_1 j_0 - k_1 k_0 \\
i &= w_1 i_0 + i_1 w_0 + j_1 k_0 - k_1 j_0 \\
j &= w_1 j_0 - i_1 k_0 + j_1 w_0 + k_1 i_0 \\
k &= w_1 k_0 + i_1 j_0 - j_1 i_0 + k_1 w_0 \\
\end{aligned}$$
If you normalize the final four to unit length before turning it into a matrix, you can chain any number of rotations without any issues, unlike when using matrices (where tiny rounding errors will compound).
To apply a rotation or orientation to 3D points, you normalize the four to unit length, then construct a 3×3 rotation matrix using
$$\mathbf{R} = \left[ \begin{matrix}
1 - 2(j^2 + k^2) & 2(i j - k w) & 2(i k + j w) \\
2(i j + k w) & 1 - 2(i^2 + k^2) & 2(j k - i w) \\
2(i k - j w) & 2(j k + i w) & 1 - 2(i^2 + j^2) \\
\end{matrix}\right]$$
and do matrix-vector multiplication with column vectors on the right side, \$\mathbf{R}\vec{v}\$.
For camera movement and similar, you can interpolate between two versors/bivectors by interpolating the four components, and normalizing the result. It approximates the minimum or maximum rotation path. Let's assume you want to interpolate between two orientations \$(w_0, i_0, j_0, k_0)\$ and \$(w_1, i_1, j_1, k_1)\$. First, choose the minimum rotation path by negating all four components if the \$w\$ component is negative, for both orientations. Then, linearly interpolate the components using \$0 \le \lambda \le 1\$,
$$w = (1-\lambda) w_0 + \lambda w_1, \quad i = (1-\lambda) i_0 + \lambda i_1, \quad j = (1-\lambda) j_0 + \lambda j_1, \quad k = (1-\lambda) k_0 + \lambda k_1$$
and normalize the result.