Header file clean-up:
[mesa.git] / src / mesa / math / m_matrix.c
1 /* $Id: m_matrix.c,v 1.14 2002/10/24 23:57:24 brianp Exp $ */
2
3 /*
4 * Mesa 3-D graphics library
5 * Version: 4.1
6 *
7 * Copyright (C) 1999-2002 Brian Paul All Rights Reserved.
8 *
9 * Permission is hereby granted, free of charge, to any person obtaining a
10 * copy of this software and associated documentation files (the "Software"),
11 * to deal in the Software without restriction, including without limitation
12 * the rights to use, copy, modify, merge, publish, distribute, sublicense,
13 * and/or sell copies of the Software, and to permit persons to whom the
14 * Software is furnished to do so, subject to the following conditions:
15 *
16 * The above copyright notice and this permission notice shall be included
17 * in all copies or substantial portions of the Software.
18 *
19 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
20 * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
21 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
22 * BRIAN PAUL BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN
23 * AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
24 * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
25 */
26
27
28 /*
29 * Matrix operations
30 *
31 * NOTES:
32 * 1. 4x4 transformation matrices are stored in memory in column major order.
33 * 2. Points/vertices are to be thought of as column vectors.
34 * 3. Transformation of a point p by a matrix M is: p' = M * p
35 */
36
37 #include "glheader.h"
38 #include "imports.h"
39 #include "macros.h"
40 #include "imports.h"
41 #include "mmath.h"
42
43 #include "m_matrix.h"
44
45
46 static const char *types[] = {
47 "MATRIX_GENERAL",
48 "MATRIX_IDENTITY",
49 "MATRIX_3D_NO_ROT",
50 "MATRIX_PERSPECTIVE",
51 "MATRIX_2D",
52 "MATRIX_2D_NO_ROT",
53 "MATRIX_3D"
54 };
55
56
57 static GLfloat Identity[16] = {
58 1.0, 0.0, 0.0, 0.0,
59 0.0, 1.0, 0.0, 0.0,
60 0.0, 0.0, 1.0, 0.0,
61 0.0, 0.0, 0.0, 1.0
62 };
63
64
65
66
67 /*
68 * This matmul was contributed by Thomas Malik
69 *
70 * Perform a 4x4 matrix multiplication (product = a x b).
71 * Input: a, b - matrices to multiply
72 * Output: product - product of a and b
73 * WARNING: (product != b) assumed
74 * NOTE: (product == a) allowed
75 *
76 * KW: 4*16 = 64 muls
77 */
78 #define A(row,col) a[(col<<2)+row]
79 #define B(row,col) b[(col<<2)+row]
80 #define P(row,col) product[(col<<2)+row]
81
82 static void matmul4( GLfloat *product, const GLfloat *a, const GLfloat *b )
83 {
84 GLint i;
85 for (i = 0; i < 4; i++) {
86 const GLfloat ai0=A(i,0), ai1=A(i,1), ai2=A(i,2), ai3=A(i,3);
87 P(i,0) = ai0 * B(0,0) + ai1 * B(1,0) + ai2 * B(2,0) + ai3 * B(3,0);
88 P(i,1) = ai0 * B(0,1) + ai1 * B(1,1) + ai2 * B(2,1) + ai3 * B(3,1);
89 P(i,2) = ai0 * B(0,2) + ai1 * B(1,2) + ai2 * B(2,2) + ai3 * B(3,2);
90 P(i,3) = ai0 * B(0,3) + ai1 * B(1,3) + ai2 * B(2,3) + ai3 * B(3,3);
91 }
92 }
93
94
95 /* Multiply two matrices known to occupy only the top three rows, such
96 * as typical model matrices, and ortho matrices.
97 */
98 static void matmul34( GLfloat *product, const GLfloat *a, const GLfloat *b )
99 {
100 GLint i;
101 for (i = 0; i < 3; i++) {
102 const GLfloat ai0=A(i,0), ai1=A(i,1), ai2=A(i,2), ai3=A(i,3);
103 P(i,0) = ai0 * B(0,0) + ai1 * B(1,0) + ai2 * B(2,0);
104 P(i,1) = ai0 * B(0,1) + ai1 * B(1,1) + ai2 * B(2,1);
105 P(i,2) = ai0 * B(0,2) + ai1 * B(1,2) + ai2 * B(2,2);
106 P(i,3) = ai0 * B(0,3) + ai1 * B(1,3) + ai2 * B(2,3) + ai3;
107 }
108 P(3,0) = 0;
109 P(3,1) = 0;
110 P(3,2) = 0;
111 P(3,3) = 1;
112 }
113
114
115 #undef A
116 #undef B
117 #undef P
118
119
120 /*
121 * Multiply a matrix by an array of floats with known properties.
122 */
123 static void matrix_multf( GLmatrix *mat, const GLfloat *m, GLuint flags )
124 {
125 mat->flags |= (flags | MAT_DIRTY_TYPE | MAT_DIRTY_INVERSE);
126
127 if (TEST_MAT_FLAGS(mat, MAT_FLAGS_3D))
128 matmul34( mat->m, mat->m, m );
129 else
130 matmul4( mat->m, mat->m, m );
131 }
132
133
134 static void print_matrix_floats( const GLfloat m[16] )
135 {
136 int i;
137 for (i=0;i<4;i++) {
138 _mesa_debug(NULL,"\t%f %f %f %f\n", m[i], m[4+i], m[8+i], m[12+i] );
139 }
140 }
141
142 void
143 _math_matrix_print( const GLmatrix *m )
144 {
145 _mesa_debug(NULL, "Matrix type: %s, flags: %x\n", types[m->type], m->flags);
146 print_matrix_floats(m->m);
147 _mesa_debug(NULL, "Inverse: \n");
148 if (m->inv) {
149 GLfloat prod[16];
150 print_matrix_floats(m->inv);
151 matmul4(prod, m->m, m->inv);
152 _mesa_debug(NULL, "Mat * Inverse:\n");
153 print_matrix_floats(prod);
154 }
155 else {
156 _mesa_debug(NULL, " - not available\n");
157 }
158 }
159
160
161
162
163 #define SWAP_ROWS(a, b) { GLfloat *_tmp = a; (a)=(b); (b)=_tmp; }
164 #define MAT(m,r,c) (m)[(c)*4+(r)]
165
166 /*
167 * Compute inverse of 4x4 transformation matrix.
168 * Code contributed by Jacques Leroy jle@star.be
169 * Return GL_TRUE for success, GL_FALSE for failure (singular matrix)
170 */
171 static GLboolean invert_matrix_general( GLmatrix *mat )
172 {
173 const GLfloat *m = mat->m;
174 GLfloat *out = mat->inv;
175 GLfloat wtmp[4][8];
176 GLfloat m0, m1, m2, m3, s;
177 GLfloat *r0, *r1, *r2, *r3;
178
179 r0 = wtmp[0], r1 = wtmp[1], r2 = wtmp[2], r3 = wtmp[3];
180
181 r0[0] = MAT(m,0,0), r0[1] = MAT(m,0,1),
182 r0[2] = MAT(m,0,2), r0[3] = MAT(m,0,3),
183 r0[4] = 1.0, r0[5] = r0[6] = r0[7] = 0.0,
184
185 r1[0] = MAT(m,1,0), r1[1] = MAT(m,1,1),
186 r1[2] = MAT(m,1,2), r1[3] = MAT(m,1,3),
187 r1[5] = 1.0, r1[4] = r1[6] = r1[7] = 0.0,
188
189 r2[0] = MAT(m,2,0), r2[1] = MAT(m,2,1),
190 r2[2] = MAT(m,2,2), r2[3] = MAT(m,2,3),
191 r2[6] = 1.0, r2[4] = r2[5] = r2[7] = 0.0,
192
193 r3[0] = MAT(m,3,0), r3[1] = MAT(m,3,1),
194 r3[2] = MAT(m,3,2), r3[3] = MAT(m,3,3),
195 r3[7] = 1.0, r3[4] = r3[5] = r3[6] = 0.0;
196
197 /* choose pivot - or die */
198 if (fabs(r3[0])>fabs(r2[0])) SWAP_ROWS(r3, r2);
199 if (fabs(r2[0])>fabs(r1[0])) SWAP_ROWS(r2, r1);
200 if (fabs(r1[0])>fabs(r0[0])) SWAP_ROWS(r1, r0);
201 if (0.0 == r0[0]) return GL_FALSE;
202
203 /* eliminate first variable */
204 m1 = r1[0]/r0[0]; m2 = r2[0]/r0[0]; m3 = r3[0]/r0[0];
205 s = r0[1]; r1[1] -= m1 * s; r2[1] -= m2 * s; r3[1] -= m3 * s;
206 s = r0[2]; r1[2] -= m1 * s; r2[2] -= m2 * s; r3[2] -= m3 * s;
207 s = r0[3]; r1[3] -= m1 * s; r2[3] -= m2 * s; r3[3] -= m3 * s;
208 s = r0[4];
209 if (s != 0.0) { r1[4] -= m1 * s; r2[4] -= m2 * s; r3[4] -= m3 * s; }
210 s = r0[5];
211 if (s != 0.0) { r1[5] -= m1 * s; r2[5] -= m2 * s; r3[5] -= m3 * s; }
212 s = r0[6];
213 if (s != 0.0) { r1[6] -= m1 * s; r2[6] -= m2 * s; r3[6] -= m3 * s; }
214 s = r0[7];
215 if (s != 0.0) { r1[7] -= m1 * s; r2[7] -= m2 * s; r3[7] -= m3 * s; }
216
217 /* choose pivot - or die */
218 if (fabs(r3[1])>fabs(r2[1])) SWAP_ROWS(r3, r2);
219 if (fabs(r2[1])>fabs(r1[1])) SWAP_ROWS(r2, r1);
220 if (0.0 == r1[1]) return GL_FALSE;
221
222 /* eliminate second variable */
223 m2 = r2[1]/r1[1]; m3 = r3[1]/r1[1];
224 r2[2] -= m2 * r1[2]; r3[2] -= m3 * r1[2];
225 r2[3] -= m2 * r1[3]; r3[3] -= m3 * r1[3];
226 s = r1[4]; if (0.0 != s) { r2[4] -= m2 * s; r3[4] -= m3 * s; }
227 s = r1[5]; if (0.0 != s) { r2[5] -= m2 * s; r3[5] -= m3 * s; }
228 s = r1[6]; if (0.0 != s) { r2[6] -= m2 * s; r3[6] -= m3 * s; }
229 s = r1[7]; if (0.0 != s) { r2[7] -= m2 * s; r3[7] -= m3 * s; }
230
231 /* choose pivot - or die */
232 if (fabs(r3[2])>fabs(r2[2])) SWAP_ROWS(r3, r2);
233 if (0.0 == r2[2]) return GL_FALSE;
234
235 /* eliminate third variable */
236 m3 = r3[2]/r2[2];
237 r3[3] -= m3 * r2[3], r3[4] -= m3 * r2[4],
238 r3[5] -= m3 * r2[5], r3[6] -= m3 * r2[6],
239 r3[7] -= m3 * r2[7];
240
241 /* last check */
242 if (0.0 == r3[3]) return GL_FALSE;
243
244 s = 1.0F/r3[3]; /* now back substitute row 3 */
245 r3[4] *= s; r3[5] *= s; r3[6] *= s; r3[7] *= s;
246
247 m2 = r2[3]; /* now back substitute row 2 */
248 s = 1.0F/r2[2];
249 r2[4] = s * (r2[4] - r3[4] * m2), r2[5] = s * (r2[5] - r3[5] * m2),
250 r2[6] = s * (r2[6] - r3[6] * m2), r2[7] = s * (r2[7] - r3[7] * m2);
251 m1 = r1[3];
252 r1[4] -= r3[4] * m1, r1[5] -= r3[5] * m1,
253 r1[6] -= r3[6] * m1, r1[7] -= r3[7] * m1;
254 m0 = r0[3];
255 r0[4] -= r3[4] * m0, r0[5] -= r3[5] * m0,
256 r0[6] -= r3[6] * m0, r0[7] -= r3[7] * m0;
257
258 m1 = r1[2]; /* now back substitute row 1 */
259 s = 1.0F/r1[1];
260 r1[4] = s * (r1[4] - r2[4] * m1), r1[5] = s * (r1[5] - r2[5] * m1),
261 r1[6] = s * (r1[6] - r2[6] * m1), r1[7] = s * (r1[7] - r2[7] * m1);
262 m0 = r0[2];
263 r0[4] -= r2[4] * m0, r0[5] -= r2[5] * m0,
264 r0[6] -= r2[6] * m0, r0[7] -= r2[7] * m0;
265
266 m0 = r0[1]; /* now back substitute row 0 */
267 s = 1.0F/r0[0];
268 r0[4] = s * (r0[4] - r1[4] * m0), r0[5] = s * (r0[5] - r1[5] * m0),
269 r0[6] = s * (r0[6] - r1[6] * m0), r0[7] = s * (r0[7] - r1[7] * m0);
270
271 MAT(out,0,0) = r0[4]; MAT(out,0,1) = r0[5],
272 MAT(out,0,2) = r0[6]; MAT(out,0,3) = r0[7],
273 MAT(out,1,0) = r1[4]; MAT(out,1,1) = r1[5],
274 MAT(out,1,2) = r1[6]; MAT(out,1,3) = r1[7],
275 MAT(out,2,0) = r2[4]; MAT(out,2,1) = r2[5],
276 MAT(out,2,2) = r2[6]; MAT(out,2,3) = r2[7],
277 MAT(out,3,0) = r3[4]; MAT(out,3,1) = r3[5],
278 MAT(out,3,2) = r3[6]; MAT(out,3,3) = r3[7];
279
280 return GL_TRUE;
281 }
282 #undef SWAP_ROWS
283
284
285 /* Adapted from graphics gems II.
286 */
287 static GLboolean invert_matrix_3d_general( GLmatrix *mat )
288 {
289 const GLfloat *in = mat->m;
290 GLfloat *out = mat->inv;
291 GLfloat pos, neg, t;
292 GLfloat det;
293
294 /* Calculate the determinant of upper left 3x3 submatrix and
295 * determine if the matrix is singular.
296 */
297 pos = neg = 0.0;
298 t = MAT(in,0,0) * MAT(in,1,1) * MAT(in,2,2);
299 if (t >= 0.0) pos += t; else neg += t;
300
301 t = MAT(in,1,0) * MAT(in,2,1) * MAT(in,0,2);
302 if (t >= 0.0) pos += t; else neg += t;
303
304 t = MAT(in,2,0) * MAT(in,0,1) * MAT(in,1,2);
305 if (t >= 0.0) pos += t; else neg += t;
306
307 t = -MAT(in,2,0) * MAT(in,1,1) * MAT(in,0,2);
308 if (t >= 0.0) pos += t; else neg += t;
309
310 t = -MAT(in,1,0) * MAT(in,0,1) * MAT(in,2,2);
311 if (t >= 0.0) pos += t; else neg += t;
312
313 t = -MAT(in,0,0) * MAT(in,2,1) * MAT(in,1,2);
314 if (t >= 0.0) pos += t; else neg += t;
315
316 det = pos + neg;
317
318 if (det*det < 1e-25)
319 return GL_FALSE;
320
321 det = 1.0F / det;
322 MAT(out,0,0) = ( (MAT(in,1,1)*MAT(in,2,2) - MAT(in,2,1)*MAT(in,1,2) )*det);
323 MAT(out,0,1) = (- (MAT(in,0,1)*MAT(in,2,2) - MAT(in,2,1)*MAT(in,0,2) )*det);
324 MAT(out,0,2) = ( (MAT(in,0,1)*MAT(in,1,2) - MAT(in,1,1)*MAT(in,0,2) )*det);
325 MAT(out,1,0) = (- (MAT(in,1,0)*MAT(in,2,2) - MAT(in,2,0)*MAT(in,1,2) )*det);
326 MAT(out,1,1) = ( (MAT(in,0,0)*MAT(in,2,2) - MAT(in,2,0)*MAT(in,0,2) )*det);
327 MAT(out,1,2) = (- (MAT(in,0,0)*MAT(in,1,2) - MAT(in,1,0)*MAT(in,0,2) )*det);
328 MAT(out,2,0) = ( (MAT(in,1,0)*MAT(in,2,1) - MAT(in,2,0)*MAT(in,1,1) )*det);
329 MAT(out,2,1) = (- (MAT(in,0,0)*MAT(in,2,1) - MAT(in,2,0)*MAT(in,0,1) )*det);
330 MAT(out,2,2) = ( (MAT(in,0,0)*MAT(in,1,1) - MAT(in,1,0)*MAT(in,0,1) )*det);
331
332 /* Do the translation part */
333 MAT(out,0,3) = - (MAT(in,0,3) * MAT(out,0,0) +
334 MAT(in,1,3) * MAT(out,0,1) +
335 MAT(in,2,3) * MAT(out,0,2) );
336 MAT(out,1,3) = - (MAT(in,0,3) * MAT(out,1,0) +
337 MAT(in,1,3) * MAT(out,1,1) +
338 MAT(in,2,3) * MAT(out,1,2) );
339 MAT(out,2,3) = - (MAT(in,0,3) * MAT(out,2,0) +
340 MAT(in,1,3) * MAT(out,2,1) +
341 MAT(in,2,3) * MAT(out,2,2) );
342
343 return GL_TRUE;
344 }
345
346
347 static GLboolean invert_matrix_3d( GLmatrix *mat )
348 {
349 const GLfloat *in = mat->m;
350 GLfloat *out = mat->inv;
351
352 if (!TEST_MAT_FLAGS(mat, MAT_FLAGS_ANGLE_PRESERVING)) {
353 return invert_matrix_3d_general( mat );
354 }
355
356 if (mat->flags & MAT_FLAG_UNIFORM_SCALE) {
357 GLfloat scale = (MAT(in,0,0) * MAT(in,0,0) +
358 MAT(in,0,1) * MAT(in,0,1) +
359 MAT(in,0,2) * MAT(in,0,2));
360
361 if (scale == 0.0)
362 return GL_FALSE;
363
364 scale = 1.0F / scale;
365
366 /* Transpose and scale the 3 by 3 upper-left submatrix. */
367 MAT(out,0,0) = scale * MAT(in,0,0);
368 MAT(out,1,0) = scale * MAT(in,0,1);
369 MAT(out,2,0) = scale * MAT(in,0,2);
370 MAT(out,0,1) = scale * MAT(in,1,0);
371 MAT(out,1,1) = scale * MAT(in,1,1);
372 MAT(out,2,1) = scale * MAT(in,1,2);
373 MAT(out,0,2) = scale * MAT(in,2,0);
374 MAT(out,1,2) = scale * MAT(in,2,1);
375 MAT(out,2,2) = scale * MAT(in,2,2);
376 }
377 else if (mat->flags & MAT_FLAG_ROTATION) {
378 /* Transpose the 3 by 3 upper-left submatrix. */
379 MAT(out,0,0) = MAT(in,0,0);
380 MAT(out,1,0) = MAT(in,0,1);
381 MAT(out,2,0) = MAT(in,0,2);
382 MAT(out,0,1) = MAT(in,1,0);
383 MAT(out,1,1) = MAT(in,1,1);
384 MAT(out,2,1) = MAT(in,1,2);
385 MAT(out,0,2) = MAT(in,2,0);
386 MAT(out,1,2) = MAT(in,2,1);
387 MAT(out,2,2) = MAT(in,2,2);
388 }
389 else {
390 /* pure translation */
391 MEMCPY( out, Identity, sizeof(Identity) );
392 MAT(out,0,3) = - MAT(in,0,3);
393 MAT(out,1,3) = - MAT(in,1,3);
394 MAT(out,2,3) = - MAT(in,2,3);
395 return GL_TRUE;
396 }
397
398 if (mat->flags & MAT_FLAG_TRANSLATION) {
399 /* Do the translation part */
400 MAT(out,0,3) = - (MAT(in,0,3) * MAT(out,0,0) +
401 MAT(in,1,3) * MAT(out,0,1) +
402 MAT(in,2,3) * MAT(out,0,2) );
403 MAT(out,1,3) = - (MAT(in,0,3) * MAT(out,1,0) +
404 MAT(in,1,3) * MAT(out,1,1) +
405 MAT(in,2,3) * MAT(out,1,2) );
406 MAT(out,2,3) = - (MAT(in,0,3) * MAT(out,2,0) +
407 MAT(in,1,3) * MAT(out,2,1) +
408 MAT(in,2,3) * MAT(out,2,2) );
409 }
410 else {
411 MAT(out,0,3) = MAT(out,1,3) = MAT(out,2,3) = 0.0;
412 }
413
414 return GL_TRUE;
415 }
416
417
418
419 static GLboolean invert_matrix_identity( GLmatrix *mat )
420 {
421 MEMCPY( mat->inv, Identity, sizeof(Identity) );
422 return GL_TRUE;
423 }
424
425
426 static GLboolean invert_matrix_3d_no_rot( GLmatrix *mat )
427 {
428 const GLfloat *in = mat->m;
429 GLfloat *out = mat->inv;
430
431 if (MAT(in,0,0) == 0 || MAT(in,1,1) == 0 || MAT(in,2,2) == 0 )
432 return GL_FALSE;
433
434 MEMCPY( out, Identity, 16 * sizeof(GLfloat) );
435 MAT(out,0,0) = 1.0F / MAT(in,0,0);
436 MAT(out,1,1) = 1.0F / MAT(in,1,1);
437 MAT(out,2,2) = 1.0F / MAT(in,2,2);
438
439 if (mat->flags & MAT_FLAG_TRANSLATION) {
440 MAT(out,0,3) = - (MAT(in,0,3) * MAT(out,0,0));
441 MAT(out,1,3) = - (MAT(in,1,3) * MAT(out,1,1));
442 MAT(out,2,3) = - (MAT(in,2,3) * MAT(out,2,2));
443 }
444
445 return GL_TRUE;
446 }
447
448
449 static GLboolean invert_matrix_2d_no_rot( GLmatrix *mat )
450 {
451 const GLfloat *in = mat->m;
452 GLfloat *out = mat->inv;
453
454 if (MAT(in,0,0) == 0 || MAT(in,1,1) == 0)
455 return GL_FALSE;
456
457 MEMCPY( out, Identity, 16 * sizeof(GLfloat) );
458 MAT(out,0,0) = 1.0F / MAT(in,0,0);
459 MAT(out,1,1) = 1.0F / MAT(in,1,1);
460
461 if (mat->flags & MAT_FLAG_TRANSLATION) {
462 MAT(out,0,3) = - (MAT(in,0,3) * MAT(out,0,0));
463 MAT(out,1,3) = - (MAT(in,1,3) * MAT(out,1,1));
464 }
465
466 return GL_TRUE;
467 }
468
469
470 #if 0
471 /* broken */
472 static GLboolean invert_matrix_perspective( GLmatrix *mat )
473 {
474 const GLfloat *in = mat->m;
475 GLfloat *out = mat->inv;
476
477 if (MAT(in,2,3) == 0)
478 return GL_FALSE;
479
480 MEMCPY( out, Identity, 16 * sizeof(GLfloat) );
481
482 MAT(out,0,0) = 1.0F / MAT(in,0,0);
483 MAT(out,1,1) = 1.0F / MAT(in,1,1);
484
485 MAT(out,0,3) = MAT(in,0,2);
486 MAT(out,1,3) = MAT(in,1,2);
487
488 MAT(out,2,2) = 0;
489 MAT(out,2,3) = -1;
490
491 MAT(out,3,2) = 1.0F / MAT(in,2,3);
492 MAT(out,3,3) = MAT(in,2,2) * MAT(out,3,2);
493
494 return GL_TRUE;
495 }
496 #endif
497
498
499 typedef GLboolean (*inv_mat_func)( GLmatrix *mat );
500
501
502 static inv_mat_func inv_mat_tab[7] = {
503 invert_matrix_general,
504 invert_matrix_identity,
505 invert_matrix_3d_no_rot,
506 #if 0
507 /* Don't use this function for now - it fails when the projection matrix
508 * is premultiplied by a translation (ala Chromium's tilesort SPU).
509 */
510 invert_matrix_perspective,
511 #else
512 invert_matrix_general,
513 #endif
514 invert_matrix_3d, /* lazy! */
515 invert_matrix_2d_no_rot,
516 invert_matrix_3d
517 };
518
519
520 static GLboolean matrix_invert( GLmatrix *mat )
521 {
522 if (inv_mat_tab[mat->type](mat)) {
523 mat->flags &= ~MAT_FLAG_SINGULAR;
524 return GL_TRUE;
525 } else {
526 mat->flags |= MAT_FLAG_SINGULAR;
527 MEMCPY( mat->inv, Identity, sizeof(Identity) );
528 return GL_FALSE;
529 }
530 }
531
532
533
534
535
536
537 /*
538 * Generate a 4x4 transformation matrix from glRotate parameters, and
539 * postmultiply the input matrix by it.
540 * This function contributed by Erich Boleyn (erich@uruk.org).
541 * Optimizatios contributed by Rudolf Opalla (rudi@khm.de).
542 */
543 void
544 _math_matrix_rotate( GLmatrix *mat,
545 GLfloat angle, GLfloat x, GLfloat y, GLfloat z )
546 {
547 GLfloat xx, yy, zz, xy, yz, zx, xs, ys, zs, one_c, s, c;
548 GLfloat m[16];
549 GLboolean optimized;
550
551 s = (GLfloat) sin( angle * DEG2RAD );
552 c = (GLfloat) cos( angle * DEG2RAD );
553
554 MEMCPY(m, Identity, sizeof(GLfloat)*16);
555 optimized = GL_FALSE;
556
557 #define M(row,col) m[col*4+row]
558
559 if (x == 0.0F) {
560 if (y == 0.0F) {
561 if (z != 0.0F) {
562 optimized = GL_TRUE;
563 /* rotate only around z-axis */
564 M(0,0) = c;
565 M(1,1) = c;
566 if (z < 0.0F) {
567 M(0,1) = s;
568 M(1,0) = -s;
569 }
570 else {
571 M(0,1) = -s;
572 M(1,0) = s;
573 }
574 }
575 }
576 else if (z == 0.0F) {
577 optimized = GL_TRUE;
578 /* rotate only around y-axis */
579 M(0,0) = c;
580 M(2,2) = c;
581 if (y < 0.0F) {
582 M(0,2) = -s;
583 M(2,0) = s;
584 }
585 else {
586 M(0,2) = s;
587 M(2,0) = -s;
588 }
589 }
590 }
591 else if (y == 0.0F) {
592 if (z == 0.0F) {
593 optimized = GL_TRUE;
594 /* rotate only around x-axis */
595 M(1,1) = c;
596 M(2,2) = c;
597 if (y < 0.0F) {
598 M(1,2) = s;
599 M(2,1) = -s;
600 }
601 else {
602 M(1,2) = -s;
603 M(2,1) = s;
604 }
605 }
606 }
607
608 if (!optimized) {
609 const GLfloat mag = (GLfloat) GL_SQRT(x * x + y * y + z * z);
610
611 if (mag <= 1.0e-4) {
612 /* no rotation, leave mat as-is */
613 return;
614 }
615
616 x /= mag;
617 y /= mag;
618 z /= mag;
619
620
621 /*
622 * Arbitrary axis rotation matrix.
623 *
624 * This is composed of 5 matrices, Rz, Ry, T, Ry', Rz', multiplied
625 * like so: Rz * Ry * T * Ry' * Rz'. T is the final rotation
626 * (which is about the X-axis), and the two composite transforms
627 * Ry' * Rz' and Rz * Ry are (respectively) the rotations necessary
628 * from the arbitrary axis to the X-axis then back. They are
629 * all elementary rotations.
630 *
631 * Rz' is a rotation about the Z-axis, to bring the axis vector
632 * into the x-z plane. Then Ry' is applied, rotating about the
633 * Y-axis to bring the axis vector parallel with the X-axis. The
634 * rotation about the X-axis is then performed. Ry and Rz are
635 * simply the respective inverse transforms to bring the arbitrary
636 * axis back to it's original orientation. The first transforms
637 * Rz' and Ry' are considered inverses, since the data from the
638 * arbitrary axis gives you info on how to get to it, not how
639 * to get away from it, and an inverse must be applied.
640 *
641 * The basic calculation used is to recognize that the arbitrary
642 * axis vector (x, y, z), since it is of unit length, actually
643 * represents the sines and cosines of the angles to rotate the
644 * X-axis to the same orientation, with theta being the angle about
645 * Z and phi the angle about Y (in the order described above)
646 * as follows:
647 *
648 * cos ( theta ) = x / sqrt ( 1 - z^2 )
649 * sin ( theta ) = y / sqrt ( 1 - z^2 )
650 *
651 * cos ( phi ) = sqrt ( 1 - z^2 )
652 * sin ( phi ) = z
653 *
654 * Note that cos ( phi ) can further be inserted to the above
655 * formulas:
656 *
657 * cos ( theta ) = x / cos ( phi )
658 * sin ( theta ) = y / sin ( phi )
659 *
660 * ...etc. Because of those relations and the standard trigonometric
661 * relations, it is pssible to reduce the transforms down to what
662 * is used below. It may be that any primary axis chosen will give the
663 * same results (modulo a sign convention) using thie method.
664 *
665 * Particularly nice is to notice that all divisions that might
666 * have caused trouble when parallel to certain planes or
667 * axis go away with care paid to reducing the expressions.
668 * After checking, it does perform correctly under all cases, since
669 * in all the cases of division where the denominator would have
670 * been zero, the numerator would have been zero as well, giving
671 * the expected result.
672 */
673
674 xx = x * x;
675 yy = y * y;
676 zz = z * z;
677 xy = x * y;
678 yz = y * z;
679 zx = z * x;
680 xs = x * s;
681 ys = y * s;
682 zs = z * s;
683 one_c = 1.0F - c;
684
685 /* We already hold the identity-matrix so we can skip some statements */
686 M(0,0) = (one_c * xx) + c;
687 M(0,1) = (one_c * xy) - zs;
688 M(0,2) = (one_c * zx) + ys;
689 /* M(0,3) = 0.0F; */
690
691 M(1,0) = (one_c * xy) + zs;
692 M(1,1) = (one_c * yy) + c;
693 M(1,2) = (one_c * yz) - xs;
694 /* M(1,3) = 0.0F; */
695
696 M(2,0) = (one_c * zx) - ys;
697 M(2,1) = (one_c * yz) + xs;
698 M(2,2) = (one_c * zz) + c;
699 /* M(2,3) = 0.0F; */
700
701 /*
702 M(3,0) = 0.0F;
703 M(3,1) = 0.0F;
704 M(3,2) = 0.0F;
705 M(3,3) = 1.0F;
706 */
707 }
708 #undef M
709
710 matrix_multf( mat, m, MAT_FLAG_ROTATION );
711 }
712
713
714
715 void
716 _math_matrix_frustum( GLmatrix *mat,
717 GLfloat left, GLfloat right,
718 GLfloat bottom, GLfloat top,
719 GLfloat nearval, GLfloat farval )
720 {
721 GLfloat x, y, a, b, c, d;
722 GLfloat m[16];
723
724 x = (2.0F*nearval) / (right-left);
725 y = (2.0F*nearval) / (top-bottom);
726 a = (right+left) / (right-left);
727 b = (top+bottom) / (top-bottom);
728 c = -(farval+nearval) / ( farval-nearval);
729 d = -(2.0F*farval*nearval) / (farval-nearval); /* error? */
730
731 #define M(row,col) m[col*4+row]
732 M(0,0) = x; M(0,1) = 0.0F; M(0,2) = a; M(0,3) = 0.0F;
733 M(1,0) = 0.0F; M(1,1) = y; M(1,2) = b; M(1,3) = 0.0F;
734 M(2,0) = 0.0F; M(2,1) = 0.0F; M(2,2) = c; M(2,3) = d;
735 M(3,0) = 0.0F; M(3,1) = 0.0F; M(3,2) = -1.0F; M(3,3) = 0.0F;
736 #undef M
737
738 matrix_multf( mat, m, MAT_FLAG_PERSPECTIVE );
739 }
740
741 void
742 _math_matrix_ortho( GLmatrix *mat,
743 GLfloat left, GLfloat right,
744 GLfloat bottom, GLfloat top,
745 GLfloat nearval, GLfloat farval )
746 {
747 GLfloat x, y, z;
748 GLfloat tx, ty, tz;
749 GLfloat m[16];
750
751 x = 2.0F / (right-left);
752 y = 2.0F / (top-bottom);
753 z = -2.0F / (farval-nearval);
754 tx = -(right+left) / (right-left);
755 ty = -(top+bottom) / (top-bottom);
756 tz = -(farval+nearval) / (farval-nearval);
757
758 #define M(row,col) m[col*4+row]
759 M(0,0) = x; M(0,1) = 0.0F; M(0,2) = 0.0F; M(0,3) = tx;
760 M(1,0) = 0.0F; M(1,1) = y; M(1,2) = 0.0F; M(1,3) = ty;
761 M(2,0) = 0.0F; M(2,1) = 0.0F; M(2,2) = z; M(2,3) = tz;
762 M(3,0) = 0.0F; M(3,1) = 0.0F; M(3,2) = 0.0F; M(3,3) = 1.0F;
763 #undef M
764
765 matrix_multf( mat, m, (MAT_FLAG_GENERAL_SCALE|MAT_FLAG_TRANSLATION));
766 }
767
768
769 #define ZERO(x) (1<<x)
770 #define ONE(x) (1<<(x+16))
771
772 #define MASK_NO_TRX (ZERO(12) | ZERO(13) | ZERO(14))
773 #define MASK_NO_2D_SCALE ( ONE(0) | ONE(5))
774
775 #define MASK_IDENTITY ( ONE(0) | ZERO(4) | ZERO(8) | ZERO(12) |\
776 ZERO(1) | ONE(5) | ZERO(9) | ZERO(13) |\
777 ZERO(2) | ZERO(6) | ONE(10) | ZERO(14) |\
778 ZERO(3) | ZERO(7) | ZERO(11) | ONE(15) )
779
780 #define MASK_2D_NO_ROT ( ZERO(4) | ZERO(8) | \
781 ZERO(1) | ZERO(9) | \
782 ZERO(2) | ZERO(6) | ONE(10) | ZERO(14) |\
783 ZERO(3) | ZERO(7) | ZERO(11) | ONE(15) )
784
785 #define MASK_2D ( ZERO(8) | \
786 ZERO(9) | \
787 ZERO(2) | ZERO(6) | ONE(10) | ZERO(14) |\
788 ZERO(3) | ZERO(7) | ZERO(11) | ONE(15) )
789
790
791 #define MASK_3D_NO_ROT ( ZERO(4) | ZERO(8) | \
792 ZERO(1) | ZERO(9) | \
793 ZERO(2) | ZERO(6) | \
794 ZERO(3) | ZERO(7) | ZERO(11) | ONE(15) )
795
796 #define MASK_3D ( \
797 \
798 \
799 ZERO(3) | ZERO(7) | ZERO(11) | ONE(15) )
800
801
802 #define MASK_PERSPECTIVE ( ZERO(4) | ZERO(12) |\
803 ZERO(1) | ZERO(13) |\
804 ZERO(2) | ZERO(6) | \
805 ZERO(3) | ZERO(7) | ZERO(15) )
806
807 #define SQ(x) ((x)*(x))
808
809 /* Determine type and flags from scratch. This is expensive enough to
810 * only want to do it once.
811 */
812 static void analyse_from_scratch( GLmatrix *mat )
813 {
814 const GLfloat *m = mat->m;
815 GLuint mask = 0;
816 GLuint i;
817
818 for (i = 0 ; i < 16 ; i++) {
819 if (m[i] == 0.0) mask |= (1<<i);
820 }
821
822 if (m[0] == 1.0F) mask |= (1<<16);
823 if (m[5] == 1.0F) mask |= (1<<21);
824 if (m[10] == 1.0F) mask |= (1<<26);
825 if (m[15] == 1.0F) mask |= (1<<31);
826
827 mat->flags &= ~MAT_FLAGS_GEOMETRY;
828
829 /* Check for translation - no-one really cares
830 */
831 if ((mask & MASK_NO_TRX) != MASK_NO_TRX)
832 mat->flags |= MAT_FLAG_TRANSLATION;
833
834 /* Do the real work
835 */
836 if (mask == (GLuint) MASK_IDENTITY) {
837 mat->type = MATRIX_IDENTITY;
838 }
839 else if ((mask & MASK_2D_NO_ROT) == (GLuint) MASK_2D_NO_ROT) {
840 mat->type = MATRIX_2D_NO_ROT;
841
842 if ((mask & MASK_NO_2D_SCALE) != MASK_NO_2D_SCALE)
843 mat->flags = MAT_FLAG_GENERAL_SCALE;
844 }
845 else if ((mask & MASK_2D) == (GLuint) MASK_2D) {
846 GLfloat mm = DOT2(m, m);
847 GLfloat m4m4 = DOT2(m+4,m+4);
848 GLfloat mm4 = DOT2(m,m+4);
849
850 mat->type = MATRIX_2D;
851
852 /* Check for scale */
853 if (SQ(mm-1) > SQ(1e-6) ||
854 SQ(m4m4-1) > SQ(1e-6))
855 mat->flags |= MAT_FLAG_GENERAL_SCALE;
856
857 /* Check for rotation */
858 if (SQ(mm4) > SQ(1e-6))
859 mat->flags |= MAT_FLAG_GENERAL_3D;
860 else
861 mat->flags |= MAT_FLAG_ROTATION;
862
863 }
864 else if ((mask & MASK_3D_NO_ROT) == (GLuint) MASK_3D_NO_ROT) {
865 mat->type = MATRIX_3D_NO_ROT;
866
867 /* Check for scale */
868 if (SQ(m[0]-m[5]) < SQ(1e-6) &&
869 SQ(m[0]-m[10]) < SQ(1e-6)) {
870 if (SQ(m[0]-1.0) > SQ(1e-6)) {
871 mat->flags |= MAT_FLAG_UNIFORM_SCALE;
872 }
873 }
874 else {
875 mat->flags |= MAT_FLAG_GENERAL_SCALE;
876 }
877 }
878 else if ((mask & MASK_3D) == (GLuint) MASK_3D) {
879 GLfloat c1 = DOT3(m,m);
880 GLfloat c2 = DOT3(m+4,m+4);
881 GLfloat c3 = DOT3(m+8,m+8);
882 GLfloat d1 = DOT3(m, m+4);
883 GLfloat cp[3];
884
885 mat->type = MATRIX_3D;
886
887 /* Check for scale */
888 if (SQ(c1-c2) < SQ(1e-6) && SQ(c1-c3) < SQ(1e-6)) {
889 if (SQ(c1-1.0) > SQ(1e-6))
890 mat->flags |= MAT_FLAG_UNIFORM_SCALE;
891 /* else no scale at all */
892 }
893 else {
894 mat->flags |= MAT_FLAG_GENERAL_SCALE;
895 }
896
897 /* Check for rotation */
898 if (SQ(d1) < SQ(1e-6)) {
899 CROSS3( cp, m, m+4 );
900 SUB_3V( cp, cp, (m+8) );
901 if (LEN_SQUARED_3FV(cp) < SQ(1e-6))
902 mat->flags |= MAT_FLAG_ROTATION;
903 else
904 mat->flags |= MAT_FLAG_GENERAL_3D;
905 }
906 else {
907 mat->flags |= MAT_FLAG_GENERAL_3D; /* shear, etc */
908 }
909 }
910 else if ((mask & MASK_PERSPECTIVE) == MASK_PERSPECTIVE && m[11]==-1.0F) {
911 mat->type = MATRIX_PERSPECTIVE;
912 mat->flags |= MAT_FLAG_GENERAL;
913 }
914 else {
915 mat->type = MATRIX_GENERAL;
916 mat->flags |= MAT_FLAG_GENERAL;
917 }
918 }
919
920
921 /* Analyse a matrix given that its flags are accurate - this is the
922 * more common operation, hopefully.
923 */
924 static void analyse_from_flags( GLmatrix *mat )
925 {
926 const GLfloat *m = mat->m;
927
928 if (TEST_MAT_FLAGS(mat, 0)) {
929 mat->type = MATRIX_IDENTITY;
930 }
931 else if (TEST_MAT_FLAGS(mat, (MAT_FLAG_TRANSLATION |
932 MAT_FLAG_UNIFORM_SCALE |
933 MAT_FLAG_GENERAL_SCALE))) {
934 if ( m[10]==1.0F && m[14]==0.0F ) {
935 mat->type = MATRIX_2D_NO_ROT;
936 }
937 else {
938 mat->type = MATRIX_3D_NO_ROT;
939 }
940 }
941 else if (TEST_MAT_FLAGS(mat, MAT_FLAGS_3D)) {
942 if ( m[ 8]==0.0F
943 && m[ 9]==0.0F
944 && m[2]==0.0F && m[6]==0.0F && m[10]==1.0F && m[14]==0.0F) {
945 mat->type = MATRIX_2D;
946 }
947 else {
948 mat->type = MATRIX_3D;
949 }
950 }
951 else if ( m[4]==0.0F && m[12]==0.0F
952 && m[1]==0.0F && m[13]==0.0F
953 && m[2]==0.0F && m[6]==0.0F
954 && m[3]==0.0F && m[7]==0.0F && m[11]==-1.0F && m[15]==0.0F) {
955 mat->type = MATRIX_PERSPECTIVE;
956 }
957 else {
958 mat->type = MATRIX_GENERAL;
959 }
960 }
961
962
963 void
964 _math_matrix_analyse( GLmatrix *mat )
965 {
966 if (mat->flags & MAT_DIRTY_TYPE) {
967 if (mat->flags & MAT_DIRTY_FLAGS)
968 analyse_from_scratch( mat );
969 else
970 analyse_from_flags( mat );
971 }
972
973 if (mat->inv && (mat->flags & MAT_DIRTY_INVERSE)) {
974 matrix_invert( mat );
975 }
976
977 mat->flags &= ~(MAT_DIRTY_FLAGS|
978 MAT_DIRTY_TYPE|
979 MAT_DIRTY_INVERSE);
980 }
981
982
983 void
984 _math_matrix_copy( GLmatrix *to, const GLmatrix *from )
985 {
986 MEMCPY( to->m, from->m, sizeof(Identity) );
987 to->flags = from->flags;
988 to->type = from->type;
989
990 if (to->inv != 0) {
991 if (from->inv == 0) {
992 matrix_invert( to );
993 }
994 else {
995 MEMCPY(to->inv, from->inv, sizeof(GLfloat)*16);
996 }
997 }
998 }
999
1000
1001 void
1002 _math_matrix_scale( GLmatrix *mat, GLfloat x, GLfloat y, GLfloat z )
1003 {
1004 GLfloat *m = mat->m;
1005 m[0] *= x; m[4] *= y; m[8] *= z;
1006 m[1] *= x; m[5] *= y; m[9] *= z;
1007 m[2] *= x; m[6] *= y; m[10] *= z;
1008 m[3] *= x; m[7] *= y; m[11] *= z;
1009
1010 if (fabs(x - y) < 1e-8 && fabs(x - z) < 1e-8)
1011 mat->flags |= MAT_FLAG_UNIFORM_SCALE;
1012 else
1013 mat->flags |= MAT_FLAG_GENERAL_SCALE;
1014
1015 mat->flags |= (MAT_DIRTY_TYPE |
1016 MAT_DIRTY_INVERSE);
1017 }
1018
1019
1020 void
1021 _math_matrix_translate( GLmatrix *mat, GLfloat x, GLfloat y, GLfloat z )
1022 {
1023 GLfloat *m = mat->m;
1024 m[12] = m[0] * x + m[4] * y + m[8] * z + m[12];
1025 m[13] = m[1] * x + m[5] * y + m[9] * z + m[13];
1026 m[14] = m[2] * x + m[6] * y + m[10] * z + m[14];
1027 m[15] = m[3] * x + m[7] * y + m[11] * z + m[15];
1028
1029 mat->flags |= (MAT_FLAG_TRANSLATION |
1030 MAT_DIRTY_TYPE |
1031 MAT_DIRTY_INVERSE);
1032 }
1033
1034
1035 void
1036 _math_matrix_loadf( GLmatrix *mat, const GLfloat *m )
1037 {
1038 MEMCPY( mat->m, m, 16*sizeof(GLfloat) );
1039 mat->flags = (MAT_FLAG_GENERAL | MAT_DIRTY);
1040 }
1041
1042 void
1043 _math_matrix_ctr( GLmatrix *m )
1044 {
1045 m->m = (GLfloat *) ALIGN_MALLOC( 16 * sizeof(GLfloat), 16 );
1046 if (m->m)
1047 MEMCPY( m->m, Identity, sizeof(Identity) );
1048 m->inv = NULL;
1049 m->type = MATRIX_IDENTITY;
1050 m->flags = 0;
1051 }
1052
1053 void
1054 _math_matrix_dtr( GLmatrix *m )
1055 {
1056 if (m->m) {
1057 ALIGN_FREE( m->m );
1058 m->m = NULL;
1059 }
1060 if (m->inv) {
1061 ALIGN_FREE( m->inv );
1062 m->inv = NULL;
1063 }
1064 }
1065
1066
1067 void
1068 _math_matrix_alloc_inv( GLmatrix *m )
1069 {
1070 if (!m->inv) {
1071 m->inv = (GLfloat *) ALIGN_MALLOC( 16 * sizeof(GLfloat), 16 );
1072 if (m->inv)
1073 MEMCPY( m->inv, Identity, 16 * sizeof(GLfloat) );
1074 }
1075 }
1076
1077
1078 void
1079 _math_matrix_mul_matrix( GLmatrix *dest, const GLmatrix *a, const GLmatrix *b )
1080 {
1081 dest->flags = (a->flags |
1082 b->flags |
1083 MAT_DIRTY_TYPE |
1084 MAT_DIRTY_INVERSE);
1085
1086 if (TEST_MAT_FLAGS(dest, MAT_FLAGS_3D))
1087 matmul34( dest->m, a->m, b->m );
1088 else
1089 matmul4( dest->m, a->m, b->m );
1090 }
1091
1092
1093 void
1094 _math_matrix_mul_floats( GLmatrix *dest, const GLfloat *m )
1095 {
1096 dest->flags |= (MAT_FLAG_GENERAL |
1097 MAT_DIRTY_TYPE |
1098 MAT_DIRTY_INVERSE);
1099
1100 matmul4( dest->m, dest->m, m );
1101 }
1102
1103 void
1104 _math_matrix_set_identity( GLmatrix *mat )
1105 {
1106 MEMCPY( mat->m, Identity, 16*sizeof(GLfloat) );
1107
1108 if (mat->inv)
1109 MEMCPY( mat->inv, Identity, 16*sizeof(GLfloat) );
1110
1111 mat->type = MATRIX_IDENTITY;
1112 mat->flags &= ~(MAT_DIRTY_FLAGS|
1113 MAT_DIRTY_TYPE|
1114 MAT_DIRTY_INVERSE);
1115 }
1116
1117
1118
1119 void
1120 _math_transposef( GLfloat to[16], const GLfloat from[16] )
1121 {
1122 to[0] = from[0];
1123 to[1] = from[4];
1124 to[2] = from[8];
1125 to[3] = from[12];
1126 to[4] = from[1];
1127 to[5] = from[5];
1128 to[6] = from[9];
1129 to[7] = from[13];
1130 to[8] = from[2];
1131 to[9] = from[6];
1132 to[10] = from[10];
1133 to[11] = from[14];
1134 to[12] = from[3];
1135 to[13] = from[7];
1136 to[14] = from[11];
1137 to[15] = from[15];
1138 }
1139
1140
1141 void
1142 _math_transposed( GLdouble to[16], const GLdouble from[16] )
1143 {
1144 to[0] = from[0];
1145 to[1] = from[4];
1146 to[2] = from[8];
1147 to[3] = from[12];
1148 to[4] = from[1];
1149 to[5] = from[5];
1150 to[6] = from[9];
1151 to[7] = from[13];
1152 to[8] = from[2];
1153 to[9] = from[6];
1154 to[10] = from[10];
1155 to[11] = from[14];
1156 to[12] = from[3];
1157 to[13] = from[7];
1158 to[14] = from[11];
1159 to[15] = from[15];
1160 }
1161
1162 void
1163 _math_transposefd( GLfloat to[16], const GLdouble from[16] )
1164 {
1165 to[0] = (GLfloat) from[0];
1166 to[1] = (GLfloat) from[4];
1167 to[2] = (GLfloat) from[8];
1168 to[3] = (GLfloat) from[12];
1169 to[4] = (GLfloat) from[1];
1170 to[5] = (GLfloat) from[5];
1171 to[6] = (GLfloat) from[9];
1172 to[7] = (GLfloat) from[13];
1173 to[8] = (GLfloat) from[2];
1174 to[9] = (GLfloat) from[6];
1175 to[10] = (GLfloat) from[10];
1176 to[11] = (GLfloat) from[14];
1177 to[12] = (GLfloat) from[3];
1178 to[13] = (GLfloat) from[7];
1179 to[14] = (GLfloat) from[11];
1180 to[15] = (GLfloat) from[15];
1181 }