1use ndarray::Dimension;
2
3use crate::*;
4
5fn minor
6<
7 E : MatEl + nd::NdFloat,
8 Descriptor : mat::Descriptor
9>
10(
11 from : &Mat4< E, Descriptor >,
12 to : &mut Mat3< E, Descriptor >,
13 i : usize,
14 j : usize
15)
16where
17Mat4< E, Descriptor > : RawSliceMut< Scalar = E > + IndexingRef< Scalar = E, Index = Ix2 >,
18Mat3< E, Descriptor > : RawSliceMut< Scalar = E >
19{
20 for( id, ( _, v ) ) in from
21 .iter_indexed_unstable()
22 .filter(
23 | ( id, _ ) |
24 {
25 let ( r, c ) = id.into_pattern();
26 r != i && c != j
27 }
28 ).enumerate()
29 {
30 to.raw_slice_mut()[ id ] = *v;
31 }
32}
33
34fn cofactor
35<
36 E : MatEl + nd::NdFloat,
37 Descriptor : mat::Descriptor
38>
39(
40 from : &Mat4< E, Descriptor >,
41 to : &mut Mat3< E, Descriptor >,
42 i : usize,
43 j : usize
44) -> E
45where
46Mat4< E, Descriptor > :
47 RawSliceMut< Scalar = E > +
48 IndexingRef< Scalar = E, Index = Ix2 >,
49Mat3< E, Descriptor > :
50 RawSliceMut< Scalar = E > +
51 ScalarRef< Scalar = E, Index = Ix2 > +
52 ConstLayout< Index = Ix2 > +
53 IndexingMut< Scalar = E, Index = Ix2 >
54{
55 let k = E::from( ( -1i32 ).pow( ( i + j ) as u32 ) ).unwrap();
56 minor( from, to, i, j );
57 k * to.determinant()
58}
59
60impl< E, Descriptor > Mat< 4, 4, E, Descriptor >
61where
62E : MatEl + nd::NdFloat,
63Descriptor : mat::Descriptor,
64Self : ScalarMut< Scalar = E, Index = Ix2 > +
65 RawSliceMut< Scalar = E > +
66 ConstLayout< Index = Ix2 > +
67 IndexingMut< Scalar = E, Index = Ix2 >
68{
69 pub fn to_array( &self ) -> [ E; 16 ]
71 {
72 self.raw_slice().try_into().unwrap()
73 }
74
75
76 pub fn determinant( &self ) -> E
78 where
79 Mat< 3, 3, E, Descriptor > :
80 RawSliceMut< Scalar = E > +
81 ScalarMut< Scalar = E, Index = Ix2 > +
82 ConstLayout< Index = Ix2 > +
83 IndexingMut< Scalar = E, Index = Ix2 >
84 {
85 let _a11 = *self.scalar_ref( Ix2( 0, 0 ) );
86 let _a12 = *self.scalar_ref( Ix2( 0, 1 ) );
87 let _a13 = *self.scalar_ref( Ix2( 0, 2 ) );
88 let _a14 = *self.scalar_ref( Ix2( 0, 3 ) );
89
90 let mut m = Mat3::< E, Descriptor >::default();
91
92 minor( self, &mut m, 0, 0 );
93 let _det11 = m.determinant();
94 minor( self, &mut m, 0, 1 );
95 let _det12 = m.determinant();
96 minor( self, &mut m, 0, 2 );
97 let _det13 = m.determinant();
98 minor( self, &mut m, 0, 3 );
99 let _det14 = m.determinant();
100
101 _a11 * _det11 - _a12 * _det12 + _a13 * _det13 - _a14 * _det14
102 }
103
104 pub fn inverse( &self ) -> Option< Self >
107 where
108 Mat< 3, 3, E, Descriptor > :
109 RawSliceMut< Scalar = E > +
110 ScalarMut< Scalar = E, Index = Ix2 > +
111 ConstLayout< Index = Ix2 > +
112 IndexingMut< Scalar = E, Index = Ix2 >
113 {
114 let det = self.determinant();
115
116 if det == E::zero() { return None; }
117
118 let mut cfm = Mat3::default();
119 let mut cf = | i, j |
120 {
121 cofactor( self, &mut cfm, i, j )
122 };
123
124 let adj = Self::from_column_major
125 ([
126 cf( 0, 0 ), cf( 0, 1 ), cf( 0, 2 ), cf( 0, 3 ),
127 cf( 1, 0 ), cf( 1, 1 ), cf( 1, 2 ), cf( 1, 3 ),
128 cf( 2, 0 ), cf( 2, 1 ), cf( 2, 2 ), cf( 2, 3 ),
129 cf( 3, 0 ), cf( 3, 1 ), cf( 3, 2 ), cf( 3, 3 ),
130 ]);
131
132 Some( adj / det )
133 }
134}