| Geant4 Cross Reference |
1 // Copyright (C) 2010, Guy Barrand. All rights 1 // Copyright (C) 2010, Guy Barrand. All rights reserved.
2 // See the file tools.license for terms. 2 // See the file tools.license for terms.
3 3
4 #ifndef tools_mat4 4 #ifndef tools_mat4
5 #define tools_mat4 5 #define tools_mat4
6 6
7 #include "mat" 7 #include "mat"
8 8
9 //#include <cmath> 9 //#include <cmath>
10 10
11 namespace tools { 11 namespace tools {
12 12
13 template <class T> 13 template <class T>
14 class mat4 : public mat<T,4> { 14 class mat4 : public mat<T,4> {
15 typedef mat<T,4> parent; 15 typedef mat<T,4> parent;
16 typedef mat<T,4> pr; 16 typedef mat<T,4> pr;
17 public: 17 public:
18 #ifdef TOOLS_MEM 18 #ifdef TOOLS_MEM
19 mat4(bool a_inc = true):parent(a_inc) {} 19 mat4(bool a_inc = true):parent(a_inc) {}
20 #else 20 #else
21 mat4():parent() {} 21 mat4():parent() {}
22 #endif 22 #endif
23 mat4(const mat<T,4>& a_from):parent(a_from){ 23 mat4(const mat<T,4>& a_from):parent(a_from){}
24 virtual ~mat4() {} 24 virtual ~mat4() {}
25 public: 25 public:
26 mat4(const mat4& a_from):parent(a_from){} 26 mat4(const mat4& a_from):parent(a_from){}
27 mat4& operator=(const mat4& a_from){ 27 mat4& operator=(const mat4& a_from){
28 parent::operator=(a_from); 28 parent::operator=(a_from);
29 return *this; 29 return *this;
30 } 30 }
31 public: 31 public:
32 mat4(const T& a_00,const T& a_01,const T& a_ 32 mat4(const T& a_00,const T& a_01,const T& a_02,const T& a_03, //first row
33 const T& a_10,const T& a_11,const T& a_ 33 const T& a_10,const T& a_11,const T& a_12,const T& a_13, //second row
34 const T& a_20,const T& a_21,const T& a_ 34 const T& a_20,const T& a_21,const T& a_22,const T& a_23, //third row
35 const T& a_30,const T& a_31,const T& a_ 35 const T& a_30,const T& a_31,const T& a_32,const T& a_33) //fourth row
36 { 36 {
37 set_matrix(a_00,a_01,a_02,a_03, 37 set_matrix(a_00,a_01,a_02,a_03,
38 a_10,a_11,a_12,a_13, 38 a_10,a_11,a_12,a_13,
39 a_20,a_21,a_22,a_23, 39 a_20,a_21,a_22,a_23,
40 a_30,a_31,a_32,a_33); 40 a_30,a_31,a_32,a_33);
41 } 41 }
42 public: 42 public:
43 void set_matrix(const mat4<T>& a_m) {parent: 43 void set_matrix(const mat4<T>& a_m) {parent::set_matrix(a_m);}
44 void set_matrix( 44 void set_matrix(
45 const T& a_00,const T& a_01,const T& a_02, 45 const T& a_00,const T& a_01,const T& a_02,const T& a_03, //1 row
46 const T& a_10,const T& a_11,const T& a_12, 46 const T& a_10,const T& a_11,const T& a_12,const T& a_13, //2 row
47 const T& a_20,const T& a_21,const T& a_22, 47 const T& a_20,const T& a_21,const T& a_22,const T& a_23, //3 row
48 const T& a_30,const T& a_31,const T& a_32, 48 const T& a_30,const T& a_31,const T& a_32,const T& a_33) //4 row
49 { 49 {
50 //a_<R><C> 50 //a_<R><C>
51 //pr::m_vec[R + C * 4]; 51 //pr::m_vec[R + C * 4];
52 pr::m_vec[0] = a_00;pr::m_vec[4] = a_01;pr 52 pr::m_vec[0] = a_00;pr::m_vec[4] = a_01;pr::m_vec[ 8] = a_02;pr::m_vec[12] = a_03;
53 pr::m_vec[1] = a_10;pr::m_vec[5] = a_11;pr 53 pr::m_vec[1] = a_10;pr::m_vec[5] = a_11;pr::m_vec[ 9] = a_12;pr::m_vec[13] = a_13;
54 pr::m_vec[2] = a_20;pr::m_vec[6] = a_21;pr 54 pr::m_vec[2] = a_20;pr::m_vec[6] = a_21;pr::m_vec[10] = a_22;pr::m_vec[14] = a_23;
55 pr::m_vec[3] = a_30;pr::m_vec[7] = a_31;pr 55 pr::m_vec[3] = a_30;pr::m_vec[7] = a_31;pr::m_vec[11] = a_32;pr::m_vec[15] = a_33;
56 } 56 }
57 57
58 void set_scale(const T& a_s) {_set_scale(a_s 58 void set_scale(const T& a_s) {_set_scale(a_s,a_s,a_s,pr::m_vec);}
59 void set_scale(const T& a_1,const T& a_2,con 59 void set_scale(const T& a_1,const T& a_2,const T& a_3) {_set_scale(a_1,a_2,a_3,pr::m_vec);}
60 void set_translate(const T& a_x,const T& a_y 60 void set_translate(const T& a_x,const T& a_y,const T& a_z) {_set_translate(a_x,a_y,a_z,pr::m_vec);}
61 void set_rotate(const T& a_x,const T& a_y,co 61 void set_rotate(const T& a_x,const T& a_y,const T& a_z,const T& a_angle,T(*a_sin)(T),T(*a_cos)(T)) {
62 _set_rotate(a_x,a_y,a_z,a_angle,pr::m_vec, 62 _set_rotate(a_x,a_y,a_z,a_angle,pr::m_vec,a_sin,a_cos);
63 } 63 }
64 void set_ortho(const T& a_l,const T& a_r, 64 void set_ortho(const T& a_l,const T& a_r, //left,right
65 const T& a_b,const T& a_t, 65 const T& a_b,const T& a_t, //bottom,top
66 const T& a_n,const T& a_f) { 66 const T& a_n,const T& a_f) { //znear,zfar
67 // from man glOrtho. 67 // from man glOrtho.
68 T tx = -(a_r+a_l)/(a_r-a_l); 68 T tx = -(a_r+a_l)/(a_r-a_l);
69 T ty = -(a_t+a_b)/(a_t-a_b); 69 T ty = -(a_t+a_b)/(a_t-a_b);
70 T tz = -(a_f+a_n)/(a_f-a_n); 70 T tz = -(a_f+a_n)/(a_f-a_n);
71 71
72 T a_00,a_01,a_02,a_03; 72 T a_00,a_01,a_02,a_03;
73 T a_10,a_11,a_12,a_13; 73 T a_10,a_11,a_12,a_13;
74 T a_20,a_21,a_22,a_23; 74 T a_20,a_21,a_22,a_23;
75 T a_30,a_31,a_32,a_33; 75 T a_30,a_31,a_32,a_33;
76 a_00 = 2/(a_r-a_l);a_01 = 0;a_02 76 a_00 = 2/(a_r-a_l);a_01 = 0;a_02 = 0;a_03 = tx;
77 a_10 = 0;a_11 = 2/(a_t-a_b);a_12 77 a_10 = 0;a_11 = 2/(a_t-a_b);a_12 = 0;a_13 = ty;
78 a_20 = 0;a_21 = 0;a_22 78 a_20 = 0;a_21 = 0;a_22 = -2/(a_f-a_n);a_23 = tz;
79 a_30 = 0;a_31 = 0;a_32 79 a_30 = 0;a_31 = 0;a_32 = 0;a_33 = 1;
80 set_matrix(a_00,a_01,a_02,a_03, //1 row 80 set_matrix(a_00,a_01,a_02,a_03, //1 row
81 a_10,a_11,a_12,a_13, //2 row 81 a_10,a_11,a_12,a_13, //2 row
82 a_20,a_21,a_22,a_23, //3 row 82 a_20,a_21,a_22,a_23, //3 row
83 a_30,a_31,a_32,a_33); //4 row 83 a_30,a_31,a_32,a_33); //4 row
84 84
85 //NOTE : Z(x,y, z,1) = -2z/(f-n)+tz = [-2z 85 //NOTE : Z(x,y, z,1) = -2z/(f-n)+tz = [-2z-f-n]/(f-n)
86 // W(x,y, z,1) = 1 -> Z/W = Z 86 // W(x,y, z,1) = 1 -> Z/W = Z
87 // Z(x,y,-n,1) = -1 87 // Z(x,y,-n,1) = -1
88 // Z(x,y,-f,1) = 1 88 // Z(x,y,-f,1) = 1
89 89
90 // Z(x,y,-(f+n)/2,1) = 0 90 // Z(x,y,-(f+n)/2,1) = 0
91 91
92 // X(x,0, z,1) = 2x/(r-l)+tx = [2x-r 92 // X(x,0, z,1) = 2x/(r-l)+tx = [2x-r-l]/(r-l)
93 // X(r,0, z,1) = 1 93 // X(r,0, z,1) = 1
94 94
95 // the view direction is then (0,0,1) in t 95 // the view direction is then (0,0,1) in the final projection.
96 } 96 }
97 void set_frustum(const T& a_l,const T& a_r, 97 void set_frustum(const T& a_l,const T& a_r, //left,right
98 const T& a_b,const T& a_t, 98 const T& a_b,const T& a_t, //bottom,top
99 const T& a_n,const T& a_f) 99 const T& a_n,const T& a_f) { //znear,zfar
100 // from man glFrustum. 100 // from man glFrustum.
101 T A = (a_r+a_l)/(a_r-a_l); 101 T A = (a_r+a_l)/(a_r-a_l);
102 T B = (a_t+a_b)/(a_t-a_b); 102 T B = (a_t+a_b)/(a_t-a_b);
103 T C = -(a_f+a_n)/(a_f-a_n); 103 T C = -(a_f+a_n)/(a_f-a_n);
104 T D = -(2*a_f*a_n)/(a_f-a_n); 104 T D = -(2*a_f*a_n)/(a_f-a_n);
105 105
106 T a_00,a_01,a_02,a_03; 106 T a_00,a_01,a_02,a_03;
107 T a_10,a_11,a_12,a_13; 107 T a_10,a_11,a_12,a_13;
108 T a_20,a_21,a_22,a_23; 108 T a_20,a_21,a_22,a_23;
109 T a_30,a_31,a_32,a_33; 109 T a_30,a_31,a_32,a_33;
110 a_00 = (2*a_n)/(a_r-a_l);a_01 = 110 a_00 = (2*a_n)/(a_r-a_l);a_01 = 0;a_02 = A;a_03 = 0;
111 a_10 = 0;a_11 = (2*a_n)/(a 111 a_10 = 0;a_11 = (2*a_n)/(a_t-a_b);a_12 = B;a_13 = 0;
112 a_20 = 0;a_21 = 112 a_20 = 0;a_21 = 0;a_22 = C;a_23 = D;
113 a_30 = 0;a_31 = 113 a_30 = 0;a_31 = 0;a_32 = -1;a_33 = 0;
114 set_matrix(a_00,a_01,a_02,a_03, //1 row 114 set_matrix(a_00,a_01,a_02,a_03, //1 row
115 a_10,a_11,a_12,a_13, //2 row 115 a_10,a_11,a_12,a_13, //2 row
116 a_20,a_21,a_22,a_23, //3 row 116 a_20,a_21,a_22,a_23, //3 row
117 a_30,a_31,a_32,a_33); //4 row 117 a_30,a_31,a_32,a_33); //4 row
118 118
119 //NOTE : Z(x,y, z,1) = C*z+D = -[(f+n)*z+2 119 //NOTE : Z(x,y, z,1) = C*z+D = -[(f+n)*z+2*f*n]/(f-n)
120 120
121 // Z(x,y,-n,1) = -[-fn-nn+2fn]/(f-n) 121 // Z(x,y,-n,1) = -[-fn-nn+2fn]/(f-n) = -[fn-nn]/(f-n) = -n
122 // W(x,y,-n,1) = n 122 // W(x,y,-n,1) = n
123 // -> Z/W(x,y,-n,1) = -1 123 // -> Z/W(x,y,-n,1) = -1
124 124
125 // Z(x,y,-2fn/(f+n),1) = 0 125 // Z(x,y,-2fn/(f+n),1) = 0
126 // -> Z/W = 0 126 // -> Z/W = 0
127 127
128 // Z(x,y,-f,1) = -[-ff-fn+2fn]/(f-n) 128 // Z(x,y,-f,1) = -[-ff-fn+2fn]/(f-n) = -[fn-ff]/(f-n) = f
129 // W(x,y,-f,1) = f 129 // W(x,y,-f,1) = f
130 // -> Z/W(x,y,-f,1) = 1 130 // -> Z/W(x,y,-f,1) = 1
131 131
132 // X(x,0, z,1) = 2nx/(r-l)+z(r+l)/(r 132 // X(x,0, z,1) = 2nx/(r-l)+z(r+l)/(r-l) = [2nx+zr+zl]/(r-l)
133 // X(r,0,-n,1) = [nr-nl]/(r-l) = n 133 // X(r,0,-n,1) = [nr-nl]/(r-l) = n
134 // W(r,0,-n,1) = n 134 // W(r,0,-n,1) = n
135 // -> X/W(r,0,-n,1) = 1 135 // -> X/W(r,0,-n,1) = 1
136 136
137 // X(l,0,-n,1) = (2nl-n(r+l))/(r-l) 137 // X(l,0,-n,1) = (2nl-n(r+l))/(r-l) = -n
138 // W(l,0,-n,1) = n 138 // W(l,0,-n,1) = n
139 // -> X/W(l,0,-n,1) = -1 139 // -> X/W(l,0,-n,1) = -1
140 140
141 // lrbt corners are in plane z=-1 at xy=+/ 141 // lrbt corners are in plane z=-1 at xy=+/-1.
142 142
143 // eye ? 143 // eye ?
144 // eye before ? (0,0,z,1) -> (zA,zB,zC+D,- 144 // eye before ? (0,0,z,1) -> (zA,zB,zC+D,-z) /W -> (-A,-B,-(C+D/z),1)
145 // z=0 -> (0,0,D=-2fn(f-n),0) 145 // z=0 -> (0,0,D=-2fn(f-n),0)
146 146
147 } 147 }
148 148
149 void get_translate(T& a_x,T& a_y,T& a_z) con 149 void get_translate(T& a_x,T& a_y,T& a_z) const {
150 a_x = pr::m_vec[12]; 150 a_x = pr::m_vec[12];
151 a_y = pr::m_vec[13]; 151 a_y = pr::m_vec[13];
152 a_z = pr::m_vec[14]; 152 a_z = pr::m_vec[14];
153 } 153 }
154 154
155 void no_translate() { 155 void no_translate() {
156 pr::m_vec[12] = 0; 156 pr::m_vec[12] = 0;
157 pr::m_vec[13] = 0; 157 pr::m_vec[13] = 0;
158 pr::m_vec[14] = 0; 158 pr::m_vec[14] = 0;
159 } 159 }
160 160
161 //void set_translate_only(const T& a_x,const T 161 //void set_translate_only(const T& a_x,const T& a_y,const T& a_z) {
162 // pr::m_vec[12] = a_x; 162 // pr::m_vec[12] = a_x;
163 // pr::m_vec[13] = a_y; 163 // pr::m_vec[13] = a_y;
164 // pr::m_vec[14] = a_z; 164 // pr::m_vec[14] = a_z;
165 //} 165 //}
166 166
167 public: 167 public:
168 void mul_4(T& a_x,T& a_y,T& a_z,T& a_p) cons 168 void mul_4(T& a_x,T& a_y,T& a_z,T& a_p) const {
169 // a_[x,y,z,p] = this * a_[x,y,z,p] 169 // a_[x,y,z,p] = this * a_[x,y,z,p]
170 //pr::m_vec[R + C * 4]; 170 //pr::m_vec[R + C * 4];
171 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m 171 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m_vec[ 8] = 02;pr::m_vec[12] = 03;
172 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m 172 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m_vec[ 9] = 12;pr::m_vec[13] = 13;
173 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m 173 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m_vec[10] = 22;pr::m_vec[14] = 23;
174 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m 174 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m_vec[11] = 32;pr::m_vec[15] = 33;
175 T x = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr 175 T x = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr::m_vec[ 8]*a_z+pr::m_vec[12]*a_p;
176 T y = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr 176 T y = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr::m_vec[ 9]*a_z+pr::m_vec[13]*a_p;
177 T z = pr::m_vec[2]*a_x+pr::m_vec[6]*a_y+pr 177 T z = pr::m_vec[2]*a_x+pr::m_vec[6]*a_y+pr::m_vec[10]*a_z+pr::m_vec[14]*a_p;
178 T p = pr::m_vec[3]*a_x+pr::m_vec[7]*a_y+pr 178 T p = pr::m_vec[3]*a_x+pr::m_vec[7]*a_y+pr::m_vec[11]*a_z+pr::m_vec[15]*a_p;
179 a_x = x; 179 a_x = x;
180 a_y = y; 180 a_y = y;
181 a_z = z; 181 a_z = z;
182 a_p = p; 182 a_p = p;
183 } 183 }
184 void mul_4_opt(T& a_x,T& a_y,T& a_z,T& a_p,T 184 void mul_4_opt(T& a_x,T& a_y,T& a_z,T& a_p,T a_tmp[4]) const {
185 // a_[x,y,z,p] = this * a_[x,y,z,p] 185 // a_[x,y,z,p] = this * a_[x,y,z,p]
186 //pr::m_vec[R + C * 4]; 186 //pr::m_vec[R + C * 4];
187 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m 187 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m_vec[ 8] = 02;pr::m_vec[12] = 03;
188 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m 188 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m_vec[ 9] = 12;pr::m_vec[13] = 13;
189 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m 189 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m_vec[10] = 22;pr::m_vec[14] = 23;
190 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m 190 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m_vec[11] = 32;pr::m_vec[15] = 33;
191 a_tmp[0] = pr::m_vec[0]*a_x+pr::m_vec[4]*a 191 a_tmp[0] = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr::m_vec[ 8]*a_z+pr::m_vec[12]*a_p;
192 a_tmp[1] = pr::m_vec[1]*a_x+pr::m_vec[5]*a 192 a_tmp[1] = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr::m_vec[ 9]*a_z+pr::m_vec[13]*a_p;
193 a_tmp[2] = pr::m_vec[2]*a_x+pr::m_vec[6]*a 193 a_tmp[2] = pr::m_vec[2]*a_x+pr::m_vec[6]*a_y+pr::m_vec[10]*a_z+pr::m_vec[14]*a_p;
194 a_tmp[3] = pr::m_vec[3]*a_x+pr::m_vec[7]*a 194 a_tmp[3] = pr::m_vec[3]*a_x+pr::m_vec[7]*a_y+pr::m_vec[11]*a_z+pr::m_vec[15]*a_p;
195 a_x = a_tmp[0]; 195 a_x = a_tmp[0];
196 a_y = a_tmp[1]; 196 a_y = a_tmp[1];
197 a_z = a_tmp[2]; 197 a_z = a_tmp[2];
198 a_p = a_tmp[3]; 198 a_p = a_tmp[3];
199 } 199 }
200 void mul_3(T& a_x,T& a_y,T& a_z) const { 200 void mul_3(T& a_x,T& a_y,T& a_z) const {
201 // a_[x,y,z] = this * a_[x,y,z] 201 // a_[x,y,z] = this * a_[x,y,z]
202 //pr::m_vec[R + C * 4]; 202 //pr::m_vec[R + C * 4];
203 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m 203 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m_vec[ 8] = 02;pr::m_vec[12] = 03;
204 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m 204 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m_vec[ 9] = 12;pr::m_vec[13] = 13;
205 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m 205 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m_vec[10] = 22;pr::m_vec[14] = 23;
206 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m 206 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m_vec[11] = 32;pr::m_vec[15] = 33;
207 T x = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr 207 T x = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr::m_vec[ 8]*a_z+pr::m_vec[12];
208 T y = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr 208 T y = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr::m_vec[ 9]*a_z+pr::m_vec[13];
209 T z = pr::m_vec[2]*a_x+pr::m_vec[6]*a_y+pr 209 T z = pr::m_vec[2]*a_x+pr::m_vec[6]*a_y+pr::m_vec[10]*a_z+pr::m_vec[14];
210 a_x = x; 210 a_x = x;
211 a_y = y; 211 a_y = y;
212 a_z = z; 212 a_z = z;
213 } 213 }
214 void mul_3_opt(T& a_x,T& a_y,T& a_z,T a_tmp[ 214 void mul_3_opt(T& a_x,T& a_y,T& a_z,T a_tmp[3]) const {
215 // a_[x,y,z] = this * a_[x,y,z] 215 // a_[x,y,z] = this * a_[x,y,z]
216 //pr::m_vec[R + C * 4]; 216 //pr::m_vec[R + C * 4];
217 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m 217 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m_vec[ 8] = 02;pr::m_vec[12] = 03;
218 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m 218 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m_vec[ 9] = 12;pr::m_vec[13] = 13;
219 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m 219 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m_vec[10] = 22;pr::m_vec[14] = 23;
220 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m 220 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m_vec[11] = 32;pr::m_vec[15] = 33;
221 a_tmp[0] = pr::m_vec[0]*a_x+pr::m_vec[4]*a 221 a_tmp[0] = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr::m_vec[ 8]*a_z+pr::m_vec[12];
222 a_tmp[1] = pr::m_vec[1]*a_x+pr::m_vec[5]*a 222 a_tmp[1] = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr::m_vec[ 9]*a_z+pr::m_vec[13];
223 a_tmp[2] = pr::m_vec[2]*a_x+pr::m_vec[6]*a 223 a_tmp[2] = pr::m_vec[2]*a_x+pr::m_vec[6]*a_y+pr::m_vec[10]*a_z+pr::m_vec[14];
224 a_x = a_tmp[0]; 224 a_x = a_tmp[0];
225 a_y = a_tmp[1]; 225 a_y = a_tmp[1];
226 a_z = a_tmp[2]; 226 a_z = a_tmp[2];
227 } 227 }
228 void mul_2(T& a_x,T& a_y) const { 228 void mul_2(T& a_x,T& a_y) const {
229 // a_[x,y] = this * a_[x,y] 229 // a_[x,y] = this * a_[x,y]
230 //pr::m_vec[R + C * 4]; 230 //pr::m_vec[R + C * 4];
231 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m 231 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m_vec[ 8] = 02;pr::m_vec[12] = 03;
232 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m 232 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m_vec[ 9] = 12;pr::m_vec[13] = 13;
233 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m 233 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m_vec[10] = 22;pr::m_vec[14] = 23;
234 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m 234 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m_vec[11] = 32;pr::m_vec[15] = 33;
235 T x = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr 235 T x = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr::m_vec[12];
236 T y = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr 236 T y = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr::m_vec[13];
237 a_x = x; 237 a_x = x;
238 a_y = y; 238 a_y = y;
239 } 239 }
240 240
241 void mul_dir_3(T& a_x,T& a_y,T& a_z) const { 241 void mul_dir_3(T& a_x,T& a_y,T& a_z) const {
242 // used to multiply normals. 242 // used to multiply normals.
243 // a_[x,y,z] = rot_part(this) * a_[x,y,z] 243 // a_[x,y,z] = rot_part(this) * a_[x,y,z]
244 244
245 T x = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr 245 T x = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr::m_vec[ 8]*a_z;
246 T y = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr 246 T y = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr::m_vec[ 9]*a_z;
247 T z = pr::m_vec[2]*a_x+pr::m_vec[6]*a_y+pr 247 T z = pr::m_vec[2]*a_x+pr::m_vec[6]*a_y+pr::m_vec[10]*a_z;
248 a_x = x; 248 a_x = x;
249 a_y = y; 249 a_y = y;
250 a_z = z; 250 a_z = z;
251 } 251 }
252 void mul_dir_3_opt(T& a_x,T& a_y,T& a_z,T a_ 252 void mul_dir_3_opt(T& a_x,T& a_y,T& a_z,T a_tmp[3]) const {
253 // used to multiply normals. 253 // used to multiply normals.
254 // a_[x,y,z] = rot_part(this) * a_[x,y,z] 254 // a_[x,y,z] = rot_part(this) * a_[x,y,z]
255 255
256 a_tmp[0] = pr::m_vec[0]*a_x+pr::m_vec[4]*a 256 a_tmp[0] = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr::m_vec[ 8]*a_z;
257 a_tmp[1] = pr::m_vec[1]*a_x+pr::m_vec[5]*a 257 a_tmp[1] = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr::m_vec[ 9]*a_z;
258 a_tmp[2] = pr::m_vec[2]*a_x+pr::m_vec[6]*a 258 a_tmp[2] = pr::m_vec[2]*a_x+pr::m_vec[6]*a_y+pr::m_vec[10]*a_z;
259 a_x = a_tmp[0]; 259 a_x = a_tmp[0];
260 a_y = a_tmp[1]; 260 a_y = a_tmp[1];
261 a_z = a_tmp[2]; 261 a_z = a_tmp[2];
262 } 262 }
263 263
264 void mul_trans_3(T& a_x,T& a_y,T& a_z) const 264 void mul_trans_3(T& a_x,T& a_y,T& a_z) const {
265 // used in sg::healpix. 265 // used in sg::healpix.
266 // a_[x,y,z] = trans_part(this) * a_[x,y,z 266 // a_[x,y,z] = trans_part(this) * a_[x,y,z]
267 a_x += pr::m_vec[12]; 267 a_x += pr::m_vec[12];
268 a_y += pr::m_vec[13]; 268 a_y += pr::m_vec[13];
269 a_z += pr::m_vec[14]; 269 a_z += pr::m_vec[14];
270 } 270 }
271 271
272 template <class VEC> 272 template <class VEC>
273 void mul_dir_3(VEC& a_dir) const {mul_dir_3( 273 void mul_dir_3(VEC& a_dir) const {mul_dir_3(a_dir[0],a_dir[1],a_dir[2]);}
274 274
275 void mul_scale(const T& a_sx,const T& a_sy,c 275 void mul_scale(const T& a_sx,const T& a_sy,const T& a_sz) {
276 // this = this * mat4_scale(a_s[x,y,z] 276 // this = this * mat4_scale(a_s[x,y,z]
277 //pr::m_vec[R + C * 4]; 277 //pr::m_vec[R + C * 4];
278 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m 278 //pr::m_vec[0]= 00;pr::m_vec[4] = 01;pr::m_vec[ 8] = 02;pr::m_vec[12] = 03;
279 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m 279 //pr::m_vec[1]= 10;pr::m_vec[5] = 11;pr::m_vec[ 9] = 12;pr::m_vec[13] = 13;
280 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m 280 //pr::m_vec[2]= 20;pr::m_vec[6] = 21;pr::m_vec[10] = 22;pr::m_vec[14] = 23;
281 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m 281 //pr::m_vec[3]= 30;pr::m_vec[7] = 31;pr::m_vec[11] = 32;pr::m_vec[15] = 33;
282 pr::m_vec[0] *= a_sx; 282 pr::m_vec[0] *= a_sx;
283 pr::m_vec[1] *= a_sx; 283 pr::m_vec[1] *= a_sx;
284 pr::m_vec[2] *= a_sx; 284 pr::m_vec[2] *= a_sx;
285 pr::m_vec[3] *= a_sx; 285 pr::m_vec[3] *= a_sx;
286 286
287 pr::m_vec[4] *= a_sy; 287 pr::m_vec[4] *= a_sy;
288 pr::m_vec[5] *= a_sy; 288 pr::m_vec[5] *= a_sy;
289 pr::m_vec[6] *= a_sy; 289 pr::m_vec[6] *= a_sy;
290 pr::m_vec[7] *= a_sy; 290 pr::m_vec[7] *= a_sy;
291 291
292 pr::m_vec[ 8] *= a_sz; 292 pr::m_vec[ 8] *= a_sz;
293 pr::m_vec[ 9] *= a_sz; 293 pr::m_vec[ 9] *= a_sz;
294 pr::m_vec[10] *= a_sz; 294 pr::m_vec[10] *= a_sz;
295 pr::m_vec[11] *= a_sz; 295 pr::m_vec[11] *= a_sz;
296 } 296 }
297 void mul_scale(const T& a_s) { 297 void mul_scale(const T& a_s) {
298 pr::m_vec[0] *= a_s; 298 pr::m_vec[0] *= a_s;
299 pr::m_vec[1] *= a_s; 299 pr::m_vec[1] *= a_s;
300 pr::m_vec[2] *= a_s; 300 pr::m_vec[2] *= a_s;
301 pr::m_vec[3] *= a_s; 301 pr::m_vec[3] *= a_s;
302 302
303 pr::m_vec[4] *= a_s; 303 pr::m_vec[4] *= a_s;
304 pr::m_vec[5] *= a_s; 304 pr::m_vec[5] *= a_s;
305 pr::m_vec[6] *= a_s; 305 pr::m_vec[6] *= a_s;
306 pr::m_vec[7] *= a_s; 306 pr::m_vec[7] *= a_s;
307 307
308 pr::m_vec[ 8] *= a_s; 308 pr::m_vec[ 8] *= a_s;
309 pr::m_vec[ 9] *= a_s; 309 pr::m_vec[ 9] *= a_s;
310 pr::m_vec[10] *= a_s; 310 pr::m_vec[10] *= a_s;
311 pr::m_vec[11] *= a_s; 311 pr::m_vec[11] *= a_s;
312 } 312 }
313 void mul_translate(const T& a_x,const T& a_y 313 void mul_translate(const T& a_x,const T& a_y,const T& a_z) {
314 pr::m_vec[12] = pr::m_vec[0]*a_x+pr::m_vec 314 pr::m_vec[12] = pr::m_vec[0]*a_x+pr::m_vec[4]*a_y+pr::m_vec[ 8]*a_z+pr::m_vec[12];
315 pr::m_vec[13] = pr::m_vec[1]*a_x+pr::m_vec 315 pr::m_vec[13] = pr::m_vec[1]*a_x+pr::m_vec[5]*a_y+pr::m_vec[ 9]*a_z+pr::m_vec[13];
316 pr::m_vec[14] = pr::m_vec[2]*a_x+pr::m_vec 316 pr::m_vec[14] = pr::m_vec[2]*a_x+pr::m_vec[6]*a_y+pr::m_vec[10]*a_z+pr::m_vec[14];
317 pr::m_vec[15] = pr::m_vec[3]*a_x+pr::m_vec 317 pr::m_vec[15] = pr::m_vec[3]*a_x+pr::m_vec[7]*a_y+pr::m_vec[11]*a_z+pr::m_vec[15];
318 } 318 }
319 319
320 void mul_rotate(const T& a_x,const T& a_y,co 320 void mul_rotate(const T& a_x,const T& a_y,const T& a_z,const T& a_angle,T(*a_sin)(T),T(*a_cos)(T)) {
321 T rot[16]; 321 T rot[16];
322 _set_rotate(a_x,a_y,a_z,a_angle,rot,a_sin, 322 _set_rotate(a_x,a_y,a_z,a_angle,rot,a_sin,a_cos);
323 parent::_mul_mtx(rot); 323 parent::_mul_mtx(rot);
324 } 324 }
325 325
326 void left_mul_rotate(const T& a_x,const T& a 326 void left_mul_rotate(const T& a_x,const T& a_y,const T& a_z,const T& a_angle,T(*a_sin)(T),T(*a_cos)(T)) {
327 T _m[16]; 327 T _m[16];
328 _set_rotate(a_x,a_y,a_z,a_angle,_m,a_sin,a 328 _set_rotate(a_x,a_y,a_z,a_angle,_m,a_sin,a_cos);
329 parent::_left_mul_mtx(_m); 329 parent::_left_mul_mtx(_m);
330 } 330 }
331 331
332 void left_mul_scale(const T& a_x,const T& a_ 332 void left_mul_scale(const T& a_x,const T& a_y,const T& a_z) {
333 T _m[16]; 333 T _m[16];
334 _set_scale(a_x,a_y,a_z,_m); 334 _set_scale(a_x,a_y,a_z,_m);
335 parent::_left_mul_mtx(_m); 335 parent::_left_mul_mtx(_m);
336 } 336 }
337 337
338 void left_mul_translate(const T& a_x,const T 338 void left_mul_translate(const T& a_x,const T& a_y,const T& a_z) {
339 T _m[16]; 339 T _m[16];
340 _set_translate(a_x,a_y,a_z,_m); 340 _set_translate(a_x,a_y,a_z,_m);
341 parent::_left_mul_mtx(_m); 341 parent::_left_mul_mtx(_m);
342 } 342 }
343 343
344 void v00(const T& a_value){pr::m_vec[0+0*4] 344 void v00(const T& a_value){pr::m_vec[0+0*4] = a_value;}
345 void v10(const T& a_value){pr::m_vec[1+0*4] 345 void v10(const T& a_value){pr::m_vec[1+0*4] = a_value;}
346 void v20(const T& a_value){pr::m_vec[2+0*4] 346 void v20(const T& a_value){pr::m_vec[2+0*4] = a_value;}
347 void v30(const T& a_value){pr::m_vec[3+0*4] 347 void v30(const T& a_value){pr::m_vec[3+0*4] = a_value;}
348 348
349 void v01(const T& a_value){pr::m_vec[0+1*4] 349 void v01(const T& a_value){pr::m_vec[0+1*4] = a_value;}
350 void v11(const T& a_value){pr::m_vec[1+1*4] 350 void v11(const T& a_value){pr::m_vec[1+1*4] = a_value;}
351 void v21(const T& a_value){pr::m_vec[2+1*4] 351 void v21(const T& a_value){pr::m_vec[2+1*4] = a_value;}
352 void v31(const T& a_value){pr::m_vec[3+1*4] 352 void v31(const T& a_value){pr::m_vec[3+1*4] = a_value;}
353 353
354 void v02(const T& a_value){pr::m_vec[0+2*4] 354 void v02(const T& a_value){pr::m_vec[0+2*4] = a_value;}
355 void v12(const T& a_value){pr::m_vec[1+2*4] 355 void v12(const T& a_value){pr::m_vec[1+2*4] = a_value;}
356 void v22(const T& a_value){pr::m_vec[2+2*4] 356 void v22(const T& a_value){pr::m_vec[2+2*4] = a_value;}
357 void v32(const T& a_value){pr::m_vec[3+2*4] 357 void v32(const T& a_value){pr::m_vec[3+2*4] = a_value;}
358 358
359 void v03(const T& a_value){pr::m_vec[0+3*4] 359 void v03(const T& a_value){pr::m_vec[0+3*4] = a_value;}
360 void v13(const T& a_value){pr::m_vec[1+3*4] 360 void v13(const T& a_value){pr::m_vec[1+3*4] = a_value;}
361 void v23(const T& a_value){pr::m_vec[2+3*4] 361 void v23(const T& a_value){pr::m_vec[2+3*4] = a_value;}
362 void v33(const T& a_value){pr::m_vec[3+3*4] 362 void v33(const T& a_value){pr::m_vec[3+3*4] = a_value;}
363 363
364 const T& v00() const {return pr::m_vec[0+0*4 364 const T& v00() const {return pr::m_vec[0+0*4];}
365 const T& v10() const {return pr::m_vec[1+0*4 365 const T& v10() const {return pr::m_vec[1+0*4];}
366 const T& v20() const {return pr::m_vec[2+0*4 366 const T& v20() const {return pr::m_vec[2+0*4];}
367 const T& v30() const {return pr::m_vec[3+0*4 367 const T& v30() const {return pr::m_vec[3+0*4];}
368 368
369 const T& v01() const {return pr::m_vec[0+1*4 369 const T& v01() const {return pr::m_vec[0+1*4];}
370 const T& v11() const {return pr::m_vec[1+1*4 370 const T& v11() const {return pr::m_vec[1+1*4];}
371 const T& v21() const {return pr::m_vec[2+1*4 371 const T& v21() const {return pr::m_vec[2+1*4];}
372 const T& v31() const {return pr::m_vec[3+1*4 372 const T& v31() const {return pr::m_vec[3+1*4];}
373 373
374 const T& v02() const {return pr::m_vec[0+2*4 374 const T& v02() const {return pr::m_vec[0+2*4];}
375 const T& v12() const {return pr::m_vec[1+2*4 375 const T& v12() const {return pr::m_vec[1+2*4];}
376 const T& v22() const {return pr::m_vec[2+2*4 376 const T& v22() const {return pr::m_vec[2+2*4];}
377 const T& v32() const {return pr::m_vec[3+2*4 377 const T& v32() const {return pr::m_vec[3+2*4];}
378 378
379 const T& v03() const {return pr::m_vec[0+3*4 379 const T& v03() const {return pr::m_vec[0+3*4];}
380 const T& v13() const {return pr::m_vec[1+3*4 380 const T& v13() const {return pr::m_vec[1+3*4];}
381 const T& v23() const {return pr::m_vec[2+3*4 381 const T& v23() const {return pr::m_vec[2+3*4];}
382 const T& v33() const {return pr::m_vec[3+3*4 382 const T& v33() const {return pr::m_vec[3+3*4];}
383 383
384 protected: 384 protected:
385 static void _set_translate(const T& a_x,cons 385 static void _set_translate(const T& a_x,const T& a_y,const T& a_z,T v[]) {
386 v[0] = T(1);v[4] = 0;v[ 8] = 0;v[12] 386 v[0] = T(1);v[4] = 0;v[ 8] = 0;v[12] = a_x;
387 v[1] = 0;v[5] = T(1);v[ 9] = 0;v[13] 387 v[1] = 0;v[5] = T(1);v[ 9] = 0;v[13] = a_y;
388 v[2] = 0;v[6] = 0;v[10] = T(1);v[14] 388 v[2] = 0;v[6] = 0;v[10] = T(1);v[14] = a_z;
389 v[3] = 0;v[7] = 0;v[11] = 0;v[15] 389 v[3] = 0;v[7] = 0;v[11] = 0;v[15] = T(1);
390 } 390 }
391 391
392 static void _set_scale(const T& a_1,const T& 392 static void _set_scale(const T& a_1,const T& a_2,const T& a_3,T v[]) {
393 v[0] = a_1;v[4] = 0;v[ 8] = 0;v[12] = 393 v[0] = a_1;v[4] = 0;v[ 8] = 0;v[12] = 0;
394 v[1] = 0;v[5] = a_2;v[ 9] = 0;v[13] = 394 v[1] = 0;v[5] = a_2;v[ 9] = 0;v[13] = 0;
395 v[2] = 0;v[6] = 0;v[10] = a_3;v[14] = 395 v[2] = 0;v[6] = 0;v[10] = a_3;v[14] = 0;
396 v[3] = 0;v[7] = 0;v[11] = 0;v[15] = 396 v[3] = 0;v[7] = 0;v[11] = 0;v[15] = T(1);
397 } 397 }
398 398
399 static void _set_rotate(const T& a_x,const T 399 static void _set_rotate(const T& a_x,const T& a_y,const T& a_z,const T& a_angle,T v[],T(*a_sin)(T),T(*a_cos)(T)) {
400 //WARNING : (a_x,a_y,a_z) must be a normal 400 //WARNING : (a_x,a_y,a_z) must be a normalized vector.
401 // v[] = exp(-a_angle*n(a_x,a_y,a_z)*Rs). 401 // v[] = exp(-a_angle*n(a_x,a_y,a_z)*Rs). It models the rotation of an object in the same coordinate system (active rotation).
402 T si = a_sin(a_angle); 402 T si = a_sin(a_angle);
403 T co = a_cos(a_angle); 403 T co = a_cos(a_angle);
404 T x = a_x; 404 T x = a_x;
405 T y = a_y; 405 T y = a_y;
406 T z = a_z; 406 T z = a_z;
407 T x2 = x*x; 407 T x2 = x*x;
408 T y2 = y*y; 408 T y2 = y*y;
409 T z2 = z*z; 409 T z2 = z*z;
410 T xy = x*y; 410 T xy = x*y;
411 T xz = x*z; 411 T xz = x*z;
412 T yz = y*z; 412 T yz = y*z;
413 v[0] = x2+(1-x2)*co;v[4] = xy*(1-co)-z*si; 413 v[0] = x2+(1-x2)*co;v[4] = xy*(1-co)-z*si;v[ 8] = xz*(1-co)+y*si;v[12] = 0;
414 v[1] = xy*(1-co)+z*si;v[5] = y2+(1-y2)*co; 414 v[1] = xy*(1-co)+z*si;v[5] = y2+(1-y2)*co;v[ 9] = yz*(1-co)-x*si;v[13] = 0;
415 v[2] = xz*(1-co)-y*si;v[6] = yz*(1-co)+x*s 415 v[2] = xz*(1-co)-y*si;v[6] = yz*(1-co)+x*si;v[10] = z2+(1-z2)*co;v[14] = 0;
416 v[3] = 0;v[7] = 0;v[ 416 v[3] = 0;v[7] = 0;v[11] = 0;v[15] = 1;
417 417
418 // If : 418 // If :
419 // n =(x,y,z) 419 // n =(x,y,z)
420 // n2 = x2+y2+z2 = 1 420 // n2 = x2+y2+z2 = 1
421 // n.E = x*E1+y*E2+z*E3 421 // n.E = x*E1+y*E2+z*E3
422 // with : 422 // with :
423 // E1 E2 E3 423 // E1 E2 E3
424 // 0 0 0 0 0 -1 0 1 0 424 // 0 0 0 0 0 -1 0 1 0
425 // 0 0 1 0 0 0 -1 0 0 425 // 0 0 1 0 0 0 -1 0 0
426 // 0 -1 0 1 0 0 0 0 0 426 // 0 -1 0 1 0 0 0 0 0
427 // 427 //
428 // R(r,c) = cos(theta)*Id(r,c)+(1-cos(thet 428 // R(r,c) = cos(theta)*Id(r,c)+(1-cos(theta))*nr*nc-sin(theta)*(n.E)(r,c)
429 // 429 //
430 // R = exp(-theta*(n.E)) // active rotatio 430 // R = exp(-theta*(n.E)) // active rotation (it models the rotation of an object by theta along n).
431 431
432 } 432 }
433 433
434 public: 434 public:
435 void mul_mtx_rot_root(const T& a_00,const T& 435 void mul_mtx_rot_root(const T& a_00,const T& a_01,const T& a_02, //1 row
436 const T& a_10,const T& 436 const T& a_10,const T& a_11,const T& a_12, //2 row
437 const T& a_20,const T& 437 const T& a_20,const T& a_21,const T& a_22) //3 row
438 { 438 {
439 T* tv = pr::m_vec; 439 T* tv = pr::m_vec;
440 //pr::m_vec[0] = 00;pr::m_vec[4] = 01;pr:: 440 //pr::m_vec[0] = 00;pr::m_vec[4] = 01;pr::m_vec[ 8] = 02;pr::m_vec[12] = 03;
441 //pr::m_vec[1] = 10;pr::m_vec[5] = 11;pr:: 441 //pr::m_vec[1] = 10;pr::m_vec[5] = 11;pr::m_vec[ 9] = 12;pr::m_vec[13] = 13;
442 //pr::m_vec[2] = 20;pr::m_vec[6] = 21;pr:: 442 //pr::m_vec[2] = 20;pr::m_vec[6] = 21;pr::m_vec[10] = 22;pr::m_vec[14] = 23;
443 //pr::m_vec[3] = 30;pr::m_vec[7] = 31;pr:: 443 //pr::m_vec[3] = 30;pr::m_vec[7] = 31;pr::m_vec[11] = 32;pr::m_vec[15] = 33;
444 444
445 T tv_0 = tv[0]; 445 T tv_0 = tv[0];
446 T tv_1 = tv[1]; 446 T tv_1 = tv[1];
447 T tv_2 = tv[2]; 447 T tv_2 = tv[2];
448 T tv_3 = tv[3]; 448 T tv_3 = tv[3];
449 T tv_4 = tv[4]; 449 T tv_4 = tv[4];
450 T tv_5 = tv[5]; 450 T tv_5 = tv[5];
451 T tv_6 = tv[6]; 451 T tv_6 = tv[6];
452 T tv_7 = tv[7]; 452 T tv_7 = tv[7];
453 T tv_8 = tv[8]; 453 T tv_8 = tv[8];
454 T tv_9 = tv[9]; 454 T tv_9 = tv[9];
455 T tv_10 = tv[10]; 455 T tv_10 = tv[10];
456 T tv_11 = tv[11]; 456 T tv_11 = tv[11];
457 T tv_12 = tv[12]; 457 T tv_12 = tv[12];
458 T tv_13 = tv[13]; 458 T tv_13 = tv[13];
459 T tv_14 = tv[14]; 459 T tv_14 = tv[14];
460 T tv_15 = tv[15]; 460 T tv_15 = tv[15];
461 461
462 T fv_0 = a_00; 462 T fv_0 = a_00;
463 T fv_1 = a_10; 463 T fv_1 = a_10;
464 T fv_2 = a_20; 464 T fv_2 = a_20;
465 //T fv_3 = 0; 465 //T fv_3 = 0;
466 466
467 T fv_4 = a_01; 467 T fv_4 = a_01;
468 T fv_5 = a_11; 468 T fv_5 = a_11;
469 T fv_6 = a_21; 469 T fv_6 = a_21;
470 //T fv_7 = 0; 470 //T fv_7 = 0;
471 471
472 T fv_8 = a_02; 472 T fv_8 = a_02;
473 T fv_9 = a_12; 473 T fv_9 = a_12;
474 T fv_10 = a_22; 474 T fv_10 = a_22;
475 //T fv_11 = 0; 475 //T fv_11 = 0;
476 476
477 //T fv_12 = 0; 477 //T fv_12 = 0;
478 //T fv_13 = 0; 478 //T fv_13 = 0;
479 //T fv_14 = 0; 479 //T fv_14 = 0;
480 //T fv_15 = 1; 480 //T fv_15 = 1;
481 481
482 tv[0] = tv_0*fv_0+tv_4*fv_1+ tv_8*fv_2; 482 tv[0] = tv_0*fv_0+tv_4*fv_1+ tv_8*fv_2;
483 tv[1] = tv_1*fv_0+tv_5*fv_1+ tv_9*fv_2; 483 tv[1] = tv_1*fv_0+tv_5*fv_1+ tv_9*fv_2;
484 tv[2] = tv_2*fv_0+tv_6*fv_1+tv_10*fv_2; 484 tv[2] = tv_2*fv_0+tv_6*fv_1+tv_10*fv_2;
485 tv[3] = tv_3*fv_0+tv_7*fv_1+tv_11*fv_2; 485 tv[3] = tv_3*fv_0+tv_7*fv_1+tv_11*fv_2;
486 486
487 tv[4] = tv_0*fv_4+tv_4*fv_5+ tv_8*fv_6; 487 tv[4] = tv_0*fv_4+tv_4*fv_5+ tv_8*fv_6;
488 tv[5] = tv_1*fv_4+tv_5*fv_5+ tv_9*fv_6; 488 tv[5] = tv_1*fv_4+tv_5*fv_5+ tv_9*fv_6;
489 tv[6] = tv_2*fv_4+tv_6*fv_5+tv_10*fv_6; 489 tv[6] = tv_2*fv_4+tv_6*fv_5+tv_10*fv_6;
490 tv[7] = tv_3*fv_4+tv_7*fv_5+tv_11*fv_6; 490 tv[7] = tv_3*fv_4+tv_7*fv_5+tv_11*fv_6;
491 491
492 tv[8] = tv_0*fv_8+tv_4*fv_9+ tv_8*fv_10; 492 tv[8] = tv_0*fv_8+tv_4*fv_9+ tv_8*fv_10;
493 tv[9] = tv_1*fv_8+tv_5*fv_9+ tv_9*fv_10; 493 tv[9] = tv_1*fv_8+tv_5*fv_9+ tv_9*fv_10;
494 tv[10] = tv_2*fv_8+tv_6*fv_9+tv_10*fv_10; 494 tv[10] = tv_2*fv_8+tv_6*fv_9+tv_10*fv_10;
495 tv[11] = tv_3*fv_8+tv_7*fv_9+tv_11*fv_10; 495 tv[11] = tv_3*fv_8+tv_7*fv_9+tv_11*fv_10;
496 496
497 tv[12] = tv_12; 497 tv[12] = tv_12;
498 tv[13] = tv_13; 498 tv[13] = tv_13;
499 tv[14] = tv_14; 499 tv[14] = tv_14;
500 tv[15] = tv_15; 500 tv[15] = tv_15;
501 } 501 }
502 502
503 template <class MAT3> 503 template <class MAT3>
504 void set_mat3(MAT3& a_m3) { 504 void set_mat3(MAT3& a_m3) {
505 a_m3.set_matrix(pr::m_vec[0],pr::m_vec[4], 505 a_m3.set_matrix(pr::m_vec[0],pr::m_vec[4],pr::m_vec[ 8], //1 row
506 pr::m_vec[1],pr::m_vec[5], 506 pr::m_vec[1],pr::m_vec[5],pr::m_vec[ 9], //2 row
507 pr::m_vec[2],pr::m_vec[6], 507 pr::m_vec[2],pr::m_vec[6],pr::m_vec[10]); //3 row
508 } 508 }
509 private:static void check_instantiation() {mat 509 private:static void check_instantiation() {mat4<float> dummy;}
510 }; 510 };
511 511
512 ////////////////////////////////////////////// 512 ////////////////////////////////////////////////
513 /// common matrices : //////////////////////// 513 /// common matrices : //////////////////////////
514 ////////////////////////////////////////////// 514 ////////////////////////////////////////////////
515 515
516 #ifdef TOOLS_MEM 516 #ifdef TOOLS_MEM
517 template <class T> inline const mat4<T>& mat4_ 517 template <class T> inline const mat4<T>& mat4_zero() {static const mat4<T> s_v(false);return s_v;}
518 #else 518 #else
519 template <class T> inline const mat4<T>& mat4_ 519 template <class T> inline const mat4<T>& mat4_zero() {static const mat4<T> s_v;return s_v;}
520 #endif 520 #endif
521 521
522 /* 522 /*
523 template <class T> 523 template <class T>
524 inline const mat4<T>& mat4_identity() { 524 inline const mat4<T>& mat4_identity() {
525 static mat4<T> s_v(false); 525 static mat4<T> s_v(false);
526 static bool s_first = true; 526 static bool s_first = true;
527 if(s_first) {s_v.set_identity();s_first=fals 527 if(s_first) {s_v.set_identity();s_first=false;}
528 return s_v; 528 return s_v;
529 } 529 }
530 */ 530 */
531 531
532 //for sf, mf : 532 //for sf, mf :
533 template <class T> 533 template <class T>
534 inline const T* get_data(const mat4<T>& a_v) { 534 inline const T* get_data(const mat4<T>& a_v) {return a_v.data();}
535 535
536 } 536 }
537 537
538 #include <ostream> 538 #include <ostream>
539 539
540 namespace tools { 540 namespace tools {
541 541
542 template <class T> 542 template <class T>
543 inline std::ostream& operator<<(std::ostream& 543 inline std::ostream& operator<<(std::ostream& a_out,const mat4<T>& a_mtx){
544 const T* v = a_mtx.data(); 544 const T* v = a_mtx.data();
545 a_out << v[0] << "," << v[4] << "," << v[ 8] 545 a_out << v[0] << "," << v[4] << "," << v[ 8] << "," << v[12] << std::endl
546 << v[1] << "," << v[5] << "," << v[ 9] 546 << v[1] << "," << v[5] << "," << v[ 9] << "," << v[13] << std::endl
547 << v[2] << "," << v[6] << "," << v[10] 547 << v[2] << "," << v[6] << "," << v[10] << "," << v[14] << std::endl
548 << v[3] << "," << v[7] << "," << v[11] 548 << v[3] << "," << v[7] << "," << v[11] << "," << v[15] << std::endl;
549 return a_out; 549 return a_out;
550 } 550 }
551 551
552 } 552 }
553 553
554 #endif 554 #endif