Skip to content

Commit

Permalink
clean up the tests, a little documentation and cleanup elsewhere
Browse files Browse the repository at this point in the history
  • Loading branch information
@froyo-np
froyo-np committed Jan 6, 2025
1 parent 6962553 commit cde3fd5
Show file tree
Hide file tree
Showing 3 changed files with 84 additions and 111 deletions.
6 changes: 3 additions & 3 deletions packages/geometry/src/axisAngle.ts
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,18 +27,18 @@ export function composeRotation(b: AxisAngle, a: AxisAngle): AxisAngle {
const cosB = Math.cos(b2)
const A = a.axis
const B = b.axis
// a2 and b2 are called half angles...
const gamma = 2 * Math.acos(cosB * cosA - (Vec3.dot(B, A) * sinB * sinA))

const D = Vec3.add(
Vec3.add(Vec3.scale(B, sinB * cosA),
Vec3.scale(A, sinA * cosB)),
Vec3.scale(Vec3.cross(B, A), sinB * sinA));
// TODO: normalization wont save us when the vector is near zero, or when the angle is k2pi.... sanitize!

const dir = Vec3.normalize(D)
if (!Vec3.finite(dir)) {
return Vec3.finite(a.axis) ? a : identity
}
return { radians: gamma, axis: Vec3.normalize(D) }
return { radians: gamma, axis: dir }
}


Expand Down
42 changes: 25 additions & 17 deletions packages/geometry/src/matrix.ts
View file
Original file line number Diff line number Diff line change
@@ -1,13 +1,11 @@
import { VectorLibFactory } from './vector'
import { Vec3, type vec3 } from './vec3'
import { type vec4 } from './vec4';
// specifically, square matricies...
type VectorConstraint = ReadonlyArray<number>;
type ColumnConstraint<V extends VectorConstraint> = ReadonlyArray<V>

// examples...

export type mat4 = readonly [vec4, vec4, vec4, vec4]

// these vec types are for internal use only

type _Vec<N extends number, E> = N extends 2 ? [E, E] :
(N extends 3 ? [E, E, E] : (
Expand All @@ -18,6 +16,11 @@ type Vec<N extends number, E> = N extends 2 ? readonly [E, E] :
N extends 4 ? readonly [E, E, E, E] : never
))

// a template for generating utils for square matricies
// note that you'll see a few 'as any' in the implementation here -
// I've been having trouble getting TS to use a literal number as a template that relates to the length of tuples
// all such cases are very short, and easy to verify as not-actually-risky

function MatrixLib<Dim extends 2 | 3 | 4>(N: Dim) {
type mColumn = _Vec<Dim, number>
type mMatrix = _Vec<Dim, mColumn>
Expand Down Expand Up @@ -45,18 +48,19 @@ function MatrixLib<Dim extends 2 | 3 | 4>(N: Dim) {
}
return arr as any;
}

