1 /*
2  * Copyright (C) 2007 The Android Open Source Project
3  *
4  * Licensed under the Apache License, Version 2.0 (the "License");
5  * you may not use this file except in compliance with the License.
6  * You may obtain a copy of the License at
7  *
8  *      http://www.apache.org/licenses/LICENSE-2.0
9  *
10  * Unless required by applicable law or agreed to in writing, software
11  * distributed under the License is distributed on an "AS IS" BASIS,
12  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  * See the License for the specific language governing permissions and
14  * limitations under the License.
15  */
16 
17 package android.opengl;
18 
19 import androidx.annotation.NonNull;
20 
21 /**
22  * Matrix math utilities. These methods operate on OpenGL ES format
23  * matrices and vectors stored in float arrays.
24  * <p>
25  * Matrices are 4 x 4 column-vector matrices stored in column-major
26  * order:
27  * <pre>
28  *  m[offset +  0] m[offset +  4] m[offset +  8] m[offset + 12]
29  *  m[offset +  1] m[offset +  5] m[offset +  9] m[offset + 13]
30  *  m[offset +  2] m[offset +  6] m[offset + 10] m[offset + 14]
31  *  m[offset +  3] m[offset +  7] m[offset + 11] m[offset + 15]</pre>
32  *
33  * Vectors are 4 x 1 column vectors stored in order:
34  * <pre>
35  * v[offset + 0]
36  * v[offset + 1]
37  * v[offset + 2]
38  * v[offset + 3]</pre>
39  */
40 public class Matrix {
41 
42     /** Temporary memory for operations that need temporary matrix data. */
43     private static final ThreadLocal<float[]> ThreadTmp = new ThreadLocal() {
44         @Override protected float[] initialValue() {
45             return new float[32];
46         }
47     };
48 
49     /**
50      * @deprecated All methods are static, do not instantiate this class.
51      */
52     @Deprecated
Matrix()53     public Matrix() {}
54 
overlap( float[] a, int aStart, int aLength, float[] b, int bStart, int bLength)55     private static boolean overlap(
56             float[] a, int aStart, int aLength, float[] b, int bStart, int bLength) {
57         if (a != b) {
58             return false;
59         }
60 
61         if (aStart == bStart) {
62             return true;
63         }
64 
65         int aEnd = aStart + aLength;
66         int bEnd = bStart + bLength;
67 
68         if (aEnd == bEnd) {
69             return true;
70         }
71 
72         if (aStart < bStart && bStart < aEnd) {
73             return true;
74         }
75         if (aStart < bEnd   && bEnd   < aEnd) {
76             return true;
77         }
78 
79         if (bStart < aStart && aStart < bEnd) {
80             return true;
81         }
82         if (bStart < aEnd   && aEnd   < bEnd) {
83             return true;
84         }
85 
86         return false;
87     }
88 
89     /**
90      * Multiplies two 4x4 matrices together and stores the result in a third 4x4
91      * matrix. In matrix notation: result = lhs x rhs. Due to the way
92      * matrix multiplication works, the result matrix will have the same
93      * effect as first multiplying by the rhs matrix, then multiplying by
94      * the lhs matrix. This is the opposite of what you might expect.
95      * <p>
96      * The same float array may be passed for result, lhs, and/or rhs. This
97      * operation is expected to do the correct thing if the result elements
98      * overlap with either of the lhs or rhs elements.
99      *
100      * @param result The float array that holds the result.
101      * @param resultOffset The offset into the result array where the result is
102      *        stored.
103      * @param lhs The float array that holds the left-hand-side matrix.
104      * @param lhsOffset The offset into the lhs array where the lhs is stored
105      * @param rhs The float array that holds the right-hand-side matrix.
106      * @param rhsOffset The offset into the rhs array where the rhs is stored.
107      *
108      * @throws IllegalArgumentException under any of the following conditions:
109      * result, lhs, or rhs are null;
110      * resultOffset + 16 > result.length
111      * or lhsOffset + 16 > lhs.length
112      * or rhsOffset + 16 > rhs.length;
113      * resultOffset < 0 or lhsOffset < 0 or rhsOffset < 0
114      */
multiplyMM(float[] result, int resultOffset, float[] lhs, int lhsOffset, float[] rhs, int rhsOffset)115     public static void multiplyMM(float[] result, int resultOffset,
116             float[] lhs, int lhsOffset, float[] rhs, int rhsOffset) {
117         // error checking
118         if (result == null) {
119             throw new IllegalArgumentException("result == null");
120         }
121         if (lhs == null) {
122             throw new IllegalArgumentException("lhs == null");
123         }
124         if (rhs == null) {
125             throw new IllegalArgumentException("rhs == null");
126         }
127         if (resultOffset < 0) {
128             throw new IllegalArgumentException("resultOffset < 0");
129         }
130         if (lhsOffset < 0) {
131             throw new IllegalArgumentException("lhsOffset < 0");
132         }
133         if (rhsOffset < 0) {
134             throw new IllegalArgumentException("rhsOffset < 0");
135         }
136         if (result.length < resultOffset + 16) {
137             throw new IllegalArgumentException("result.length < resultOffset + 16");
138         }
139         if (lhs.length < lhsOffset + 16) {
140             throw new IllegalArgumentException("lhs.length < lhsOffset + 16");
141         }
142         if (rhs.length < rhsOffset + 16) {
143             throw new IllegalArgumentException("rhs.length < rhsOffset + 16");
144         }
145 
146         // Check for overlap between rhs and result or lhs and result
147         if ( overlap(result, resultOffset, 16, lhs, lhsOffset, 16)
148                 || overlap(result, resultOffset, 16, rhs, rhsOffset, 16) ) {
149             float[] tmp = ThreadTmp.get();
150             for (int i=0; i<4; i++) {
151                 final float rhs_i0 = rhs[ 4*i + 0 + rhsOffset ];
152                 float ri0 = lhs[ 0 + lhsOffset ] * rhs_i0;
153                 float ri1 = lhs[ 1 + lhsOffset ] * rhs_i0;
154                 float ri2 = lhs[ 2 + lhsOffset ] * rhs_i0;
155                 float ri3 = lhs[ 3 + lhsOffset ] * rhs_i0;
156                 for (int j=1; j<4; j++) {
157                     final float rhs_ij = rhs[ 4*i + j + rhsOffset];
158                     ri0 += lhs[ 4*j + 0 + lhsOffset ] * rhs_ij;
159                     ri1 += lhs[ 4*j + 1 + lhsOffset ] * rhs_ij;
160                     ri2 += lhs[ 4*j + 2 + lhsOffset ] * rhs_ij;
161                     ri3 += lhs[ 4*j + 3 + lhsOffset ] * rhs_ij;
162                 }
163                 tmp[ 4*i + 0 ] = ri0;
164                 tmp[ 4*i + 1 ] = ri1;
165                 tmp[ 4*i + 2 ] = ri2;
166                 tmp[ 4*i + 3 ] = ri3;
167             }
168 
169             // copy from tmp to result
170             for (int i=0; i < 16; i++) {
171                 result[ i + resultOffset ] = tmp[ i ];
172             }
173 
174         } else {
175             for (int i=0; i<4; i++) {
176                 final float rhs_i0 = rhs[ 4*i + 0 + rhsOffset ];
177                 float ri0 = lhs[ 0 + lhsOffset ] * rhs_i0;
178                 float ri1 = lhs[ 1 + lhsOffset ] * rhs_i0;
179                 float ri2 = lhs[ 2 + lhsOffset ] * rhs_i0;
180                 float ri3 = lhs[ 3 + lhsOffset ] * rhs_i0;
181                 for (int j=1; j<4; j++) {
182                     final float rhs_ij = rhs[ 4*i + j + rhsOffset];
183                     ri0 += lhs[ 4*j + 0 + lhsOffset ] * rhs_ij;
184                     ri1 += lhs[ 4*j + 1 + lhsOffset ] * rhs_ij;
185                     ri2 += lhs[ 4*j + 2 + lhsOffset ] * rhs_ij;
186                     ri3 += lhs[ 4*j + 3 + lhsOffset ] * rhs_ij;
187                 }
188                 result[ 4*i + 0 + resultOffset ] = ri0;
189                 result[ 4*i + 1 + resultOffset ] = ri1;
190                 result[ 4*i + 2 + resultOffset ] = ri2;
191                 result[ 4*i + 3 + resultOffset ] = ri3;
192             }
193         }
194     }
195 
196     /**
197      * Multiplies a 4 element vector by a 4x4 matrix and stores the result in a
198      * 4-element column vector. In matrix notation: result = lhs x rhs
199      * <p>
200      * The same float array may be passed for resultVec, lhsMat, and/or rhsVec.
201      * This operation is expected to do the correct thing if the result elements
202      * overlap with either of the lhs or rhs elements.
203      *
204      * @param resultVec The float array that holds the result vector.
205      * @param resultVecOffset The offset into the result array where the result
206      *        vector is stored.
207      * @param lhsMat The float array that holds the left-hand-side matrix.
208      * @param lhsMatOffset The offset into the lhs array where the lhs is stored
209      * @param rhsVec The float array that holds the right-hand-side vector.
210      * @param rhsVecOffset The offset into the rhs vector where the rhs vector
211      *        is stored.
212      *
213      * @throws IllegalArgumentException under any of the following conditions:
214      * resultVec, lhsMat, or rhsVec are null;
215      * resultVecOffset + 4  > resultVec.length
216      * or lhsMatOffset + 16 > lhsMat.length
217      * or rhsVecOffset + 4  > rhsVec.length;
218      * resultVecOffset < 0 or lhsMatOffset < 0 or rhsVecOffset < 0
219      */
multiplyMV(float[] resultVec, int resultVecOffset, float[] lhsMat, int lhsMatOffset, float[] rhsVec, int rhsVecOffset)220     public static void multiplyMV(float[] resultVec,
221             int resultVecOffset, float[] lhsMat, int lhsMatOffset,
222             float[] rhsVec, int rhsVecOffset) {
223         // error checking
224         if (resultVec == null) {
225             throw new IllegalArgumentException("resultVec == null");
226         }
227         if (lhsMat == null) {
228             throw new IllegalArgumentException("lhsMat == null");
229         }
230         if (rhsVec == null) {
231             throw new IllegalArgumentException("rhsVec == null");
232         }
233         if (resultVecOffset < 0) {
234             throw new IllegalArgumentException("resultVecOffset < 0");
235         }
236         if (lhsMatOffset < 0) {
237             throw new IllegalArgumentException("lhsMatOffset < 0");
238         }
239         if (rhsVecOffset < 0) {
240             throw new IllegalArgumentException("rhsVecOffset < 0");
241         }
242         if (resultVec.length < resultVecOffset + 4) {
243             throw new IllegalArgumentException("resultVec.length < resultVecOffset + 4");
244         }
245         if (lhsMat.length < lhsMatOffset + 16) {
246             throw new IllegalArgumentException("lhsMat.length < lhsMatOffset + 16");
247         }
248         if (rhsVec.length < rhsVecOffset + 4) {
249             throw new IllegalArgumentException("rhsVec.length < rhsVecOffset + 4");
250         }
251 
252         float tmp0 = lhsMat[0 + 4 * 0 + lhsMatOffset] * rhsVec[0 + rhsVecOffset] +
253                      lhsMat[0 + 4 * 1 + lhsMatOffset] * rhsVec[1 + rhsVecOffset] +
254                      lhsMat[0 + 4 * 2 + lhsMatOffset] * rhsVec[2 + rhsVecOffset] +
255                      lhsMat[0 + 4 * 3 + lhsMatOffset] * rhsVec[3 + rhsVecOffset];
256         float tmp1 = lhsMat[1 + 4 * 0 + lhsMatOffset] * rhsVec[0 + rhsVecOffset] +
257                      lhsMat[1 + 4 * 1 + lhsMatOffset] * rhsVec[1 + rhsVecOffset] +
258                      lhsMat[1 + 4 * 2 + lhsMatOffset] * rhsVec[2 + rhsVecOffset] +
259                      lhsMat[1 + 4 * 3 + lhsMatOffset] * rhsVec[3 + rhsVecOffset];
260         float tmp2 = lhsMat[2 + 4 * 0 + lhsMatOffset] * rhsVec[0 + rhsVecOffset] +
261                      lhsMat[2 + 4 * 1 + lhsMatOffset] * rhsVec[1 + rhsVecOffset] +
262                      lhsMat[2 + 4 * 2 + lhsMatOffset] * rhsVec[2 + rhsVecOffset] +
263                      lhsMat[2 + 4 * 3 + lhsMatOffset] * rhsVec[3 + rhsVecOffset];
264         float tmp3 = lhsMat[3 + 4 * 0 + lhsMatOffset] * rhsVec[0 + rhsVecOffset] +
265                      lhsMat[3 + 4 * 1 + lhsMatOffset] * rhsVec[1 + rhsVecOffset] +
266                      lhsMat[3 + 4 * 2 + lhsMatOffset] * rhsVec[2 + rhsVecOffset] +
267                      lhsMat[3 + 4 * 3 + lhsMatOffset] * rhsVec[3 + rhsVecOffset];
268 
269         resultVec[ 0 + resultVecOffset ] = tmp0;
270         resultVec[ 1 + resultVecOffset ] = tmp1;
271         resultVec[ 2 + resultVecOffset ] = tmp2;
272         resultVec[ 3 + resultVecOffset ] = tmp3;
273     }
274 
275     /**
276      * Transposes a 4 x 4 matrix.
277      * <p>
278      * mTrans and m must not overlap.
279      *
280      * @param mTrans the array that holds the output transposed matrix
281      * @param mTransOffset an offset into mTrans where the transposed matrix is
282      *        stored.
283      * @param m the input array
284      * @param mOffset an offset into m where the input matrix is stored.
285      */
transposeM(float[] mTrans, int mTransOffset, float[] m, int mOffset)286     public static void transposeM(float[] mTrans, int mTransOffset, float[] m,
287             int mOffset) {
288         for (int i = 0; i < 4; i++) {
289             int mBase = i * 4 + mOffset;
290             mTrans[i + mTransOffset] = m[mBase];
291             mTrans[i + 4 + mTransOffset] = m[mBase + 1];
292             mTrans[i + 8 + mTransOffset] = m[mBase + 2];
293             mTrans[i + 12 + mTransOffset] = m[mBase + 3];
294         }
295     }
296 
297     /**
298      * Inverts a 4 x 4 matrix.
299      * <p>
300      * mInv and m must not overlap.
301      *
302      * @param mInv the array that holds the output inverted matrix
303      * @param mInvOffset an offset into mInv where the inverted matrix is
304      *        stored.
305      * @param m the input array
306      * @param mOffset an offset into m where the input matrix is stored.
307      * @return true if the matrix could be inverted, false if it could not.
308      */
invertM(float[] mInv, int mInvOffset, float[] m, int mOffset)309     public static boolean invertM(float[] mInv, int mInvOffset, float[] m,
310             int mOffset) {
311         // Invert a 4 x 4 matrix using Cramer's Rule
312 
313         // transpose matrix
314         final float src0  = m[mOffset +  0];
315         final float src4  = m[mOffset +  1];
316         final float src8  = m[mOffset +  2];
317         final float src12 = m[mOffset +  3];
318 
319         final float src1  = m[mOffset +  4];
320         final float src5  = m[mOffset +  5];
321         final float src9  = m[mOffset +  6];
322         final float src13 = m[mOffset +  7];
323 
324         final float src2  = m[mOffset +  8];
325         final float src6  = m[mOffset +  9];
326         final float src10 = m[mOffset + 10];
327         final float src14 = m[mOffset + 11];
328 
329         final float src3  = m[mOffset + 12];
330         final float src7  = m[mOffset + 13];
331         final float src11 = m[mOffset + 14];
332         final float src15 = m[mOffset + 15];
333 
334         // calculate pairs for first 8 elements (cofactors)
335         final float atmp0  = src10 * src15;
336         final float atmp1  = src11 * src14;
337         final float atmp2  = src9  * src15;
338         final float atmp3  = src11 * src13;
339         final float atmp4  = src9  * src14;
340         final float atmp5  = src10 * src13;
341         final float atmp6  = src8  * src15;
342         final float atmp7  = src11 * src12;
343         final float atmp8  = src8  * src14;
344         final float atmp9  = src10 * src12;
345         final float atmp10 = src8  * src13;
346         final float atmp11 = src9  * src12;
347 
348         // calculate first 8 elements (cofactors)
349         final float dst0  = (atmp0 * src5 + atmp3 * src6 + atmp4  * src7)
350                           - (atmp1 * src5 + atmp2 * src6 + atmp5  * src7);
351         final float dst1  = (atmp1 * src4 + atmp6 * src6 + atmp9  * src7)
352                           - (atmp0 * src4 + atmp7 * src6 + atmp8  * src7);
353         final float dst2  = (atmp2 * src4 + atmp7 * src5 + atmp10 * src7)
354                           - (atmp3 * src4 + atmp6 * src5 + atmp11 * src7);
355         final float dst3  = (atmp5 * src4 + atmp8 * src5 + atmp11 * src6)
356                           - (atmp4 * src4 + atmp9 * src5 + atmp10 * src6);
357         final float dst4  = (atmp1 * src1 + atmp2 * src2 + atmp5  * src3)
358                           - (atmp0 * src1 + atmp3 * src2 + atmp4  * src3);
359         final float dst5  = (atmp0 * src0 + atmp7 * src2 + atmp8  * src3)
360                           - (atmp1 * src0 + atmp6 * src2 + atmp9  * src3);
361         final float dst6  = (atmp3 * src0 + atmp6 * src1 + atmp11 * src3)
362                           - (atmp2 * src0 + atmp7 * src1 + atmp10 * src3);
363         final float dst7  = (atmp4 * src0 + atmp9 * src1 + atmp10 * src2)
364                           - (atmp5 * src0 + atmp8 * src1 + atmp11 * src2);
365 
366         // calculate pairs for second 8 elements (cofactors)
367         final float btmp0  = src2 * src7;
368         final float btmp1  = src3 * src6;
369         final float btmp2  = src1 * src7;
370         final float btmp3  = src3 * src5;
371         final float btmp4  = src1 * src6;
372         final float btmp5  = src2 * src5;
373         final float btmp6  = src0 * src7;
374         final float btmp7  = src3 * src4;
375         final float btmp8  = src0 * src6;
376         final float btmp9  = src2 * src4;
377         final float btmp10 = src0 * src5;
378         final float btmp11 = src1 * src4;
379 
380         // calculate second 8 elements (cofactors)
381         final float dst8  = (btmp0  * src13 + btmp3  * src14 + btmp4  * src15)
382                           - (btmp1  * src13 + btmp2  * src14 + btmp5  * src15);
383         final float dst9  = (btmp1  * src12 + btmp6  * src14 + btmp9  * src15)
384                           - (btmp0  * src12 + btmp7  * src14 + btmp8  * src15);
385         final float dst10 = (btmp2  * src12 + btmp7  * src13 + btmp10 * src15)
386                           - (btmp3  * src12 + btmp6  * src13 + btmp11 * src15);
387         final float dst11 = (btmp5  * src12 + btmp8  * src13 + btmp11 * src14)
388                           - (btmp4  * src12 + btmp9  * src13 + btmp10 * src14);
389         final float dst12 = (btmp2  * src10 + btmp5  * src11 + btmp1  * src9 )
390                           - (btmp4  * src11 + btmp0  * src9  + btmp3  * src10);
391         final float dst13 = (btmp8  * src11 + btmp0  * src8  + btmp7  * src10)
392                           - (btmp6  * src10 + btmp9  * src11 + btmp1  * src8 );
393         final float dst14 = (btmp6  * src9  + btmp11 * src11 + btmp3  * src8 )
394                           - (btmp10 * src11 + btmp2  * src8  + btmp7  * src9 );
395         final float dst15 = (btmp10 * src10 + btmp4  * src8  + btmp9  * src9 )
396                           - (btmp8  * src9  + btmp11 * src10 + btmp5  * src8 );
397 
398         // calculate determinant
399         final float det =
400                 src0 * dst0 + src1 * dst1 + src2 * dst2 + src3 * dst3;
401 
402         if (det == 0.0f) {
403             return false;
404         }
405 
406         // calculate matrix inverse
407         final float invdet = 1.0f / det;
408         mInv[     mInvOffset] = dst0  * invdet;
409         mInv[ 1 + mInvOffset] = dst1  * invdet;
410         mInv[ 2 + mInvOffset] = dst2  * invdet;
411         mInv[ 3 + mInvOffset] = dst3  * invdet;
412 
413         mInv[ 4 + mInvOffset] = dst4  * invdet;
414         mInv[ 5 + mInvOffset] = dst5  * invdet;
415         mInv[ 6 + mInvOffset] = dst6  * invdet;
416         mInv[ 7 + mInvOffset] = dst7  * invdet;
417 
418         mInv[ 8 + mInvOffset] = dst8  * invdet;
419         mInv[ 9 + mInvOffset] = dst9  * invdet;
420         mInv[10 + mInvOffset] = dst10 * invdet;
421         mInv[11 + mInvOffset] = dst11 * invdet;
422 
423         mInv[12 + mInvOffset] = dst12 * invdet;
424         mInv[13 + mInvOffset] = dst13 * invdet;
425         mInv[14 + mInvOffset] = dst14 * invdet;
426         mInv[15 + mInvOffset] = dst15 * invdet;
427 
428         return true;
429     }
430 
431     /**
432      * Computes an orthographic projection matrix.
433      *
434      * @param m returns the result
435      * @param mOffset
436      * @param left
437      * @param right
438      * @param bottom
439      * @param top
440      * @param near
441      * @param far
442      */
orthoM(float[] m, int mOffset, float left, float right, float bottom, float top, float near, float far)443     public static void orthoM(float[] m, int mOffset,
444         float left, float right, float bottom, float top,
445         float near, float far) {
446         if (left == right) {
447             throw new IllegalArgumentException("left == right");
448         }
449         if (bottom == top) {
450             throw new IllegalArgumentException("bottom == top");
451         }
452         if (near == far) {
453             throw new IllegalArgumentException("near == far");
454         }
455 
456         final float r_width  = 1.0f / (right - left);
457         final float r_height = 1.0f / (top - bottom);
458         final float r_depth  = 1.0f / (far - near);
459         final float x =  2.0f * (r_width);
460         final float y =  2.0f * (r_height);
461         final float z = -2.0f * (r_depth);
462         final float tx = -(right + left) * r_width;
463         final float ty = -(top + bottom) * r_height;
464         final float tz = -(far + near) * r_depth;
465         m[mOffset + 0] = x;
466         m[mOffset + 5] = y;
467         m[mOffset +10] = z;
468         m[mOffset +12] = tx;
469         m[mOffset +13] = ty;
470         m[mOffset +14] = tz;
471         m[mOffset +15] = 1.0f;
472         m[mOffset + 1] = 0.0f;
473         m[mOffset + 2] = 0.0f;
474         m[mOffset + 3] = 0.0f;
475         m[mOffset + 4] = 0.0f;
476         m[mOffset + 6] = 0.0f;
477         m[mOffset + 7] = 0.0f;
478         m[mOffset + 8] = 0.0f;
479         m[mOffset + 9] = 0.0f;
480         m[mOffset + 11] = 0.0f;
481     }
482 
483 
484     /**
485      * Defines a projection matrix in terms of six clip planes.
486      *
487      * @param m the float array that holds the output perspective matrix
488      * @param offset the offset into float array m where the perspective
489      *        matrix data is written
490      * @param left
491      * @param right
492      * @param bottom
493      * @param top
494      * @param near
495      * @param far
496      */
frustumM(float[] m, int offset, float left, float right, float bottom, float top, float near, float far)497     public static void frustumM(float[] m, int offset,
498             float left, float right, float bottom, float top,
499             float near, float far) {
500         if (left == right) {
501             throw new IllegalArgumentException("left == right");
502         }
503         if (top == bottom) {
504             throw new IllegalArgumentException("top == bottom");
505         }
506         if (near == far) {
507             throw new IllegalArgumentException("near == far");
508         }
509         if (near <= 0.0f) {
510             throw new IllegalArgumentException("near <= 0.0f");
511         }
512         if (far <= 0.0f) {
513             throw new IllegalArgumentException("far <= 0.0f");
514         }
515         final float r_width  = 1.0f / (right - left);
516         final float r_height = 1.0f / (top - bottom);
517         final float r_depth  = 1.0f / (near - far);
518         final float x = 2.0f * (near * r_width);
519         final float y = 2.0f * (near * r_height);
520         final float A = (right + left) * r_width;
521         final float B = (top + bottom) * r_height;
522         final float C = (far + near) * r_depth;
523         final float D = 2.0f * (far * near * r_depth);
524         m[offset + 0] = x;
525         m[offset + 5] = y;
526         m[offset + 8] = A;
527         m[offset +  9] = B;
528         m[offset + 10] = C;
529         m[offset + 14] = D;
530         m[offset + 11] = -1.0f;
531         m[offset +  1] = 0.0f;
532         m[offset +  2] = 0.0f;
533         m[offset +  3] = 0.0f;
534         m[offset +  4] = 0.0f;
535         m[offset +  6] = 0.0f;
536         m[offset +  7] = 0.0f;
537         m[offset + 12] = 0.0f;
538         m[offset + 13] = 0.0f;
539         m[offset + 15] = 0.0f;
540     }
541 
542     /**
543      * Defines a projection matrix in terms of a field of view angle, an
544      * aspect ratio, and z clip planes.
545      *
546      * @param m the float array that holds the perspective matrix
547      * @param offset the offset into float array m where the perspective
548      *        matrix data is written
549      * @param fovy field of view in y direction, in degrees
550      * @param aspect width to height aspect ratio of the viewport
551      * @param zNear
552      * @param zFar
553      */
perspectiveM(float[] m, int offset, float fovy, float aspect, float zNear, float zFar)554     public static void perspectiveM(float[] m, int offset,
555           float fovy, float aspect, float zNear, float zFar) {
556         float f = 1.0f / (float) Math.tan(fovy * (Math.PI / 360.0));
557         float rangeReciprocal = 1.0f / (zNear - zFar);
558 
559         m[offset + 0] = f / aspect;
560         m[offset + 1] = 0.0f;
561         m[offset + 2] = 0.0f;
562         m[offset + 3] = 0.0f;
563 
564         m[offset + 4] = 0.0f;
565         m[offset + 5] = f;
566         m[offset + 6] = 0.0f;
567         m[offset + 7] = 0.0f;
568 
569         m[offset + 8] = 0.0f;
570         m[offset + 9] = 0.0f;
571         m[offset + 10] = (zFar + zNear) * rangeReciprocal;
572         m[offset + 11] = -1.0f;
573 
574         m[offset + 12] = 0.0f;
575         m[offset + 13] = 0.0f;
576         m[offset + 14] = 2.0f * zFar * zNear * rangeReciprocal;
577         m[offset + 15] = 0.0f;
578     }
579 
580     /**
581      * Computes the length of a vector.
582      *
583      * @param x x coordinate of a vector
584      * @param y y coordinate of a vector
585      * @param z z coordinate of a vector
586      * @return the length of a vector
587      */
length(float x, float y, float z)588     public static float length(float x, float y, float z) {
589         return (float) Math.sqrt(x * x + y * y + z * z);
590     }
591 
592     /**
593      * Sets matrix m to the identity matrix.
594      *
595      * @param sm returns the result
596      * @param smOffset index into sm where the result matrix starts
597      */
setIdentityM(float[] sm, int smOffset)598     public static void setIdentityM(float[] sm, int smOffset) {
599         for (int i=0 ; i<16 ; i++) {
600             sm[smOffset + i] = 0;
601         }
602         for(int i = 0; i < 16; i += 5) {
603             sm[smOffset + i] = 1.0f;
604         }
605     }
606 
607     /**
608      * Scales matrix m by x, y, and z, putting the result in sm.
609      * <p>
610      * m and sm must not overlap.
611      *
612      * @param sm returns the result
613      * @param smOffset index into sm where the result matrix starts
614      * @param m source matrix
615      * @param mOffset index into m where the source matrix starts
616      * @param x scale factor x
617      * @param y scale factor y
618      * @param z scale factor z
619      */
scaleM(float[] sm, int smOffset, float[] m, int mOffset, float x, float y, float z)620     public static void scaleM(float[] sm, int smOffset,
621             float[] m, int mOffset,
622             float x, float y, float z) {
623         for (int i=0 ; i<4 ; i++) {
624             int smi = smOffset + i;
625             int mi = mOffset + i;
626             sm[     smi] = m[     mi] * x;
627             sm[ 4 + smi] = m[ 4 + mi] * y;
628             sm[ 8 + smi] = m[ 8 + mi] * z;
629             sm[12 + smi] = m[12 + mi];
630         }
631     }
632 
633     /**
634      * Scales matrix m in place by sx, sy, and sz.
635      *
636      * @param m matrix to scale
637      * @param mOffset index into m where the matrix starts
638      * @param x scale factor x
639      * @param y scale factor y
640      * @param z scale factor z
641      */
scaleM(float[] m, int mOffset, float x, float y, float z)642     public static void scaleM(float[] m, int mOffset,
643             float x, float y, float z) {
644         for (int i=0 ; i<4 ; i++) {
645             int mi = mOffset + i;
646             m[     mi] *= x;
647             m[ 4 + mi] *= y;
648             m[ 8 + mi] *= z;
649         }
650     }
651 
652     /**
653      * Translates matrix m by x, y, and z, putting the result in tm.
654      * <p>
655      * m and tm must not overlap.
656      *
657      * @param tm returns the result
658      * @param tmOffset index into sm where the result matrix starts
659      * @param m source matrix
660      * @param mOffset index into m where the source matrix starts
661      * @param x translation factor x
662      * @param y translation factor y
663      * @param z translation factor z
664      */
translateM(float[] tm, int tmOffset, float[] m, int mOffset, float x, float y, float z)665     public static void translateM(float[] tm, int tmOffset,
666             float[] m, int mOffset,
667             float x, float y, float z) {
668         for (int i=0 ; i<12 ; i++) {
669             tm[tmOffset + i] = m[mOffset + i];
670         }
671         for (int i=0 ; i<4 ; i++) {
672             int tmi = tmOffset + i;
673             int mi = mOffset + i;
674             tm[12 + tmi] = m[mi] * x + m[4 + mi] * y + m[8 + mi] * z +
675                 m[12 + mi];
676         }
677     }
678 
679     /**
680      * Translates matrix m by x, y, and z in place.
681      *
682      * @param m matrix
683      * @param mOffset index into m where the matrix starts
684      * @param x translation factor x
685      * @param y translation factor y
686      * @param z translation factor z
687      */
translateM( float[] m, int mOffset, float x, float y, float z)688     public static void translateM(
689             float[] m, int mOffset,
690             float x, float y, float z) {
691         for (int i=0 ; i<4 ; i++) {
692             int mi = mOffset + i;
693             m[12 + mi] += m[mi] * x + m[4 + mi] * y + m[8 + mi] * z;
694         }
695     }
696 
697     /**
698      * Rotates matrix m by angle a (in degrees) around the axis (x, y, z).
699      * <p>
700      * m and rm must not overlap.
701      *
702      * @param rm returns the result
703      * @param rmOffset index into rm where the result matrix starts
704      * @param m source matrix
705      * @param mOffset index into m where the source matrix starts
706      * @param a angle to rotate in degrees
707      * @param x X axis component
708      * @param y Y axis component
709      * @param z Z axis component
710      */
rotateM(float[] rm, int rmOffset, float[] m, int mOffset, float a, float x, float y, float z)711     public static void rotateM(float[] rm, int rmOffset,
712             float[] m, int mOffset,
713             float a, float x, float y, float z) {
714         float[] tmp = ThreadTmp.get();
715         setRotateM(tmp, 16, a, x, y, z);
716         multiplyMM(rm, rmOffset, m, mOffset, tmp, 16);
717     }
718 
719     /**
720      * Rotates matrix m in place by angle a (in degrees)
721      * around the axis (x, y, z).
722      *
723      * @param m source matrix
724      * @param mOffset index into m where the matrix starts
725      * @param a angle to rotate in degrees
726      * @param x X axis component
727      * @param y Y axis component
728      * @param z Z axis component
729      */
rotateM(float[] m, int mOffset, float a, float x, float y, float z)730     public static void rotateM(float[] m, int mOffset,
731             float a, float x, float y, float z) {
732         rotateM(m, mOffset, m, mOffset, a, x, y, z);
733     }
734 
735     /**
736      * Creates a matrix for rotation by angle a (in degrees)
737      * around the axis (x, y, z).
738      * <p>
739      * An optimized path will be used for rotation about a major axis
740      * (e.g. x=1.0f y=0.0f z=0.0f).
741      *
742      * @param rm returns the result
743      * @param rmOffset index into rm where the result matrix starts
744      * @param a angle to rotate in degrees
745      * @param x X axis component
746      * @param y Y axis component
747      * @param z Z axis component
748      */
setRotateM(float[] rm, int rmOffset, float a, float x, float y, float z)749     public static void setRotateM(float[] rm, int rmOffset,
750             float a, float x, float y, float z) {
751         rm[rmOffset + 3] = 0;
752         rm[rmOffset + 7] = 0;
753         rm[rmOffset + 11]= 0;
754         rm[rmOffset + 12]= 0;
755         rm[rmOffset + 13]= 0;
756         rm[rmOffset + 14]= 0;
757         rm[rmOffset + 15]= 1;
758         a *= (float) (Math.PI / 180.0f);
759         float s = (float) Math.sin(a);
760         float c = (float) Math.cos(a);
761         if (1.0f == x && 0.0f == y && 0.0f == z) {
762             rm[rmOffset + 5] = c;   rm[rmOffset + 10]= c;
763             rm[rmOffset + 6] = s;   rm[rmOffset + 9] = -s;
764             rm[rmOffset + 1] = 0;   rm[rmOffset + 2] = 0;
765             rm[rmOffset + 4] = 0;   rm[rmOffset + 8] = 0;
766             rm[rmOffset + 0] = 1;
767         } else if (0.0f == x && 1.0f == y && 0.0f == z) {
768             rm[rmOffset + 0] = c;   rm[rmOffset + 10]= c;
769             rm[rmOffset + 8] = s;   rm[rmOffset + 2] = -s;
770             rm[rmOffset + 1] = 0;   rm[rmOffset + 4] = 0;
771             rm[rmOffset + 6] = 0;   rm[rmOffset + 9] = 0;
772             rm[rmOffset + 5] = 1;
773         } else if (0.0f == x && 0.0f == y && 1.0f == z) {
774             rm[rmOffset + 0] = c;   rm[rmOffset + 5] = c;
775             rm[rmOffset + 1] = s;   rm[rmOffset + 4] = -s;
776             rm[rmOffset + 2] = 0;   rm[rmOffset + 6] = 0;
777             rm[rmOffset + 8] = 0;   rm[rmOffset + 9] = 0;
778             rm[rmOffset + 10]= 1;
779         } else {
780             float len = length(x, y, z);
781             if (1.0f != len) {
782                 float recipLen = 1.0f / len;
783                 x *= recipLen;
784                 y *= recipLen;
785                 z *= recipLen;
786             }
787             float nc = 1.0f - c;
788             float xy = x * y;
789             float yz = y * z;
790             float zx = z * x;
791             float xs = x * s;
792             float ys = y * s;
793             float zs = z * s;
794             rm[rmOffset +  0] = x*x*nc +  c;
795             rm[rmOffset +  4] =  xy*nc - zs;
796             rm[rmOffset +  8] =  zx*nc + ys;
797             rm[rmOffset +  1] =  xy*nc + zs;
798             rm[rmOffset +  5] = y*y*nc +  c;
799             rm[rmOffset +  9] =  yz*nc - xs;
800             rm[rmOffset +  2] =  zx*nc - ys;
801             rm[rmOffset +  6] =  yz*nc + xs;
802             rm[rmOffset + 10] = z*z*nc +  c;
803         }
804     }
805 
806     /**
807      * Converts Euler angles to a rotation matrix.
808      *
809      * @param rm returns the result
810      * @param rmOffset index into rm where the result matrix starts
811      * @param x angle of rotation, in degrees
812      * @param y is broken, do not use
813      * @param z angle of rotation, in degrees
814      *
815      * @deprecated This method is incorrect around the y axis. This method is
816      *             deprecated and replaced (below) by setRotateEulerM2 which
817      *             behaves correctly
818      */
819     @Deprecated
setRotateEulerM(float[] rm, int rmOffset, float x, float y, float z)820     public static void setRotateEulerM(float[] rm, int rmOffset,
821             float x, float y, float z) {
822         x *= (float) (Math.PI / 180.0f);
823         y *= (float) (Math.PI / 180.0f);
824         z *= (float) (Math.PI / 180.0f);
825         float cx = (float) Math.cos(x);
826         float sx = (float) Math.sin(x);
827         float cy = (float) Math.cos(y);
828         float sy = (float) Math.sin(y);
829         float cz = (float) Math.cos(z);
830         float sz = (float) Math.sin(z);
831         float cxsy = cx * sy;
832         float sxsy = sx * sy;
833 
834         rm[rmOffset + 0]  =   cy * cz;
835         rm[rmOffset + 1]  =  -cy * sz;
836         rm[rmOffset + 2]  =   sy;
837         rm[rmOffset + 3]  =  0.0f;
838 
839         rm[rmOffset + 4]  =  cxsy * cz + cx * sz;
840         rm[rmOffset + 5]  = -cxsy * sz + cx * cz;
841         rm[rmOffset + 6]  =  -sx * cy;
842         rm[rmOffset + 7]  =  0.0f;
843 
844         rm[rmOffset + 8]  = -sxsy * cz + sx * sz;
845         rm[rmOffset + 9]  =  sxsy * sz + sx * cz;
846         rm[rmOffset + 10] =  cx * cy;
847         rm[rmOffset + 11] =  0.0f;
848 
849         rm[rmOffset + 12] =  0.0f;
850         rm[rmOffset + 13] =  0.0f;
851         rm[rmOffset + 14] =  0.0f;
852         rm[rmOffset + 15] =  1.0f;
853     }
854 
855     /**
856      * Converts Euler angles to a rotation matrix.
857      *
858      * @param rm returns the result
859      * @param rmOffset index into rm where the result matrix starts
860      * @param x angle of rotation, in degrees
861      * @param y angle of rotation, in degrees
862      * @param z angle of rotation, in degrees
863      *
864      * @throws IllegalArgumentException if rm is null;
865      * or if rmOffset + 16 > rm.length;
866      * rmOffset < 0
867      */
setRotateEulerM2(@onNull float[] rm, int rmOffset, float x, float y, float z)868     public static void setRotateEulerM2(@NonNull float[] rm, int rmOffset,
869             float x, float y, float z) {
870         if (rm == null) {
871             throw new IllegalArgumentException("rm == null");
872         }
873         if (rmOffset < 0) {
874             throw new IllegalArgumentException("rmOffset < 0");
875         }
876         if (rm.length < rmOffset + 16) {
877             throw new IllegalArgumentException("rm.length < rmOffset + 16");
878         }
879 
880         x *= (float) (Math.PI / 180.0f);
881         y *= (float) (Math.PI / 180.0f);
882         z *= (float) (Math.PI / 180.0f);
883         float cx = (float) Math.cos(x);
884         float sx = (float) Math.sin(x);
885         float cy = (float) Math.cos(y);
886         float sy = (float) Math.sin(y);
887         float cz = (float) Math.cos(z);
888         float sz = (float) Math.sin(z);
889         float cxsy = cx * sy;
890         float sxsy = sx * sy;
891 
892         rm[rmOffset + 0]  =  cy * cz;
893         rm[rmOffset + 1]  = -cy * sz;
894         rm[rmOffset + 2]  =  sy;
895         rm[rmOffset + 3]  =  0.0f;
896 
897         rm[rmOffset + 4]  =  sxsy * cz + cx * sz;
898         rm[rmOffset + 5]  = -sxsy * sz + cx * cz;
899         rm[rmOffset + 6]  = -sx * cy;
900         rm[rmOffset + 7]  =  0.0f;
901 
902         rm[rmOffset + 8]  = -cxsy * cz + sx * sz;
903         rm[rmOffset + 9]  =  cxsy * sz + sx * cz;
904         rm[rmOffset + 10] =  cx * cy;
905         rm[rmOffset + 11] =  0.0f;
906 
907         rm[rmOffset + 12] =  0.0f;
908         rm[rmOffset + 13] =  0.0f;
909         rm[rmOffset + 14] =  0.0f;
910         rm[rmOffset + 15] =  1.0f;
911     }
912 
913     /**
914      * Defines a viewing transformation in terms of an eye point, a center of
915      * view, and an up vector.
916      *
917      * @param rm returns the result
918      * @param rmOffset index into rm where the result matrix starts
919      * @param eyeX eye point X
920      * @param eyeY eye point Y
921      * @param eyeZ eye point Z
922      * @param centerX center of view X
923      * @param centerY center of view Y
924      * @param centerZ center of view Z
925      * @param upX up vector X
926      * @param upY up vector Y
927      * @param upZ up vector Z
928      */
setLookAtM(float[] rm, int rmOffset, float eyeX, float eyeY, float eyeZ, float centerX, float centerY, float centerZ, float upX, float upY, float upZ)929     public static void setLookAtM(float[] rm, int rmOffset,
930             float eyeX, float eyeY, float eyeZ,
931             float centerX, float centerY, float centerZ, float upX, float upY,
932             float upZ) {
933 
934         // See the OpenGL GLUT documentation for gluLookAt for a description
935         // of the algorithm. We implement it in a straightforward way:
936 
937         float fx = centerX - eyeX;
938         float fy = centerY - eyeY;
939         float fz = centerZ - eyeZ;
940 
941         // Normalize f
942         float rlf = 1.0f / Matrix.length(fx, fy, fz);
943         fx *= rlf;
944         fy *= rlf;
945         fz *= rlf;
946 
947         // compute s = f x up (x means "cross product")
948         float sx = fy * upZ - fz * upY;
949         float sy = fz * upX - fx * upZ;
950         float sz = fx * upY - fy * upX;
951 
952         // and normalize s
953         float rls = 1.0f / Matrix.length(sx, sy, sz);
954         sx *= rls;
955         sy *= rls;
956         sz *= rls;
957 
958         // compute u = s x f
959         float ux = sy * fz - sz * fy;
960         float uy = sz * fx - sx * fz;
961         float uz = sx * fy - sy * fx;
962 
963         rm[rmOffset + 0] = sx;
964         rm[rmOffset + 1] = ux;
965         rm[rmOffset + 2] = -fx;
966         rm[rmOffset + 3] = 0.0f;
967 
968         rm[rmOffset + 4] = sy;
969         rm[rmOffset + 5] = uy;
970         rm[rmOffset + 6] = -fy;
971         rm[rmOffset + 7] = 0.0f;
972 
973         rm[rmOffset + 8] = sz;
974         rm[rmOffset + 9] = uz;
975         rm[rmOffset + 10] = -fz;
976         rm[rmOffset + 11] = 0.0f;
977 
978         rm[rmOffset + 12] = 0.0f;
979         rm[rmOffset + 13] = 0.0f;
980         rm[rmOffset + 14] = 0.0f;
981         rm[rmOffset + 15] = 1.0f;
982 
983         translateM(rm, rmOffset, -eyeX, -eyeY, -eyeZ);
984     }
985 }
986