M-Labs/nalgebra
3
0
Fork 1
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge branch 'dev' into complex

Browse Source 
# Conflicts:
#	src/base/ops.rs
#	src/geometry/isometry.rs
#	src/geometry/quaternion.rs
#	src/geometry/quaternion_construction.rs
#	src/geometry/rotation.rs
#	src/geometry/similarity.rs
#	src/geometry/transform.rs
#	src/geometry/translation.rs
#	src/geometry/unit_complex.rs
This commit is contained in:
sebcrozet 2019-03-31 10:48:59 +02:00
commit 38ef0cbf7b
25 changed files with 1227 additions and 74 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
Cargo.toml
View File

@ -48,6 +48,9 @@ quickcheck = { version = "0.8", optional = true }
pest = { version = "2.0", optional = true }
pest_derive = { version = "2.0", optional = true }
[patch.crates-io]
alga = { git = "https://github.com/rustsim/alga", branch = "dev" }
[dev-dependencies]
serde_json = "1.0"
rand_xorshift = "0.1"

2
src/base/blas.rs
View File

@ -494,7 +494,7 @@ where
{
/// Computes `self = a * x + b * self`.
///
/// If be is zero, `self` is never read from.
/// If `b` is zero, `self` is never read from.
///
/// # Examples:
///

62
src/base/ops.rs
View File

@ -11,7 +11,7 @@ use crate::base::allocator::{Allocator, SameShapeAllocator, SameShapeC, SameShap
use crate::base::constraint::{
AreMultipliable, DimEq, SameNumberOfColumns, SameNumberOfRows, ShapeConstraint,
};
use crate::base::dimension::{Dim, DimMul, DimName, DimProd};
use crate::base::dimension::{Dim, DimMul, DimName, DimProd, Dynamic};
use crate::base::storage::{ContiguousStorageMut, Storage, StorageMut};
use crate::base::{DefaultAllocator, Matrix, MatrixMN, MatrixN, MatrixSum, Scalar, VectorSliceN};
@ -384,6 +384,36 @@ where
}
}
impl<N, C: Dim> iter::Sum for MatrixMN<N, Dynamic, C>
where
N: Scalar + ClosedAdd + Zero,
DefaultAllocator: Allocator<N, Dynamic, C>,
{
/// # Example
/// ```
/// # use nalgebra::DVector;
/// assert_eq!(vec![DVector::repeat(3, 1.0f64),
/// DVector::repeat(3, 1.0f64),
/// DVector::repeat(3, 1.0f64)].into_iter().sum::<DVector<f64>>(),
/// DVector::repeat(3, 1.0f64) + DVector::repeat(3, 1.0f64) + DVector::repeat(3, 1.0f64));
/// ```
///
/// # Panics
/// Panics if the iterator is empty:
/// ```should_panic
/// # use std::iter;
/// # use nalgebra::DMatrix;
/// iter::empty::<DMatrix<f64>>().sum::<DMatrix<f64>>(); // panics!
/// ```
fn sum<I: Iterator<Item = MatrixMN<N, Dynamic, C>>>(mut iter: I) -> MatrixMN<N, Dynamic, C> {
if let Some(first) = iter.next() {
iter.fold(first, |acc, x| acc + x)
} else {
panic!("Cannot compute `sum` of empty iterator.")
}
}
}
impl<'a, N, R: DimName, C: DimName> iter::Sum<&'a MatrixMN<N, R, C>> for MatrixMN<N, R, C>
where
N: Scalar + ClosedAdd + Zero,
@ -394,6 +424,36 @@ where
}
}
impl<'a, N, C: Dim> iter::Sum<&'a MatrixMN<N, Dynamic, C>> for MatrixMN<N, Dynamic, C>
where
N: Scalar + ClosedAdd + Zero,
DefaultAllocator: Allocator<N, Dynamic, C>,
{
/// # Example
/// ```
/// # use nalgebra::DVector;
/// let v = &DVector::repeat(3, 1.0f64);
///
/// assert_eq!(vec![v, v, v].into_iter().sum::<DVector<f64>>(),
/// v + v + v);
/// ```
///
/// # Panics
/// Panics if the iterator is empty:
/// ```should_panic
/// # use std::iter;
/// # use nalgebra::DMatrix;
/// iter::empty::<&DMatrix<f64>>().sum::<DMatrix<f64>>(); // panics!
/// ```
fn sum<I: Iterator<Item = &'a MatrixMN<N, Dynamic, C>>>(mut iter: I) -> MatrixMN<N, Dynamic, C> {
if let Some(first) = iter.next() {
iter.fold(first.clone(), |acc, x| acc + x)
} else {
panic!("Cannot compute `sum` of empty iterator.")
}
}
}
/*
*
* Multiplication

2
src/base/unit.rs
View File

@ -222,7 +222,7 @@ impl<T: Neg> Neg for Unit<T> {
#[inline]
fn neg(self) -> Self::Output {
Unit::new_unchecked(-self.value)
Self::Output::new_unchecked(-self.value)
}
}

93
src/geometry/isometry.rs Normal file → Executable file
View File

@ -17,7 +17,7 @@ use alga::linear::Rotation;
use crate::base::allocator::Allocator;
use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1};
use crate::base::storage::Owned;
use crate::base::{DefaultAllocator, MatrixN};
use crate::base::{DefaultAllocator, MatrixN, VectorN};
use crate::geometry::{Point, Translation};
/// A direct isometry, i.e., a rotation followed by a translation, aka. a rigid-body motion, aka. an element of a Special Euclidean (SE) group.
@ -254,6 +254,97 @@ where DefaultAllocator: Allocator<N, D>
pub fn append_rotation_wrt_center_mut(&mut self, r: &R) {
self.rotation = self.rotation.append_rotation(r);
}
/// Transform the given point by this isometry.
///
/// This is the same as the multiplication `self * pt`.
///
/// # Example
///
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Isometry3, Translation3, UnitQuaternion, Vector3, Point3};
/// let tra = Translation3::new(0.0, 0.0, 3.0);
/// let rot = UnitQuaternion::from_scaled_axis(Vector3::y() * f32::consts::FRAC_PI_2);
/// let iso = Isometry3::from_parts(tra, rot);
///
/// let transformed_point = iso.transform_point(&Point3::new(1.0, 2.0, 3.0));
/// assert_relative_eq!(transformed_point, Point3::new(3.0, 2.0, 2.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self * pt
}
/// Transform the given vector by this isometry, ignoring the translation
/// component of the isometry.
///
/// This is the same as the multiplication `self * v`.
///
/// # Example
///
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Isometry3, Translation3, UnitQuaternion, Vector3};
/// let tra = Translation3::new(0.0, 0.0, 3.0);
/// let rot = UnitQuaternion::from_scaled_axis(Vector3::y() * f32::consts::FRAC_PI_2);
/// let iso = Isometry3::from_parts(tra, rot);
///
/// let transformed_point = iso.transform_vector(&Vector3::new(1.0, 2.0, 3.0));
/// assert_relative_eq!(transformed_point, Vector3::new(3.0, 2.0, -1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self * v
}
/// Transform the given point by the inverse of this isometry. This may be
/// less expensive than computing the entire isometry inverse and then
/// transforming the point.
///
/// # Example
///
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Isometry3, Translation3, UnitQuaternion, Vector3, Point3};
/// let tra = Translation3::new(0.0, 0.0, 3.0);
/// let rot = UnitQuaternion::from_scaled_axis(Vector3::y() * f32::consts::FRAC_PI_2);
/// let iso = Isometry3::from_parts(tra, rot);
///
/// let transformed_point = iso.inverse_transform_point(&Point3::new(1.0, 2.0, 3.0));
/// assert_relative_eq!(transformed_point, Point3::new(0.0, 2.0, 1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self.rotation
.inverse_transform_point(&(pt - &self.translation.vector))
}
/// Transform the given vector by the inverse of this isometry, ignoring the
/// translation component of the isometry. This may be
/// less expensive than computing the entire isometry inverse and then
/// transforming the point.
///
/// # Example
///
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Isometry3, Translation3, UnitQuaternion, Vector3};
/// let tra = Translation3::new(0.0, 0.0, 3.0);
/// let rot = UnitQuaternion::from_scaled_axis(Vector3::y() * f32::consts::FRAC_PI_2);
/// let iso = Isometry3::from_parts(tra, rot);
///
/// let transformed_point = iso.inverse_transform_vector(&Vector3::new(1.0, 2.0, 3.0));
/// assert_relative_eq!(transformed_point, Vector3::new(-3.0, 2.0, 1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self.rotation.inverse_transform_vector(v)
}
}
// NOTE: we don't require `R: Rotation<...>` here because this is not useful for the implementation

9
src/geometry/isometry_alga.rs Normal file → Executable file
View File

@ -85,12 +85,12 @@ where
{
#[inline]
fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self * pt
self.transform_point(pt)
}
#[inline]
fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self * v
self.transform_vector(v)
}
}
@ -101,13 +101,12 @@ where
{
#[inline]
fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self.rotation
.inverse_transform_point(&(pt - &self.translation.vector))
self.inverse_transform_point(pt)
}
#[inline]
fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self.rotation.inverse_transform_vector(v)
self.inverse_transform_vector(v)
}
}

457
src/geometry/quaternion.rs Normal file → Executable file
View File

@ -17,9 +17,9 @@ use alga::general::RealField;
use crate::base::dimension::{U1, U3, U4};
use crate::base::storage::{CStride, RStride};
use crate::base::{Matrix3, MatrixN, MatrixSlice, MatrixSliceMut, Unit, Vector3, Vector4};
use crate::base::{Matrix3, Matrix4, MatrixSlice, MatrixSliceMut, Unit, Vector3, Vector4};
use crate::geometry::Rotation;
use crate::geometry::{Point3, Rotation};
/// A quaternion. See the type alias `UnitQuaternion = Unit<Quaternion>` for a quaternion
/// that may be used as a rotation.
@ -124,6 +124,12 @@ impl<N: RealField> Quaternion<N> {
Self::from(self.coords.normalize())
}
/// The imaginary part of this quaternion.
#[inline]
pub fn imag(&self) -> Vector3<N> {
self.coords.xyz()
}
/// The conjugate of this quaternion.
///
/// # Example
@ -135,13 +141,7 @@ impl<N: RealField> Quaternion<N> {
/// ```
#[inline]
pub fn conjugate(&self) -> Self {
let v = Vector4::new(
-self.coords[0],
-self.coords[1],
-self.coords[2],
self.coords[3],
);
Self::from(v)
Self::from_parts(self.w, -self.imag())
}
/// Inverts this quaternion if it is not zero.
@ -307,6 +307,81 @@ impl<N: RealField> Quaternion<N> {
self.coords.dot(&rhs.coords)
}
/// Calculates the inner product (also known as the dot product).
/// See "Foundations of Game Engine Development, Volume 1: Mathematics" by Lengyel
/// Formula 4.89.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let a = Quaternion::new(0.0, 2.0, 3.0, 4.0);
/// let b = Quaternion::new(0.0, 5.0, 2.0, 1.0);
/// let expected = Quaternion::new(-20.0, 0.0, 0.0, 0.0);
/// let result = a.inner(&b);
/// assert_relative_eq!(expected, result, epsilon = 1.0e-5);
#[inline]
pub fn inner(&self, other: &Self) -> Self {
(self * other + other * self).half()
}
/// Calculates the outer product (also known as the wedge product).
/// See "Foundations of Game Engine Development, Volume 1: Mathematics" by Lengyel
/// Formula 4.89.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let a = Quaternion::new(0.0, 2.0, 3.0, 4.0);
/// let b = Quaternion::new(0.0, 5.0, 2.0, 1.0);
/// let expected = Quaternion::new(0.0, -5.0, 18.0, -11.0);
/// let result = a.outer(&b);
/// assert_relative_eq!(expected, result, epsilon = 1.0e-5);
/// ```
#[inline]
pub fn outer(&self, other: &Self) -> Self {
(self * other - other * self).half()
}
/// Calculates the projection of `self` onto `other` (also known as the parallel).
/// See "Foundations of Game Engine Development, Volume 1: Mathematics" by Lengyel
/// Formula 4.94.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let a = Quaternion::new(0.0, 2.0, 3.0, 4.0);
/// let b = Quaternion::new(0.0, 5.0, 2.0, 1.0);
/// let expected = Quaternion::new(0.0, 3.333333333333333, 1.3333333333333333, 0.6666666666666666);
/// let result = a.project(&b).unwrap();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-5);
/// ```
#[inline]
pub fn project(&self, other: &Self) -> Option<Self> {
self.inner(other).right_div(other)
}
/// Calculates the rejection of `self` from `other` (also known as the perpendicular).
/// See "Foundations of Game Engine Development, Volume 1: Mathematics" by Lengyel
/// Formula 4.94.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let a = Quaternion::new(0.0, 2.0, 3.0, 4.0);
/// let b = Quaternion::new(0.0, 5.0, 2.0, 1.0);
/// let expected = Quaternion::new(0.0, -1.3333333333333333, 1.6666666666666665, 3.3333333333333335);
/// let result = a.reject(&b).unwrap();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-5);
/// ```
#[inline]
pub fn reject(&self, other: &Self) -> Option<Self> {
self.outer(other).right_div(other)
}
/// The polar decomposition of this quaternion.
///
/// Returns, from left to right: the quaternion norm, the half rotation angle, the rotation
@ -506,6 +581,274 @@ impl<N: RealField> Quaternion<N> {
pub fn normalize_mut(&mut self) -> N {
self.coords.normalize_mut()
}
/// Calculates square of a quaternion.
#[inline]
pub fn squared(&self) -> Self {
self * self
}
/// Divides quaternion into two.
#[inline]
pub fn half(&self) -> Self {
self / ::convert(2.0f64)
}
/// Calculates square root.
#[inline]
pub fn sqrt(&self) -> Self {
self.powf(::convert(0.5))
}
/// Check if the quaternion is pure.
#[inline]
pub fn is_pure(&self) -> bool {
self.w.is_zero()
}
/// Convert quaternion to pure quaternion.
#[inline]
pub fn pure(&self) -> Self {
Self::from_imag(self.imag())
}
/// Left quaternionic division.
///
/// Calculates B<sup>-1</sup> * A where A = self, B = other.
#[inline]
pub fn left_div(&self, other: &Self) -> Option<Self> {
other.try_inverse().map(|inv| inv * self)
}
/// Right quaternionic division.
///
/// Calculates A * B<sup>-1</sup> where A = self, B = other.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let a = Quaternion::new(0.0, 1.0, 2.0, 3.0);
/// let b = Quaternion::new(0.0, 5.0, 2.0, 1.0);
/// let result = a.right_div(&b).unwrap();
/// let expected = Quaternion::new(0.4, 0.13333333333333336, -0.4666666666666667, 0.26666666666666666);
/// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn right_div(&self, other: &Self) -> Option<Self> {
other.try_inverse().map(|inv| self * inv)
}
/// Calculates the quaternionic cosinus.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let expected = Quaternion::new(58.93364616794395, -34.086183690465596, -51.1292755356984, -68.17236738093119);
/// let result = input.cos();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn cos(&self) -> Self {
let z = self.imag().magnitude();
let w = -self.w.sin() * z.sinhc();
Self::from_parts(self.w.cos() * z.cosh(), self.imag() * w)
}
/// Calculates the quaternionic arccosinus.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let result = input.cos().acos();
/// assert_relative_eq!(input, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn acos(&self) -> Self {
let u = Self::from_imag(self.imag().normalize());
let identity = Self::identity();
let z = (self + (self.squared() - identity).sqrt()).ln();
-(u * z)
}
/// Calculates the quaternionic sinus.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let expected = Quaternion::new(91.78371578403467, 21.886486853029176, 32.82973027954377, 43.77297370605835);
/// let result = input.sin();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn sin(&self) -> Self {
let z = self.imag().magnitude();
let w = self.w.cos() * z.sinhc();
Self::from_parts(self.w.sin() * z.cosh(), self.imag() * w)
}
/// Calculates the quaternionic arcsinus.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let result = input.sin().asin();
/// assert_relative_eq!(input, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn asin(&self) -> Self {
let u = Self::from_imag(self.imag().normalize());
let identity = Self::identity();
let z = ((u * self) + (identity - self.squared()).sqrt()).ln();
-(u * z)
}
/// Calculates the quaternionic tangent.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let expected = Quaternion::new(0.00003821631725009489, 0.3713971716439371, 0.5570957574659058, 0.7427943432878743);
/// let result = input.tan();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn tan(&self) -> Self {
self.sin().right_div(&self.cos()).unwrap()
}
/// Calculates the quaternionic arctangent.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let result = input.tan().atan();
/// assert_relative_eq!(input, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn atan(&self) -> Self {
let u = Self::from_imag(self.imag().normalize());
let num = u + self;
let den = u - self;
let fr = num.right_div(&den).unwrap();
let ln = fr.ln();
(u.half()) * ln
}
/// Calculates the hyperbolic quaternionic sinus.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let expected = Quaternion::new(0.7323376060463428, -0.4482074499805421, -0.6723111749708133, -0.8964148999610843);
/// let result = input.sinh();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn sinh(&self) -> Self {
(self.exp() - (-self).exp()).half()
}
/// Calculates the hyperbolic quaternionic arcsinus.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let expected = Quaternion::new(2.385889902585242, 0.514052600662788, 0.7710789009941821, 1.028105201325576);
/// let result = input.asinh();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn asinh(&self) -> Self {
let identity = Self::identity();
(self + (identity + self.squared()).sqrt()).ln()
}
/// Calculates the hyperbolic quaternionic cosinus.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let expected = Quaternion::new(0.9615851176369566, -0.3413521745610167, -0.5120282618415251, -0.6827043491220334);
/// let result = input.cosh();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn cosh(&self) -> Self {
(self.exp() + (-self).exp()).half()
}
/// Calculates the hyperbolic quaternionic arccosinus.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let expected = Quaternion::new(2.4014472020074007, 0.5162761016176176, 0.7744141524264264, 1.0325522032352352);
/// let result = input.acosh();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn acosh(&self) -> Self {
let identity = Self::identity();
(self + (self + identity).sqrt() * (self - identity).sqrt()).ln()
}
/// Calculates the hyperbolic quaternionic tangent.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let expected = Quaternion::new(1.0248695360556623, -0.10229568178876419, -0.1534435226831464, -0.20459136357752844);
/// let result = input.tanh();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn tanh(&self) -> Self {
self.sinh().right_div(&self.cosh()).unwrap()
}
/// Calculates the hyperbolic quaternionic arctangent.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::Quaternion;
/// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
/// let expected = Quaternion::new(0.03230293287000163, 0.5173453683196951, 0.7760180524795426, 1.0346907366393903);
/// let result = input.atanh();
/// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
/// ```
#[inline]
pub fn atanh(&self) -> Self {
let identity = Self::identity();
((identity + self).ln() - (identity - self).ln()).half()
}
}
impl<N: RealField + AbsDiffEq<Epsilon = N>> AbsDiffEq for Quaternion<N> {
@ -847,16 +1190,12 @@ impl<N: RealField> UnitQuaternion<N> {
/// ```
#[inline]
pub fn axis_angle(&self) -> Option<(Unit<Vector3<N>>, N)> {
if let Some(axis) = self.axis() {
Some((axis, self.angle()))
} else {
None
}
self.axis().map(|axis| (axis, self.angle()))
}
/// Compute the exponential of a quaternion.
///
/// Note that this function yields a `Quaternion<N>` because it looses the unit property.
/// Note that this function yields a `Quaternion<N>` because it loses the unit property.
#[inline]
pub fn exp(&self) -> Quaternion<N> {
self.as_ref().exp()
@ -864,7 +1203,7 @@ impl<N: RealField> UnitQuaternion<N> {
/// Compute the natural logarithm of a quaternion.
///
/// Note that this function yields a `Quaternion<N>` because it looses the unit property.
/// Note that this function yields a `Quaternion<N>` because it loses the unit property.
/// The vector part of the return value corresponds to the axis-angle representation (divided
/// by 2.0) of this unit quaternion.
///
@ -879,7 +1218,7 @@ impl<N: RealField> UnitQuaternion<N> {
#[inline]
pub fn ln(&self) -> Quaternion<N> {
if let Some(v) = self.axis() {
Quaternion::from_parts(N::zero(), v.into_inner() * self.angle())
Quaternion::from_imag(v.into_inner() * self.angle())
} else {
Quaternion::zero()
}
@ -1002,9 +1341,91 @@ impl<N: RealField> UnitQuaternion<N> {
/// assert_relative_eq!(rot.to_homogeneous(), expected, epsilon = 1.0e-6);
/// ```
#[inline]
pub fn to_homogeneous(&self) -> MatrixN<N, U4> {
pub fn to_homogeneous(&self) -> Matrix4<N> {
self.to_rotation_matrix().to_homogeneous()
}
/// Rotate a point by this unit quaternion.
///
/// This is the same as the multiplication `self * pt`.
///
/// # Example
///
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{UnitQuaternion, Vector3, Point3};
/// let rot = UnitQuaternion::from_axis_angle(&Vector3::y_axis(), f32::consts::FRAC_PI_2);
/// let transformed_point = rot.transform_point(&Point3::new(1.0, 2.0, 3.0));
///
/// assert_relative_eq!(transformed_point, Point3::new(3.0, 2.0, -1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn transform_point(&self, pt: &Point3<N>) -> Point3<N> {
self * pt
}
/// Rotate a vector by this unit quaternion.
///
/// This is the same as the multiplication `self * v`.
///
/// # Example
///
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{UnitQuaternion, Vector3};
/// let rot = UnitQuaternion::from_axis_angle(&Vector3::y_axis(), f32::consts::FRAC_PI_2);
/// let transformed_vector = rot.transform_vector(&Vector3::new(1.0, 2.0, 3.0));
///
/// assert_relative_eq!(transformed_vector, Vector3::new(3.0, 2.0, -1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn transform_vector(&self, v: &Vector3<N>) -> Vector3<N> {
self * v
}
/// Rotate a point by the inverse of this unit quaternion. This may be
/// cheaper than inverting the unit quaternion and transforming the
/// point.
///
/// # Example
///
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{UnitQuaternion, Vector3, Point3};
/// let rot = UnitQuaternion::from_axis_angle(&Vector3::y_axis(), f32::consts::FRAC_PI_2);
/// let transformed_point = rot.inverse_transform_point(&Point3::new(1.0, 2.0, 3.0));
///
/// assert_relative_eq!(transformed_point, Point3::new(-3.0, 2.0, 1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn inverse_transform_point(&self, pt: &Point3<N>) -> Point3<N> {
// FIXME: would it be useful performancewise not to call inverse explicitly (i-e. implement
// the inverse transformation explicitly here) ?
self.inverse() * pt
}
/// Rotate a vector by the inverse of this unit quaternion. This may be
/// cheaper than inverting the unit quaternion and transforming the
/// vector.
///
/// # Example
///
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{UnitQuaternion, Vector3};
/// let rot = UnitQuaternion::from_axis_angle(&Vector3::y_axis(), f32::consts::FRAC_PI_2);
/// let transformed_vector = rot.inverse_transform_vector(&Vector3::new(1.0, 2.0, 3.0));
///
/// assert_relative_eq!(transformed_vector, Vector3::new(-3.0, 2.0, 1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn inverse_transform_vector(&self, v: &Vector3<N>) -> Vector3<N> {
self.inverse() * v
}
}
impl<N: RealField + fmt::Display> fmt::Display for UnitQuaternion<N> {

10
src/geometry/quaternion_alga.rs Normal file → Executable file
View File

@ -200,26 +200,24 @@ impl_structures!(
impl<N: RealField> Transformation<Point3<N>> for UnitQuaternion<N> {
#[inline]
fn transform_point(&self, pt: &Point3<N>) -> Point3<N> {
self * pt
self.transform_point(pt)
}
#[inline]
fn transform_vector(&self, v: &Vector3<N>) -> Vector3<N> {
self * v
self.transform_vector(v)
}
}
impl<N: RealField> ProjectiveTransformation<Point3<N>> for UnitQuaternion<N> {
#[inline]
fn inverse_transform_point(&self, pt: &Point3<N>) -> Point3<N> {
// FIXME: would it be useful performancewise not to call inverse explicitly (i-e. implement
// the inverse transformation explicitly here) ?
self.inverse() * pt
self.inverse_transform_point(pt)
}
#[inline]
fn inverse_transform_vector(&self, v: &Vector3<N>) -> Vector3<N> {
self.inverse() * v
self.inverse_transform_vector(v)
}
}

27
src/geometry/quaternion_construction.rs
View File

@ -13,9 +13,7 @@ use alga::general::RealField;
use crate::base::dimension::U3;
use crate::base::storage::Storage;
#[cfg(feature = "arbitrary")]
use crate::base::Vector3;
use crate::base::{Unit, Vector, Vector4, Matrix3};
use crate::base::{Unit, Vector, Vector3, Vector4, Matrix3};
use crate::geometry::{Quaternion, Rotation3, UnitQuaternion};
@ -43,8 +41,13 @@ impl<N: RealField> Quaternion<N> {
/// ```
#[inline]
pub fn new(w: N, i: N, j: N, k: N) -> Self {
let v = Vector4::<N>::new(i, j, k, w);
Self::from(v)
Self::from(Vector4::new(i, j, k, w))
}
/// Constructs a pure quaternion.
#[inline]
pub fn from_imag(vector: Vector3<N>) -> Self {
Self::from_parts(N::zero(), vector)
}
/// Creates a new quaternion from its scalar and vector parts. Note that the arguments order does
@ -68,6 +71,12 @@ impl<N: RealField> Quaternion<N> {
Self::new(scalar, vector[0], vector[1], vector[2])
}
/// Constructs a real quaternion.
#[inline]
pub fn from_real(r: N) -> Self {
Self::from_parts(r, Vector3::zero())
}
/// Creates a new quaternion from its polar decomposition.
///
/// Note that `axis` is assumed to be a unit vector.
@ -92,7 +101,7 @@ impl<N: RealField> Quaternion<N> {
/// ```
#[inline]
pub fn identity() -> Self {
Self::new(N::one(), N::zero(), N::zero(), N::zero())
Self::from_real(N::one())
}
}
@ -106,7 +115,7 @@ impl<N: RealField> One for Quaternion<N> {
impl<N: RealField> Zero for Quaternion<N> {
#[inline]
fn zero() -> Self {
Self::new(N::zero(), N::zero(), N::zero(), N::zero())
Self::from(Vector4::zero())
}
#[inline]
@ -579,7 +588,7 @@ impl<N: RealField> UnitQuaternion<N> {
pub fn new<SB>(axisangle: Vector<N, U3, SB>) -> Self
where SB: Storage<N, U3> {
let two: N = crate::convert(2.0f64);
let q = Quaternion::<N>::from_parts(N::zero(), axisangle / two).exp();
let q = Quaternion::<N>::from_imag(axisangle / two).exp();
Self::new_unchecked(q)
}
@ -608,7 +617,7 @@ impl<N: RealField> UnitQuaternion<N> {
pub fn new_eps<SB>(axisangle: Vector<N, U3, SB>, eps: N) -> Self
where SB: Storage<N, U3> {
let two: N = crate::convert(2.0f64);
let q = Quaternion::<N>::from_parts(N::zero(), axisangle / two).exp_eps(eps);
let q = Quaternion::<N>::from_imag(axisangle / two).exp_eps(eps);
Self::new_unchecked(q)
}

4
src/geometry/quaternion_ops.rs
View File

@ -552,7 +552,7 @@ impl<N: RealField> Neg for Quaternion<N> {
#[inline]
fn neg(self) -> Self::Output {
Quaternion::from(-self.coords)
Self::Output::from(-self.coords)
}
}
@ -561,7 +561,7 @@ impl<'a, N: RealField> Neg for &'a Quaternion<N> {
#[inline]
fn neg(self) -> Self::Output {
Quaternion::from(-&self.coords)
Self::Output::from(-&self.coords)
}
}

83
src/geometry/rotation.rs Normal file → Executable file
View File

@ -18,7 +18,8 @@ use alga::general::RealField;
use crate::base::allocator::Allocator;
use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1};
use crate::base::{DefaultAllocator, MatrixN, Scalar};
use crate::base::{DefaultAllocator, MatrixN, Scalar, VectorN};
use crate::geometry::Point;
/// A rotation matrix.
#[repr(C)]
@ -351,6 +352,86 @@ where DefaultAllocator: Allocator<N, D, D>
}
}
impl<N: Real, D: DimName> Rotation<N, D>
where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>
{
/// Rotate the given point.
///
/// This is the same as the multiplication `self * pt`.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Point3, Rotation2, Rotation3, UnitQuaternion, Vector3};
/// let rot = Rotation3::new(Vector3::y() * f32::consts::FRAC_PI_2);
/// let transformed_point = rot.transform_point(&Point3::new(1.0, 2.0, 3.0));
///
/// assert_relative_eq!(transformed_point, Point3::new(3.0, 2.0, -1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self * pt
}
/// Rotate the given vector.
///
/// This is the same as the multiplication `self * v`.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Rotation2, Rotation3, UnitQuaternion, Vector3};
/// let rot = Rotation3::new(Vector3::y() * f32::consts::FRAC_PI_2);
/// let transformed_vector = rot.transform_vector(&Vector3::new(1.0, 2.0, 3.0));
///
/// assert_relative_eq!(transformed_vector, Vector3::new(3.0, 2.0, -1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self * v
}
/// Rotate the given point by the inverse of this rotation. This may be
/// cheaper than inverting the rotation and then transforming the given
/// point.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Point3, Rotation2, Rotation3, UnitQuaternion, Vector3};
/// let rot = Rotation3::new(Vector3::y() * f32::consts::FRAC_PI_2);
/// let transformed_point = rot.inverse_transform_point(&Point3::new(1.0, 2.0, 3.0));
///
/// assert_relative_eq!(transformed_point, Point3::new(-3.0, 2.0, 1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
Point::from(self.inverse_transform_vector(&pt.coords))
}
/// Rotate the given vector by the inverse of this rotation. This may be
/// cheaper than inverting the rotation and then transforming the given
/// vector.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Rotation2, Rotation3, UnitQuaternion, Vector3};
/// let rot = Rotation3::new(Vector3::y() * f32::consts::FRAC_PI_2);
/// let transformed_vector = rot.inverse_transform_vector(&Vector3::new(1.0, 2.0, 3.0));
///
/// assert_relative_eq!(transformed_vector, Vector3::new(-3.0, 2.0, 1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self.matrix().tr_mul(v)
}
}
impl<N: Scalar + Eq, D: DimName> Eq for Rotation<N, D> where DefaultAllocator: Allocator<N, D, D> {}
impl<N: Scalar + PartialEq, D: DimName> PartialEq for Rotation<N, D>

8
src/geometry/rotation_alga.rs Normal file → Executable file
View File

@ -75,12 +75,12 @@ where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>
{
#[inline]
fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self * pt
self.transform_point(pt)
}
#[inline]
fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self * v
self.transform_vector(v)
}
}
@ -89,12 +89,12 @@ where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>
{
#[inline]
fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
Point::from(self.inverse_transform_vector(&pt.coords))
self.inverse_transform_point(pt)
}
#[inline]
fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self.matrix().tr_mul(v)
self.inverse_transform_vector(v)
}
}

83
src/geometry/similarity.rs Normal file → Executable file
View File

@ -16,7 +16,7 @@ use alga::linear::Rotation;
use crate::base::allocator::Allocator;
use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1};
use crate::base::storage::Owned;
use crate::base::{DefaultAllocator, MatrixN};
use crate::base::{DefaultAllocator, MatrixN, VectorN};
use crate::geometry::{Isometry, Point, Translation};
/// A similarity, i.e., an uniform scaling, followed by a rotation, followed by a translation.
@ -238,6 +238,87 @@ where
pub fn append_rotation_wrt_center_mut(&mut self, r: &R) {
self.isometry.append_rotation_wrt_center_mut(r)
}
/// Transform the given point by this similarity.
///
/// This is the same as the multiplication `self * pt`.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Point3, Similarity3, Vector3};
/// let axisangle = Vector3::y() * f32::consts::FRAC_PI_2;
/// let translation = Vector3::new(1.0, 2.0, 3.0);
/// let sim = Similarity3::new(translation, axisangle, 3.0);
/// let transformed_point = sim.transform_point(&Point3::new(4.0, 5.0, 6.0));
/// assert_relative_eq!(transformed_point, Point3::new(19.0, 17.0, -9.0), epsilon = 1.0e-5);
/// ```
#[inline]
pub fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self * pt
}
/// Transform the given vector by this similarity, ignoring the translational
/// component.
///
/// This is the same as the multiplication `self * t`.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Similarity3, Vector3};
/// let axisangle = Vector3::y() * f32::consts::FRAC_PI_2;
/// let translation = Vector3::new(1.0, 2.0, 3.0);
/// let sim = Similarity3::new(translation, axisangle, 3.0);
/// let transformed_vector = sim.transform_vector(&Vector3::new(4.0, 5.0, 6.0));
/// assert_relative_eq!(transformed_vector, Vector3::new(18.0, 15.0, -12.0), epsilon = 1.0e-5);
/// ```
#[inline]
pub fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self * v
}
/// Transform the given point by the inverse of this similarity. This may
/// be cheaper than inverting the similarity and then transforming the
/// given point.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Point3, Similarity3, Vector3};
/// let axisangle = Vector3::y() * f32::consts::FRAC_PI_2;
/// let translation = Vector3::new(1.0, 2.0, 3.0);
/// let sim = Similarity3::new(translation, axisangle, 2.0);
/// let transformed_point = sim.inverse_transform_point(&Point3::new(4.0, 5.0, 6.0));
/// assert_relative_eq!(transformed_point, Point3::new(-1.5, 1.5, 1.5), epsilon = 1.0e-5);
/// ```
#[inline]
pub fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self.isometry.inverse_transform_point(pt) / self.scaling()
}
/// Transform the given vector by the inverse of this similarity,
/// ignoring the translational component. This may be cheaper than
/// inverting the similarity and then transforming the given vector.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use std::f32;
/// # use nalgebra::{Similarity3, Vector3};
/// let axisangle = Vector3::y() * f32::consts::FRAC_PI_2;
/// let translation = Vector3::new(1.0, 2.0, 3.0);
/// let sim = Similarity3::new(translation, axisangle, 2.0);
/// let transformed_vector = sim.inverse_transform_vector(&Vector3::new(4.0, 5.0, 6.0));
/// assert_relative_eq!(transformed_vector, Vector3::new(-3.0, 2.5, 2.0), epsilon = 1.0e-5);
/// ```
#[inline]
pub fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self.isometry.inverse_transform_vector(v) / self.scaling()
}
}
// NOTE: we don't require `R: Rotation<...>` here because this is not useful for the implementation

8
src/geometry/similarity_alga.rs Normal file → Executable file
View File

@ -82,12 +82,12 @@ where
{
#[inline]
fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self * pt
self.transform_point(pt)
}
#[inline]
fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self * v
self.transform_vector(v)
}
}
@ -98,12 +98,12 @@ where
{
#[inline]
fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self.isometry.inverse_transform_point(pt) / self.scaling()
self.inverse_transform_point(pt)
}
#[inline]
fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self.isometry.inverse_transform_vector(v) / self.scaling()
self.inverse_transform_vector(v)
}
}

57
src/geometry/transform.rs Normal file → Executable file
View File

@ -6,12 +6,14 @@ use std::marker::PhantomData;
#[cfg(feature = "serde-serialize")]
use serde::{Deserialize, Deserializer, Serialize, Serializer};
use alga::general::RealField;
use alga::general::{RealField, TwoSidedInverse};
use crate::base::allocator::Allocator;
use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1};
use crate::base::storage::Owned;
use crate::base::{DefaultAllocator, MatrixN};
use crate::base::{DefaultAllocator, MatrixN, VectorN};
use crate::geometry::Point;
/// Trait implemented by phantom types identifying the projective transformation type.
///
@ -452,6 +454,57 @@ where DefaultAllocator: Allocator<N, DimNameSum<D, U1>, DimNameSum<D, U1>>
}
}
impl<N, D: DimNameAdd<U1>, C> Transform<N, D, C>
where
N: RealField,
C: TCategory,
DefaultAllocator: Allocator<N, DimNameSum<D, U1>, DimNameSum<D, U1>>
+ Allocator<N, DimNameSum<D, U1>>
+ Allocator<N, D, D>
+ Allocator<N, D>,
{
/// Transform the given point by this transformation.
///
/// This is the same as the multiplication `self * pt`.
#[inline]
pub fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self * pt
}
/// Transform the given vector by this transformation, ignoring the
/// translational component of the transformation.
///
/// This is the same as the multiplication `self * v`.
#[inline]
pub fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self * v
}
}
impl<N: RealField, D: DimNameAdd<U1>, C: TCategory> Transform<N, D, C>
where C: SubTCategoryOf<TProjective>,
DefaultAllocator: Allocator<N, DimNameSum<D, U1>, DimNameSum<D, U1>>
+ Allocator<N, DimNameSum<D, U1>>
+ Allocator<N, D, D>
+ Allocator<N, D>,
{
/// Transform the given point by the inverse of this transformation.
/// This may be cheaper than inverting the transformation and transforming
/// the point.
#[inline]
pub fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self.two_sided_inverse() * pt
}
/// Transform the given vector by the inverse of this transformation.
/// This may be cheaper than inverting the transformation and transforming
/// the vector.
#[inline]
pub fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self.two_sided_inverse() * v
}
}
impl<N: RealField, D: DimNameAdd<U1>> Transform<N, D, TGeneral>
where DefaultAllocator: Allocator<N, DimNameSum<D, U1>, DimNameSum<D, U1>>
{

8
src/geometry/transform_alga.rs Normal file → Executable file
View File

@ -96,12 +96,12 @@ where
{
#[inline]
fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self * pt
self.transform_point(pt)
}
#[inline]
fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self * v
self.transform_vector(v)
}
}
@ -116,12 +116,12 @@ where
{
#[inline]
fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
self.two_sided_inverse() * pt
self.inverse_transform_point(pt)
}
#[inline]
fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
self.two_sided_inverse() * v
self.inverse_transform_vector(v)
}
}

40
src/geometry/translation.rs Normal file → Executable file
View File

@ -11,13 +11,15 @@ use serde::{Deserialize, Deserializer, Serialize, Serializer};
#[cfg(feature = "abomonation-serialize")]
use abomonation::Abomonation;
use alga::general::{ClosedNeg, RealField};
use alga::general::{ClosedAdd, ClosedNeg, ClosedSub, RealField};
use crate::base::allocator::Allocator;
use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1};
use crate::base::storage::Owned;
use crate::base::{DefaultAllocator, MatrixN, Scalar, VectorN};
use crate::geometry::Point;
/// A translation.
#[repr(C)]
#[derive(Debug)]
@ -190,6 +192,42 @@ where DefaultAllocator: Allocator<N, D>
}
}
impl<N: Scalar + ClosedAdd, D: DimName> Translation<N, D>
where DefaultAllocator: Allocator<N, D>
{
/// Translate the given point.
///
/// This is the same as the multiplication `self * pt`.
///
/// # Example
/// ```
/// # use nalgebra::{Translation3, Point3};
/// let t = Translation3::new(1.0, 2.0, 3.0);
/// let transformed_point = t.transform_point(&Point3::new(4.0, 5.0, 6.0));
/// assert_eq!(transformed_point, Point3::new(5.0, 7.0, 9.0));
#[inline]
pub fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
pt + &self.vector
}
}
impl<N: Scalar + ClosedSub, D: DimName> Translation<N, D>
where DefaultAllocator: Allocator<N, D>
{
/// Translate the given point by the inverse of this translation.
///
/// # Example
/// ```
/// # use nalgebra::{Translation3, Point3};
/// let t = Translation3::new(1.0, 2.0, 3.0);
/// let transformed_point = t.inverse_transform_point(&Point3::new(4.0, 5.0, 6.0));
/// assert_eq!(transformed_point, Point3::new(3.0, 3.0, 3.0));
#[inline]
pub fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
pt - &self.vector
}
}
impl<N: Scalar + Eq, D: DimName> Eq for Translation<N, D> where DefaultAllocator: Allocator<N, D> {}
impl<N: Scalar + PartialEq, D: DimName> PartialEq for Translation<N, D>

4
src/geometry/translation_alga.rs Normal file → Executable file
View File

@ -76,7 +76,7 @@ where DefaultAllocator: Allocator<N, D>
{
#[inline]
fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
pt + &self.vector
self.transform_point(pt)
}
#[inline]
@ -90,7 +90,7 @@ where DefaultAllocator: Allocator<N, D>
{
#[inline]
fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
pt - &self.vector
self.inverse_transform_point(pt)
}
#[inline]

74
src/geometry/unit_complex.rs Normal file → Executable file
View File

@ -3,8 +3,8 @@ use num_complex::Complex;
use std::fmt;
use alga::general::RealField;
use crate::base::{Matrix2, Matrix3, Unit, Vector1};
use crate::geometry::Rotation2;
use crate::base::{Matrix2, Matrix3, Unit, Vector1, Vector2};
use crate::geometry::{Rotation2, Point2};
/// A complex number with a norm equal to 1.
pub type UnitComplex<N> = Unit<Complex<N>>;
@ -251,6 +251,76 @@ impl<N: RealField> UnitComplex<N> {
pub fn to_homogeneous(&self) -> Matrix3<N> {
self.to_rotation_matrix().to_homogeneous()
}
/// Rotate the given point by this unit complex number.
///
/// This is the same as the multiplication `self * pt`.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::{UnitComplex, Point2};
/// # use std::f32;
/// let rot = UnitComplex::new(f32::consts::FRAC_PI_2);
/// let transformed_point = rot.transform_point(&Point2::new(1.0, 2.0));
/// assert_relative_eq!(transformed_point, Point2::new(-2.0, 1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn transform_point(&self, pt: &Point2<N>) -> Point2<N> {
self * pt
}
/// Rotate the given vector by this unit complex number.
///
/// This is the same as the multiplication `self * v`.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::{UnitComplex, Vector2};
/// # use std::f32;
/// let rot = UnitComplex::new(f32::consts::FRAC_PI_2);
/// let transformed_vector = rot.transform_vector(&Vector2::new(1.0, 2.0));
/// assert_relative_eq!(transformed_vector, Vector2::new(-2.0, 1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn transform_vector(&self, v: &Vector2<N>) -> Vector2<N> {
self * v
}
/// Rotate the given point by the inverse of this unit complex number.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::{UnitComplex, Point2};
/// # use std::f32;
/// let rot = UnitComplex::new(f32::consts::FRAC_PI_2);
/// let transformed_point = rot.inverse_transform_point(&Point2::new(1.0, 2.0));
/// assert_relative_eq!(transformed_point, Point2::new(2.0, -1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn inverse_transform_point(&self, pt: &Point2<N>) -> Point2<N> {
// FIXME: would it be useful performancewise not to call inverse explicitly (i-e. implement
// the inverse transformation explicitly here) ?
self.inverse() * pt
}
/// Rotate the given vector by the inverse of this unit complex number.
///
/// # Example
/// ```
/// # #[macro_use] extern crate approx;
/// # use nalgebra::{UnitComplex, Vector2};
/// # use std::f32;
/// let rot = UnitComplex::new(f32::consts::FRAC_PI_2);
/// let transformed_vector = rot.inverse_transform_vector(&Vector2::new(1.0, 2.0));
/// assert_relative_eq!(transformed_vector, Vector2::new(2.0, -1.0), epsilon = 1.0e-6);
/// ```
#[inline]
pub fn inverse_transform_vector(&self, v: &Vector2<N>) -> Vector2<N> {
self.inverse() * v
}
}
impl<N: RealField + fmt::Display> fmt::Display for UnitComplex<N> {

10
src/geometry/unit_complex_alga.rs Normal file → Executable file
View File

@ -63,12 +63,12 @@ where DefaultAllocator: Allocator<N, U2>
{
#[inline]
fn transform_point(&self, pt: &Point2<N>) -> Point2<N> {
self * pt
self.transform_point(pt)
}
#[inline]
fn transform_vector(&self, v: &Vector2<N>) -> Vector2<N> {
self * v
self.transform_vector(v)
}
}
@ -77,14 +77,12 @@ where DefaultAllocator: Allocator<N, U2>
{
#[inline]
fn inverse_transform_point(&self, pt: &Point2<N>) -> Point2<N> {
// FIXME: would it be useful performancewise not to call inverse explicitly (i-e. implement
// the inverse transformation explicitly here) ?
self.inverse() * pt
self.inverse_transform_point(pt)
}
#[inline]
fn inverse_transform_vector(&self, v: &Vector2<N>) -> Vector2<N> {
self.inverse() * v
self.inverse_transform_vector(v)
}
}

129
src/linalg/convolution.rs Normal file
View File

@ -0,0 +1,129 @@
use base::allocator::Allocator;
use base::default_allocator::DefaultAllocator;
use base::dimension::{Dim, DimAdd, DimDiff, DimMax, DimMaximum, DimSub, DimSum};
use std::cmp;
use storage::Storage;
use {zero, Real, Vector, VectorN, U1};
impl<N: Real, D1: Dim, S1: Storage<N, D1>> Vector<N, D1, S1> {
/// Returns the convolution of the target vector and a kernel
///
/// # Arguments
///
/// * `kernel` - A Vector with size > 0
///
/// # Errors
/// Inputs must statisfy `vector.len() >= kernel.len() > 0`.
///
pub fn convolve_full<D2, S2>(
&self,
kernel: Vector<N, D2, S2>,
) -> VectorN<N, DimDiff<DimSum<D1, D2>, U1>>
where
D1: DimAdd<D2>,
D2: DimAdd<D1, Output = DimSum<D1, D2>>,
DimSum<D1, D2>: DimSub<U1>,
S2: Storage<N, D2>,
DefaultAllocator: Allocator<N, DimDiff<DimSum<D1, D2>, U1>>,
{
let vec = self.len();
let ker = kernel.len();
if ker == 0 || ker > vec {
panic!("convolve_full expects `self.len() >= kernel.len() > 0`, received {} and {} respectively.",vec,ker);
}
let result_len = self.data.shape().0.add(kernel.data.shape().0).sub(U1);
let mut conv = VectorN::zeros_generic(result_len, U1);
for i in 0..(vec + ker - 1) {
let u_i = if i > vec { i - ker } else { 0 };
let u_f = cmp::min(i, vec - 1);
if u_i == u_f {
conv[i] += self[u_i] * kernel[(i - u_i)];
} else {
for u in u_i..(u_f + 1) {
if i - u < ker {
conv[i] += self[u] * kernel[(i - u)];
}
}
}
}
conv
}
/// Returns the convolution of the target vector and a kernel
/// The output convolution consists only of those elements that do not rely on the zero-padding.
/// # Arguments
///
/// * `kernel` - A Vector with size > 0
///
///
/// # Errors
/// Inputs must statisfy `self.len() >= kernel.len() > 0`.
///
pub fn convolve_valid<D2, S2>(&self, kernel: Vector<N, D2, S2>,
) -> VectorN<N, DimDiff<DimSum<D1, U1>, D2>>
where
D1: DimAdd<U1>,
D2: Dim,
DimSum<D1, U1>: DimSub<D2>,
S2: Storage<N, D2>,
DefaultAllocator: Allocator<N, DimDiff<DimSum<D1, U1>, D2>>,
{
let vec = self.len();
let ker = kernel.len();
if ker == 0 || ker > vec {
panic!("convolve_valid expects `self.len() >= kernel.len() > 0`, received {} and {} respectively.",vec,ker);
}
let result_len = self.data.shape().0.add(U1).sub(kernel.data.shape().0);
let mut conv = VectorN::zeros_generic(result_len, U1);
for i in 0..(vec - ker + 1) {
for j in 0..ker {
conv[i] += self[i + j] * kernel[ker - j - 1];
}
}
conv
}
/// Returns the convolution of the targetvector and a kernel
/// The output convolution is the same size as vector, centered with respect to the ‘full’ output.
/// # Arguments
///
/// * `kernel` - A Vector with size > 0
///
/// # Errors
/// Inputs must statisfy `self.len() >= kernel.len() > 0`.
pub fn convolve_same<D2, S2>(&self, kernel: Vector<N, D2, S2>) -> VectorN<N, DimMaximum<D1, D2>>
where
D1: DimMax<D2>,
D2: DimMax<D1, Output = DimMaximum<D1, D2>>,
S2: Storage<N, D2>,
DefaultAllocator: Allocator<N, DimMaximum<D1, D2>>,
{
let vec = self.len();
let ker = kernel.len();
if ker == 0 || ker > vec {
panic!("convolve_same expects `self.len() >= kernel.len() > 0`, received {} and {} respectively.",vec,ker);
}
let result_len = self.data.shape().0.max(kernel.data.shape().0);
let mut conv = VectorN::zeros_generic(result_len, U1);
for i in 0..vec {
for j in 0..ker {
let val = if i + j < 1 || i + j >= vec + 1 {
zero::<N>()
} else {
self[i + j - 1]
};
conv[i] += val * kernel[ker - j - 1];
}
}
conv
}
}

6
src/linalg/inverse.rs
View File

@ -36,7 +36,7 @@ impl<N: ComplexField, D: Dim, S: StorageMut<N, D, D>> SquareMatrix<N, D, S> {
0 => true,
1 => {
let determinant = self.get_unchecked((0, 0)).clone();
if determinant == N::zero() {
if determinant.is_zero() {
false
} else {
*self.get_unchecked_mut((0, 0)) = N::one() / determinant;
@ -51,7 +51,7 @@ impl<N: ComplexField, D: Dim, S: StorageMut<N, D, D>> SquareMatrix<N, D, S> {
let determinant = m11 * m22 - m21 * m12;
if determinant == N::zero() {
if determinant.is_zero() {
false
} else {
*self.get_unchecked_mut((0, 0)) = m22 / determinant;
@ -83,7 +83,7 @@ impl<N: ComplexField, D: Dim, S: StorageMut<N, D, D>> SquareMatrix<N, D, S> {
let determinant =
m11 * minor_m12_m23 - m12 * minor_m11_m23 + m13 * minor_m11_m22;
if determinant == N::zero() {
if determinant.is_zero() {
false
} else {
*self.get_unchecked_mut((0, 0)) = minor_m12_m23 / determinant;

2
src/linalg/mod.rs
View File

@ -17,6 +17,7 @@ mod solve;
mod svd;
mod symmetric_eigen;
mod symmetric_tridiagonal;
mod convolution;
//// FIXME: Not complete enough for publishing.
//// This handles only cases where each eigenvalue has multiplicity one.
@ -33,3 +34,4 @@ pub use self::schur::*;
pub use self::svd::*;
pub use self::symmetric_eigen::*;
pub use self::symmetric_tridiagonal::*;
pub use self::convolution::*;

119
tests/linalg/convolution.rs Normal file
View File

@ -0,0 +1,119 @@
use na::{Vector2,Vector3,Vector4,Vector5,DVector};
use std::panic;
//
// Should mimic calculations in Python's scipy library
// >>>from scipy.signal import convolve
//
// >>> convolve([1,2,3,4],[1,2],"same")
// array([ 1, 4, 7, 10])
#[test]
fn convolve_same_check(){
// Static Tests
let actual_s = Vector4::from_vec(vec![1.0,4.0,7.0,10.0]);
let expected_s = Vector4::new(1.0,2.0,3.0,4.0).convolve_same(Vector2::new(1.0,2.0));
assert!(relative_eq!(actual_s, expected_s, epsilon = 1.0e-7));
// Dynamic Tests
let actual_d = DVector::from_vec(vec![1.0,4.0,7.0,10.0]);
let expected_d = DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_same(DVector::from_vec(vec![1.0,2.0]));
assert!(relative_eq!(actual_d, expected_d, epsilon = 1.0e-7));
// Panic Tests
// These really only apply to dynamic sized vectors
assert!(
panic::catch_unwind(|| {
DVector::from_vec(vec![1.0,2.0]).convolve_same(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
}).is_err()
);
assert!(
panic::catch_unwind(|| {
DVector::<f32>::from_vec(vec![]).convolve_same(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
}).is_err()
);
assert!(
panic::catch_unwind(|| {
DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_same(DVector::<f32>::from_vec(vec![]));
}).is_err()
);
}
// >>> convolve([1,2,3,4],[1,2],"full")
// array([ 1, 4, 7, 10, 8])
#[test]
fn convolve_full_check(){
// Static Tests
let actual_s = Vector5::new(1.0,4.0,7.0,10.0,8.0);
let expected_s = Vector4::new(1.0,2.0,3.0,4.0).convolve_full(Vector2::new(1.0,2.0));
assert!(relative_eq!(actual_s, expected_s, epsilon = 1.0e-7));
// Dynamic Tests
let actual_d = DVector::from_vec(vec![1.0,4.0,7.0,10.0,8.0]);
let expected_d = DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_full(DVector::from_vec(vec![1.0,2.0]));
assert!(relative_eq!(actual_d, expected_d, epsilon = 1.0e-7));
// Panic Tests
// These really only apply to dynamic sized vectors
assert!(
panic::catch_unwind(|| {
DVector::from_vec(vec![1.0,2.0]).convolve_full(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
}).is_err()
);
assert!(
panic::catch_unwind(|| {
DVector::<f32>::from_vec(vec![]).convolve_full(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
}).is_err()
);
assert!(
panic::catch_unwind(|| {
DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_full(DVector::<f32>::from_vec(vec![]));
}).is_err()
);
}
// >>> convolve([1,2,3,4],[1,2],"valid")
// array([ 4, 7, 10])
#[test]
fn convolve_valid_check(){
// Static Tests
let actual_s = Vector3::from_vec(vec![4.0,7.0,10.0]);
let expected_s = Vector4::new(1.0,2.0,3.0,4.0).convolve_valid( Vector2::new(1.0,2.0));
assert!(relative_eq!(actual_s, expected_s, epsilon = 1.0e-7));
// Dynamic Tests
let actual_d = DVector::from_vec(vec![4.0,7.0,10.0]);
let expected_d = DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_valid(DVector::from_vec(vec![1.0,2.0]));
assert!(relative_eq!(actual_d, expected_d, epsilon = 1.0e-7));
// Panic Tests
// These really only apply to dynamic sized vectors
assert!(
panic::catch_unwind(|| {
DVector::from_vec(vec![1.0,2.0]).convolve_valid(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
}).is_err()
);
assert!(
panic::catch_unwind(|| {
DVector::<f32>::from_vec(vec![]).convolve_valid(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
}).is_err()
);
assert!(
panic::catch_unwind(|| {
DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_valid(DVector::<f32>::from_vec(vec![]));
}).is_err()
);
}

1
tests/linalg/mod.rs
View File

@ -11,3 +11,4 @@ mod schur;
mod solve;
mod svd;
mod tridiagonal;
mod convolution;