const identity = (): Matrix => {
const _identity = (): mMatrix => {
const mat: mMatrix = blank();
for (let i = 0; i < N; i++) {
mat[i][i] = 1;
}
return mat;
}
const identity = (): Matrix => {
return _identity();
}
const translate = (v: Column): Matrix => {
const _mat: Matrix = identity();
const mat: mMatrix = [..._mat] as any;
mat[N - 1] = ([...v] as any)
const mat: mMatrix = _identity();
mat[N - 1] = v as mColumn;
return mat;
}
const transpose = (m: Matrix): Matrix => {
Expand All @@ -76,17 +80,21 @@ function MatrixLib<Dim extends 2 | 3 | 4>(N: Dim) {
const T = transpose(a)
return map(v, (_, i) => lib.dot(v, T[i]))
}
const normalize = (a: Matrix): Matrix => {
return map(a, (c) => lib.normalize(c))
}

const toColumnMajorArray = (m: mat4): number[] => m.flatMap(x => x)
return {
identity, mul, transpose, translate, transform, normalize, toColumnMajorArray
identity, mul, transpose, translate, transform, toColumnMajorArray
}

}
type Mat4Lib = ReturnType<typeof MatrixLib<4>>;
export function rotateAboutAxis(axis: vec3, radians: number): mat4 {

/**
* @param axis
* @param radians
* @returns a 4x4 matrix, expressing a rotation about the @param axis through the origin
*/
function rotateAboutAxis(axis: vec3, radians: number): mat4 {
const sin = Math.sin(radians);
const cos = Math.cos(radians);
const icos = 1 - cos;
Expand Down Expand Up @@ -118,8 +126,7 @@ export function rotateAboutAxis(axis: vec3, radians: number): mat4 {
const W: vec4 = [0, 0, 0, 1];
return [X, Y, Z, W];
}

export function rotateAboutPoint(lib: Mat4Lib, axis: vec3, radians: number, point: vec3): mat4 {
function rotateAboutPoint(lib: Mat4Lib, axis: vec3, radians: number, point: vec3): mat4 {
const mul = lib.mul;
const back = lib.translate([...point, 1])
const rot = rotateAboutAxis(axis, radians);
Expand All @@ -128,6 +135,7 @@ export function rotateAboutPoint(lib: Mat4Lib, axis: vec3, radians: number, poin
return mul(mul(move, rot), back);

}
const moverotmove = (lib: Mat4Lib) => (axis: vec3, radians: number, point: vec3) => rotateAboutPoint(lib, axis, radians, point)
const lib = MatrixLib(4)

export const Mat4 = { ...lib, rotateAboutAxis, rotateAboutPoint }
export const Mat4 = { ...lib, rotateAboutAxis, rotateAboutPoint: moverotmove(lib) }
147 changes: 56 additions & 91 deletions packages/geometry/src/tests/rotation.test.ts
View file
Original file line number Diff line number Diff line change
@@ -1,58 +1,9 @@
import { describe, expect, test } from "vitest";
import { Mat4, rotateAboutAxis } from "../matrix";
import { Mat4 } from "../matrix";
import { Vec3, vec3 } from "../vec3"
import { Vec4, type vec4 } from "../vec4";
import { AxisAngle, composeRotation, rotateVector } from "../axisAngle";

// TODO: delete the quaternion stuff
// A versor is a (Unit) Quaternion, for representing a rotation in 3D
export type Versor = Readonly<{
qi: number;
qj: number;
qk: number;
qr: number; // the scalar term
}>

export function versorFromAxisAngle(rotation: AxisAngle): Versor {
const sin = Math.sin(rotation.radians / 2)
const cos = Math.cos(rotation.radians / 2)
const e = Vec3.normalize(rotation.axis)
const [x, y, z] = e
return {
qi: x * sin,
qj: y * sin,
qk: z * sin,
qr: cos
}
}
export function axisAngleFromVersor(rotation: Versor): AxisAngle {
const { qi, qj, qk, qr } = rotation
const q: vec3 = [qi, qj, qk]
const D = Math.sqrt(Vec3.dot(q, q))
const theta = Math.atan2(D, qr)
const axis = Vec3.scale(q, 1 / D)
return { axis, radians: theta }
}
// rotate by q2 "after" q1
// if you read q2 x q1, then you'd call compose(q2,q1)
function compose(q2: Versor, q1: Versor): Versor {
const { qr: r2, qi: i2, qj: j2, qk: k2 } = q2
const { qr: r1, qi: i1, qj: j1, qk: k1 } = q1
// source:
const r = (r2 * r1 - i2 * i1 - j2 * j1 - k2 * k1)
const i = (r2 * i1 + i2 * r1 + j2 * k1 - k2 * j1)
const j = (r2 * j1 - i2 * k1 + j2 * r1 + k2 * i1)
const k = (r2 * k1 + i2 * j1 - j2 * i1 + k2 * r1)

return {
qi: i,
qj: j,
qk: k,
qr: r
}
}


function nearly(actual: vec3, expected: vec3) {
const dst = Vec3.length(Vec3.sub(actual, expected))
if (dst > 0.0001) {
Expand Down Expand Up @@ -131,48 +82,62 @@ describe('rotation in various ways', () => {
nearly(rotateVector(total, v), rotateVector(expectedRotation, v))
nearly(rotateVector(composeRotation(total, total), v), v)
})
describe('matrix works the same', () => {
const randomAxis = (): vec3 => {
const theta = Math.PI * 100 * Math.random()
const sigma = Math.PI * 100 * Math.random()
const x = Math.cos(theta) * Math.sin(sigma)
const y = Math.sin(theta) * Math.sin(sigma)
const z = Math.cos(sigma);
// always has length 1... I think?
return [x, y, z]
})
describe('matrix works the same', () => {
const randomAxis = (): vec3 => {
const theta = Math.PI * 100 * Math.random()
const sigma = Math.PI * 100 * Math.random()
const x = Math.cos(theta) * Math.sin(sigma)
const y = Math.sin(theta) * Math.sin(sigma)
const z = Math.cos(sigma);
// always has length 1... I think?
return Vec3.normalize([x, y, z])
}
test('rotateAboutAxis and axis angle agree (right hand rule... right?)', () => {
const axis: vec3 = Vec3.normalize([0.5904, -0.6193, -0.5175])
expect(Vec3.length(axis)).toBeCloseTo(1, 8)
const v: vec3 = [0.4998, 0.0530, 0.8645]
expect(Vec3.length(v)).toBeCloseTo(1, 3)
const angle = -Math.PI / 4
const mat = Mat4.rotateAboutAxis(axis, angle)
const aa: AxisAngle = { axis, radians: angle }
const a = rotateVector(aa, v)
const b = Vec4.xyz(Mat4.transform(mat, [...v, 0]))
nearly(b, a)
})
test('repeated rotations about random axes match the equivalent matrix rotateVector...', () => {
let v: vec3 = [1, 0, 0]
for (let i = 0; i < 300; i++) {
const axis = randomAxis()
expect(Vec3.length(axis)).toBeCloseTo(1)
const angle = (Math.PI / 360) + (Math.random() * Math.PI)
const dir = Math.random() > 0.5 ? -1 : 1
const mat = Mat4.rotateAboutAxis(axis, angle * dir)
const aa: AxisAngle = { axis, radians: dir * angle }
// rotateVector v by each
const v4: vec4 = [...v, 0]
const mResult = Mat4.transform(mat, v4)
const aaResult = rotateVector(aa, v)
nearly(Vec4.xyz(mResult), aaResult)
v = aaResult
}
test('rotateAboutAxis and axis angle agree (right hand rule... right?)', () => {
const axis: vec3 = Vec3.normalize([0.5904, -0.6193, -0.5175])
expect(Vec3.length(axis)).toBeCloseTo(1, 8)
const v: vec3 = [0.4998, 0.0530, 0.8645]
expect(Vec3.length(v)).toBeCloseTo(1, 3)
const angle = -Math.PI / 4
const mat = rotateAboutAxis(axis, angle)
const aa: AxisAngle = { axis, radians: angle }
const a = rotateVector(aa, v)
const b = Vec4.xyz(Mat4.transform(mat, [...v, 0]))
nearly(b, a)
})
test('repeated rotations about random axes match the equivalent matrix rotateVector...', () => {
let v: vec3 = [1, 0, 0]
for (let i = 0; i < 300; i++) {
const axis = randomAxis()
expect(Vec3.length(axis)).toBeCloseTo(1)
const angle = (Math.PI / 360) + (Math.random() * Math.PI)
const dir = Math.random() > 0.5 ? -1 : 1
const mat = rotateAboutAxis(axis, angle * dir)
const aa: AxisAngle = { axis, radians: dir * angle }
// rotateVector v by each
// console.log('-------------------------------------')
// console.log(`rotate (${v[0].toFixed(4)}, ${v[1].toFixed(4)}, ${v[2].toFixed(4)}) about`)
// console.log(`<${axis[0].toFixed(4)}, ${axis[1].toFixed(4)}, ${axis[2].toFixed(4)}> by ${(angle * dir).toFixed(5)} `)
const v4: vec4 = [...v, 0]
const mResult = Mat4.transform(mat, v4)
const aaResult = rotateVector(aa, v)
nearly(Vec4.xyz(mResult), aaResult)
v = aaResult
}
})
})
})
describe('rotation about a point which is not the origin', () => {
test('an easy to understand example', () => {
let v: vec4 = [1, 0, 0, 1]
let yAxis: vec3 = [0, 1, 0]
let origin: vec3 = [2, 0, 0]

// o----v---|---x
// 0----1---2---3---4---

// o = the true origin
// v = the point v
// | = the y-axis through the point [2,0,0]
// x = imagine rotating v around the | by 180, you get x
const rot = Mat4.rotateAboutPoint(yAxis, Math.PI, origin)
nearly(Vec4.xyz(Mat4.transform(rot, v)), [3, 0, 0])
})
})

Expand Down

0 comments on commit cde3fd5

Please sign in to comment.