use crate::types::{Vector2, Vector3};
use std::ops::Mul;
#[derive(Debug, Clone, Copy, PartialEq)]
#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
pub struct Matrix3 {
pub m: [[f64; 3]; 3],
}
impl Matrix3 {
pub fn identity() -> Self {
Self {
m: [
[1.0, 0.0, 0.0],
[0.0, 1.0, 0.0],
[0.0, 0.0, 1.0],
],
}
}
pub fn zero() -> Self {
Self {
m: [[0.0; 3]; 3],
}
}
pub fn from_rows(row0: [f64; 3], row1: [f64; 3], row2: [f64; 3]) -> Self {
Self {
m: [row0, row1, row2],
}
}
pub fn rotation_z(angle: f64) -> Self {
let cos = angle.cos();
let sin = angle.sin();
Self {
m: [
[cos, -sin, 0.0],
[sin, cos, 0.0],
[0.0, 0.0, 1.0],
],
}
}
pub fn scaling(sx: f64, sy: f64, sz: f64) -> Self {
Self {
m: [
[sx, 0.0, 0.0],
[0.0, sy, 0.0],
[0.0, 0.0, sz],
],
}
}
pub fn arbitrary_axis(normal: Vector3) -> Self {
const ARBITRARY_AXIS_THRESHOLD: f64 = 1.0 / 64.0;
let normal = normal.normalize();
let ax = if normal.x.abs() < ARBITRARY_AXIS_THRESHOLD
&& normal.y.abs() < ARBITRARY_AXIS_THRESHOLD
{
Vector3::new(0.0, 1.0, 0.0) } else {
Vector3::new(0.0, 0.0, 1.0) };
let x_dir = ax.cross(&normal).normalize();
let y_dir = normal.cross(&x_dir).normalize();
Self {
m: [
[x_dir.x, y_dir.x, normal.x],
[x_dir.y, y_dir.y, normal.y],
[x_dir.z, y_dir.z, normal.z],
],
}
}
pub fn transpose(&self) -> Self {
Self {
m: [
[self.m[0][0], self.m[1][0], self.m[2][0]],
[self.m[0][1], self.m[1][1], self.m[2][1]],
[self.m[0][2], self.m[1][2], self.m[2][2]],
],
}
}
pub fn determinant(&self) -> f64 {
self.m[0][0] * (self.m[1][1] * self.m[2][2] - self.m[1][2] * self.m[2][1])
- self.m[0][1] * (self.m[1][0] * self.m[2][2] - self.m[1][2] * self.m[2][0])
+ self.m[0][2] * (self.m[1][0] * self.m[2][1] - self.m[1][1] * self.m[2][0])
}
pub fn inverse(&self) -> Option<Self> {
let det = self.determinant();
if det.abs() < 1e-10 {
return None;
}
let inv_det = 1.0 / det;
Some(Self {
m: [
[
(self.m[1][1] * self.m[2][2] - self.m[1][2] * self.m[2][1]) * inv_det,
(self.m[0][2] * self.m[2][1] - self.m[0][1] * self.m[2][2]) * inv_det,
(self.m[0][1] * self.m[1][2] - self.m[0][2] * self.m[1][1]) * inv_det,
],
[
(self.m[1][2] * self.m[2][0] - self.m[1][0] * self.m[2][2]) * inv_det,
(self.m[0][0] * self.m[2][2] - self.m[0][2] * self.m[2][0]) * inv_det,
(self.m[0][2] * self.m[1][0] - self.m[0][0] * self.m[1][2]) * inv_det,
],
[
(self.m[1][0] * self.m[2][1] - self.m[1][1] * self.m[2][0]) * inv_det,
(self.m[0][1] * self.m[2][0] - self.m[0][0] * self.m[2][1]) * inv_det,
(self.m[0][0] * self.m[1][1] - self.m[0][1] * self.m[1][0]) * inv_det,
],
],
})
}
pub fn transform_point(&self, v: Vector3) -> Vector3 {
Vector3::new(
self.m[0][0] * v.x + self.m[0][1] * v.y + self.m[0][2] * v.z,
self.m[1][0] * v.x + self.m[1][1] * v.y + self.m[1][2] * v.z,
self.m[2][0] * v.x + self.m[2][1] * v.y + self.m[2][2] * v.z,
)
}
}
impl Mul for Matrix3 {
type Output = Self;
fn mul(self, rhs: Self) -> Self::Output {
let mut result = Self::zero();
for i in 0..3 {
for j in 0..3 {
for k in 0..3 {
result.m[i][j] += self.m[i][k] * rhs.m[k][j];
}
}
}
result
}
}
impl Mul<Vector3> for Matrix3 {
type Output = Vector3;
fn mul(self, v: Vector3) -> Self::Output {
self.transform_point(v)
}
}
impl Default for Matrix3 {
fn default() -> Self {
Self::identity()
}
}
#[derive(Debug, Clone, Copy, PartialEq)]
#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
pub struct Matrix4 {
pub m: [[f64; 4]; 4],
}
impl Matrix4 {
pub fn identity() -> Self {
Self {
m: [
[1.0, 0.0, 0.0, 0.0],
[0.0, 1.0, 0.0, 0.0],
[0.0, 0.0, 1.0, 0.0],
[0.0, 0.0, 0.0, 1.0],
],
}
}
pub fn zero() -> Self {
Self { m: [[0.0; 4]; 4] }
}
pub fn translation(tx: f64, ty: f64, tz: f64) -> Self {
Self {
m: [
[1.0, 0.0, 0.0, tx],
[0.0, 1.0, 0.0, ty],
[0.0, 0.0, 1.0, tz],
[0.0, 0.0, 0.0, 1.0],
],
}
}
pub fn scaling(sx: f64, sy: f64, sz: f64) -> Self {
Self {
m: [
[sx, 0.0, 0.0, 0.0],
[0.0, sy, 0.0, 0.0],
[0.0, 0.0, sz, 0.0],
[0.0, 0.0, 0.0, 1.0],
],
}
}
pub fn scaling_with_origin(sx: f64, sy: f64, sz: f64, origin: Vector3) -> Self {
Self::translation(origin.x, origin.y, origin.z)
* Self::scaling(sx, sy, sz)
* Self::translation(-origin.x, -origin.y, -origin.z)
}
pub fn rotation_x(angle: f64) -> Self {
let cos = angle.cos();
let sin = angle.sin();
Self {
m: [
[1.0, 0.0, 0.0, 0.0],
[0.0, cos, -sin, 0.0],
[0.0, sin, cos, 0.0],
[0.0, 0.0, 0.0, 1.0],
],
}
}
pub fn rotation_y(angle: f64) -> Self {
let cos = angle.cos();
let sin = angle.sin();
Self {
m: [
[cos, 0.0, sin, 0.0],
[0.0, 1.0, 0.0, 0.0],
[-sin, 0.0, cos, 0.0],
[0.0, 0.0, 0.0, 1.0],
],
}
}
pub fn rotation_z(angle: f64) -> Self {
let cos = angle.cos();
let sin = angle.sin();
Self {
m: [
[cos, -sin, 0.0, 0.0],
[sin, cos, 0.0, 0.0],
[0.0, 0.0, 1.0, 0.0],
[0.0, 0.0, 0.0, 1.0],
],
}
}
pub fn rotation(axis: Vector3, angle: f64) -> Self {
let axis = axis.normalize();
let cos = angle.cos();
let sin = angle.sin();
let one_minus_cos = 1.0 - cos;
let x = axis.x;
let y = axis.y;
let z = axis.z;
Self {
m: [
[
cos + x * x * one_minus_cos,
x * y * one_minus_cos - z * sin,
x * z * one_minus_cos + y * sin,
0.0,
],
[
y * x * one_minus_cos + z * sin,
cos + y * y * one_minus_cos,
y * z * one_minus_cos - x * sin,
0.0,
],
[
z * x * one_minus_cos - y * sin,
z * y * one_minus_cos + x * sin,
cos + z * z * one_minus_cos,
0.0,
],
[0.0, 0.0, 0.0, 1.0],
],
}
}
pub fn from_matrix3(m3: Matrix3) -> Self {
Self {
m: [
[m3.m[0][0], m3.m[0][1], m3.m[0][2], 0.0],
[m3.m[1][0], m3.m[1][1], m3.m[1][2], 0.0],
[m3.m[2][0], m3.m[2][1], m3.m[2][2], 0.0],
[0.0, 0.0, 0.0, 1.0],
],
}
}
pub fn to_matrix3(&self) -> Matrix3 {
Matrix3 {
m: [
[self.m[0][0], self.m[0][1], self.m[0][2]],
[self.m[1][0], self.m[1][1], self.m[1][2]],
[self.m[2][0], self.m[2][1], self.m[2][2]],
],
}
}
pub fn transpose(&self) -> Self {
let mut result = Self::zero();
for i in 0..4 {
for j in 0..4 {
result.m[i][j] = self.m[j][i];
}
}
result
}
pub fn transform_point(&self, v: Vector3) -> Vector3 {
let w = self.m[3][0] * v.x + self.m[3][1] * v.y + self.m[3][2] * v.z + self.m[3][3];
let w = if w.abs() < 1e-10 { 1.0 } else { w };
Vector3::new(
(self.m[0][0] * v.x + self.m[0][1] * v.y + self.m[0][2] * v.z + self.m[0][3]) / w,
(self.m[1][0] * v.x + self.m[1][1] * v.y + self.m[1][2] * v.z + self.m[1][3]) / w,
(self.m[2][0] * v.x + self.m[2][1] * v.y + self.m[2][2] * v.z + self.m[2][3]) / w,
)
}
pub fn transform_direction(&self, v: Vector3) -> Vector3 {
Vector3::new(
self.m[0][0] * v.x + self.m[0][1] * v.y + self.m[0][2] * v.z,
self.m[1][0] * v.x + self.m[1][1] * v.y + self.m[1][2] * v.z,
self.m[2][0] * v.x + self.m[2][1] * v.y + self.m[2][2] * v.z,
)
}
}
impl Mul for Matrix4 {
type Output = Self;
fn mul(self, rhs: Self) -> Self::Output {
let mut result = Self::zero();
for i in 0..4 {
for j in 0..4 {
for k in 0..4 {
result.m[i][j] += self.m[i][k] * rhs.m[k][j];
}
}
}
result
}
}
impl Default for Matrix4 {
fn default() -> Self {
Self::identity()
}
}
#[derive(Debug, Clone, Copy, PartialEq)]
#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
pub struct Transform {
pub matrix: Matrix4,
}
impl Transform {
pub fn identity() -> Self {
Self {
matrix: Matrix4::identity(),
}
}
pub fn from_matrix(matrix: Matrix4) -> Self {
Self { matrix }
}
pub fn from_rotation(axis: Vector3, angle: f64) -> Self {
Self {
matrix: Matrix4::rotation(axis, angle),
}
}
pub fn from_translation(translation: Vector3) -> Self {
Self {
matrix: Matrix4::translation(translation.x, translation.y, translation.z),
}
}
pub fn from_scale(scale: f64) -> Self {
Self {
matrix: Matrix4::scaling(scale, scale, scale),
}
}
pub fn from_scaling(scale: Vector3) -> Self {
Self {
matrix: Matrix4::scaling(scale.x, scale.y, scale.z),
}
}
pub fn from_scaling_with_origin(scale: Vector3, origin: Vector3) -> Self {
Self {
matrix: Matrix4::scaling_with_origin(scale.x, scale.y, scale.z, origin),
}
}
pub fn apply(&self, point: Vector3) -> Vector3 {
self.matrix.transform_point(point)
}
pub fn apply_rotation(&self, direction: Vector3) -> Vector3 {
self.matrix.transform_direction(direction)
}
pub fn then(&self, other: &Transform) -> Transform {
Transform {
matrix: other.matrix * self.matrix,
}
}
pub fn compose(&self, other: &Transform) -> Transform {
Transform {
matrix: self.matrix * other.matrix,
}
}
pub fn from_mirror_x() -> Self {
Self {
matrix: Matrix4::scaling(-1.0, 1.0, 1.0),
}
}
pub fn from_mirror_y() -> Self {
Self {
matrix: Matrix4::scaling(1.0, -1.0, 1.0),
}
}
pub fn from_mirror_z() -> Self {
Self {
matrix: Matrix4::scaling(1.0, 1.0, -1.0),
}
}
pub fn from_mirror_plane(point: Vector3, normal: Vector3) -> Self {
let n = normal.normalize();
let d = -(n.x * point.x + n.y * point.y + n.z * point.z);
Self {
matrix: Matrix4 {
m: [
[1.0 - 2.0 * n.x * n.x, -2.0 * n.x * n.y, -2.0 * n.x * n.z, -2.0 * n.x * d],
[-2.0 * n.y * n.x, 1.0 - 2.0 * n.y * n.y, -2.0 * n.y * n.z, -2.0 * n.y * d],
[-2.0 * n.z * n.x, -2.0 * n.z * n.y, 1.0 - 2.0 * n.z * n.z, -2.0 * n.z * d],
[0.0, 0.0, 0.0, 1.0],
],
},
}
}
pub fn from_mirror_line(p1: Vector3, p2: Vector3) -> Self {
let dir = (p2 - p1).normalize();
let normal = Vector3::new(-dir.y, dir.x, 0.0);
Self::from_mirror_plane(p1, normal)
}
pub fn is_identity(&self) -> bool {
*self == Self::identity()
}
}
impl Default for Transform {
fn default() -> Self {
Self::identity()
}
}
impl Mul for Transform {
type Output = Self;
fn mul(self, rhs: Self) -> Self::Output {
self.compose(&rhs)
}
}
pub fn rotate_point_2d(point: Vector2, center: Vector2, angle: f64) -> Vector2 {
let cos = angle.cos();
let sin = angle.sin();
let dx = point.x - center.x;
let dy = point.y - center.y;
Vector2::new(
center.x + dx * cos - dy * sin,
center.y + dx * sin + dy * cos,
)
}
pub fn is_zero_angle(angle: f64) -> bool {
angle.abs() < 1e-10
}
pub fn normalize_angle(angle: f64) -> f64 {
let two_pi = 2.0 * std::f64::consts::PI;
let mut a = angle % two_pi;
if a < 0.0 {
a += two_pi;
}
a
}
#[cfg(test)]
mod tests {
use super::*;
use std::f64::consts::PI;
#[test]
fn test_matrix3_identity() {
let m = Matrix3::identity();
let v = Vector3::new(1.0, 2.0, 3.0);
let result = m * v;
assert!((result.x - v.x).abs() < 1e-10);
assert!((result.y - v.y).abs() < 1e-10);
assert!((result.z - v.z).abs() < 1e-10);
}
#[test]
fn test_matrix3_rotation_z() {
let m = Matrix3::rotation_z(PI / 2.0);
let v = Vector3::new(1.0, 0.0, 0.0);
let result = m * v;
assert!(result.x.abs() < 1e-10);
assert!((result.y - 1.0).abs() < 1e-10);
}
#[test]
fn test_matrix4_translation() {
let m = Matrix4::translation(1.0, 2.0, 3.0);
let v = Vector3::new(0.0, 0.0, 0.0);
let result = m.transform_point(v);
assert!((result.x - 1.0).abs() < 1e-10);
assert!((result.y - 2.0).abs() < 1e-10);
assert!((result.z - 3.0).abs() < 1e-10);
}
#[test]
fn test_matrix4_rotation() {
let m = Matrix4::rotation(Vector3::new(0.0, 0.0, 1.0), PI / 2.0);
let v = Vector3::new(1.0, 0.0, 0.0);
let result = m.transform_point(v);
assert!(result.x.abs() < 1e-10);
assert!((result.y - 1.0).abs() < 1e-10);
}
#[test]
fn test_transform_composition() {
let t1 = Transform::from_translation(Vector3::new(1.0, 0.0, 0.0));
let t2 = Transform::from_translation(Vector3::new(0.0, 1.0, 0.0));
let combined = t1.then(&t2);
let origin = Vector3::new(0.0, 0.0, 0.0);
let result = combined.apply(origin);
assert!((result.x - 1.0).abs() < 1e-10);
assert!((result.y - 1.0).abs() < 1e-10);
}
#[test]
fn test_arbitrary_axis() {
let m = Matrix3::arbitrary_axis(Vector3::new(0.0, 0.0, 1.0));
let det = m.determinant();
assert!((det - 1.0).abs() < 1e-10); }
#[test]
fn test_mirror_x() {
let t = Transform::from_mirror_x();
let p = Vector3::new(3.0, 5.0, 7.0);
let r = t.apply(p);
assert!((r.x - (-3.0)).abs() < 1e-10);
assert!((r.y - 5.0).abs() < 1e-10);
assert!((r.z - 7.0).abs() < 1e-10);
}
#[test]
fn test_mirror_y() {
let t = Transform::from_mirror_y();
let p = Vector3::new(3.0, 5.0, 7.0);
let r = t.apply(p);
assert!((r.x - 3.0).abs() < 1e-10);
assert!((r.y - (-5.0)).abs() < 1e-10);
assert!((r.z - 7.0).abs() < 1e-10);
}
#[test]
fn test_mirror_z() {
let t = Transform::from_mirror_z();
let p = Vector3::new(3.0, 5.0, 7.0);
let r = t.apply(p);
assert!((r.x - 3.0).abs() < 1e-10);
assert!((r.y - 5.0).abs() < 1e-10);
assert!((r.z - (-7.0)).abs() < 1e-10);
}
#[test]
fn test_mirror_plane_through_origin() {
let t = Transform::from_mirror_plane(Vector3::ZERO, Vector3::UNIT_X);
let p = Vector3::new(3.0, 5.0, 7.0);
let r = t.apply(p);
assert!((r.x - (-3.0)).abs() < 1e-10);
assert!((r.y - 5.0).abs() < 1e-10);
assert!((r.z - 7.0).abs() < 1e-10);
}
#[test]
fn test_mirror_plane_offset() {
let t = Transform::from_mirror_plane(Vector3::new(2.0, 0.0, 0.0), Vector3::UNIT_X);
let p = Vector3::new(5.0, 3.0, 1.0);
let r = t.apply(p);
assert!((r.x - (-1.0)).abs() < 1e-10); assert!((r.y - 3.0).abs() < 1e-10);
assert!((r.z - 1.0).abs() < 1e-10);
}
#[test]
fn test_mirror_line() {
let t = Transform::from_mirror_line(
Vector3::new(0.0, 0.0, 0.0),
Vector3::new(1.0, 0.0, 0.0),
);
let p = Vector3::new(3.0, 5.0, 0.0);
let r = t.apply(p);
assert!((r.x - 3.0).abs() < 1e-10);
assert!((r.y - (-5.0)).abs() < 1e-10);
}
#[test]
fn test_mirror_line_diagonal() {
let t = Transform::from_mirror_line(
Vector3::ZERO,
Vector3::new(1.0, 1.0, 0.0),
);
let p = Vector3::new(3.0, 0.0, 0.0);
let r = t.apply(p);
assert!((r.x - 0.0).abs() < 1e-10);
assert!((r.y - 3.0).abs() < 1e-10);
}
#[test]
fn test_mirror_is_involution() {
let t = Transform::from_mirror_plane(
Vector3::new(1.0, 2.0, 3.0),
Vector3::new(1.0, 1.0, 0.0),
);
let p = Vector3::new(7.0, -2.0, 4.0);
let r = t.apply(t.apply(p));
assert!((r.x - p.x).abs() < 1e-10);
assert!((r.y - p.y).abs() < 1e-10);
assert!((r.z - p.z).abs() < 1e-10);
}
#[test]
fn test_normalize_angle() {
assert!((normalize_angle(0.0)).abs() < 1e-10);
assert!((normalize_angle(PI) - PI).abs() < 1e-10);
assert!((normalize_angle(-PI / 2.0) - 3.0 * PI / 2.0).abs() < 1e-10);
assert!((normalize_angle(3.0 * PI) - PI).abs() < 1e-10);
}
}