home *** CD-ROM | disk | FTP | other *** search
/ Magazyn Amiga 5 / MA_Cover_5.iso / ppc / mesa / src / matrix.c < prev    next >
Encoding:
C/C++ Source or Header  |  1998-01-31  |  28.1 KB  |  1,025 lines

  1. /* $Id: matrix.c,v 1.22 1997/10/16 23:37:23 brianp Exp $ */
  2.  
  3. /*
  4.  * Mesa 3-D graphics library
  5.  * Version:  2.5
  6.  * Copyright (C) 1995-1997  Brian Paul
  7.  *
  8.  * This library is free software; you can redistribute it and/or
  9.  * modify it under the terms of the GNU Library General Public
  10.  * License as published by the Free Software Foundation; either
  11.  * version 2 of the License, or (at your option) any later version.
  12.  *
  13.  * This library is distributed in the hope that it will be useful,
  14.  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  15.  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  16.  * Library General Public License for more details.
  17.  *
  18.  * You should have received a copy of the GNU Library General Public
  19.  * License along with this library; if not, write to the Free
  20.  * Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
  21.  */
  22.  
  23.  
  24. /*
  25.  * $Log: matrix.c,v $
  26.  * Revision 1.22  1997/10/16 23:37:23  brianp
  27.  * fixed scotter's email address
  28.  *
  29.  * Revision 1.21  1997/08/13 01:54:34  brianp
  30.  * new matrix invert code from Scott McCaskill
  31.  *
  32.  * Revision 1.20  1997/07/24 01:23:16  brianp
  33.  * changed precompiled header symbol from PCH to PC_HEADER
  34.  *
  35.  * Revision 1.19  1997/05/30 02:21:43  brianp
  36.  * gl_PopMatrix() set ctx->New*Matrix flag incorrectly
  37.  *
  38.  * Revision 1.18  1997/05/28 04:06:03  brianp
  39.  * implemented projection near/far value stack for Driver.NearFar() function
  40.  *
  41.  * Revision 1.17  1997/05/28 03:25:43  brianp
  42.  * added precompiled header (PCH) support
  43.  *
  44.  * Revision 1.16  1997/05/01 01:39:40  brianp
  45.  * replace sqrt() with GL_SQRT()
  46.  *
  47.  * Revision 1.15  1997/04/21 01:20:41  brianp
  48.  * added MATRIX_2D_NO_ROT
  49.  *
  50.  * Revision 1.14  1997/04/20 20:28:49  brianp
  51.  * replaced abort() with gl_problem()
  52.  *
  53.  * Revision 1.13  1997/04/20 16:31:08  brianp
  54.  * added NearFar device driver function
  55.  *
  56.  * Revision 1.12  1997/04/20 16:18:15  brianp
  57.  * added glOrtho and glFrustum API pointers
  58.  *
  59.  * Revision 1.11  1997/04/01 04:23:53  brianp
  60.  * added gl_analyze_*_matrix() functions
  61.  *
  62.  * Revision 1.10  1997/02/10 19:47:53  brianp
  63.  * moved buffer resize code out of gl_Viewport() into gl_ResizeBuffersMESA()
  64.  *
  65.  * Revision 1.9  1997/01/31 23:32:40  brianp
  66.  * now clear depth buffer after reallocation due to window resize
  67.  *
  68.  * Revision 1.8  1997/01/29 19:06:04  brianp
  69.  * removed extra, local definition of Identity[] matrix
  70.  *
  71.  * Revision 1.7  1997/01/28 22:19:17  brianp
  72.  * new matrix inversion code from Stephane Rehel
  73.  *
  74.  * Revision 1.6  1996/12/22 17:53:11  brianp
  75.  * faster invert_matrix() function from scotter@iname.com
  76.  *
  77.  * Revision 1.5  1996/12/02 18:58:34  brianp
  78.  * gl_rotation_matrix() now returns identity matrix if given a 0 rotation axis
  79.  *
  80.  * Revision 1.4  1996/09/27 01:29:05  brianp
  81.  * added missing default cases to switches
  82.  *
  83.  * Revision 1.3  1996/09/15 14:18:37  brianp
  84.  * now use GLframebuffer and GLvisual
  85.  *
  86.  * Revision 1.2  1996/09/14 06:46:04  brianp
  87.  * better matmul() from Jacques Leroy
  88.  *
  89.  * Revision 1.1  1996/09/13 01:38:16  brianp
  90.  * Initial revision
  91.  *
  92.  */
  93.  
  94.  
  95. /*
  96.  * Matrix operations
  97.  *
  98.  *
  99.  * NOTES:
  100.  * 1. 4x4 transformation matrices are stored in memory in column major order.
  101.  * 2. Points/vertices are to be thought of as column vectors.
  102.  * 3. Transformation of a point p by a matrix M is: p' = M * p
  103.  *
  104.  */
  105.  
  106.  
  107. #ifdef PC_HEADER
  108. #include "all.h"
  109. #else
  110. #include <math.h>
  111. #include <stdio.h>
  112. #include <stdlib.h>
  113. #include <string.h>
  114. #include "context.h"
  115. #include "dlist.h"
  116. #include "macros.h"
  117. #include "matrix.h"
  118. #include "mmath.h"
  119. #include "types.h"
  120. #endif
  121.  
  122.  
  123.  
  124. static GLfloat Identity[16] = {
  125.    1.0, 0.0, 0.0, 0.0,
  126.    0.0, 1.0, 0.0, 0.0,
  127.    0.0, 0.0, 1.0, 0.0,
  128.    0.0, 0.0, 0.0, 1.0
  129. };
  130.  
  131.  
  132.  
  133.  
  134. static void print_matrix( const GLfloat m[16] )
  135. {
  136.    int i;
  137.  
  138.    for (i=0;i<4;i++) {
  139.       printf("%f %f %f %f\n", m[i], m[4+i], m[8+i], m[12+i] );
  140.    }
  141. }
  142.  
  143.  
  144.  
  145. /*
  146.  * Perform a 4x4 matrix multiplication  (product = a x b).
  147.  * Input:  a, b - matrices to multiply
  148.  * Output:  product - product of a and b
  149.  * WARNING: (product != b) assumed
  150.  * NOTE:    (product == a) allowed    
  151.  */
  152. static void matmul( GLfloat *product, const GLfloat *a, const GLfloat *b )
  153. {
  154.    /* This matmul was contributed by Thomas Malik */
  155.    GLint i;
  156.  
  157. #define A(row,col)  a[(col<<2)+row]
  158. #define B(row,col)  b[(col<<2)+row]
  159. #define P(row,col)  product[(col<<2)+row]
  160.  
  161.    /* i-te Zeile */
  162.    for (i = 0; i < 4; i++) {
  163.       GLfloat ai0=A(i,0),  ai1=A(i,1),  ai2=A(i,2),  ai3=A(i,3);
  164.       P(i,0) = ai0 * B(0,0) + ai1 * B(1,0) + ai2 * B(2,0) + ai3 * B(3,0);
  165.       P(i,1) = ai0 * B(0,1) + ai1 * B(1,1) + ai2 * B(2,1) + ai3 * B(3,1);
  166.       P(i,2) = ai0 * B(0,2) + ai1 * B(1,2) + ai2 * B(2,2) + ai3 * B(3,2);
  167.       P(i,3) = ai0 * B(0,3) + ai1 * B(1,3) + ai2 * B(2,3) + ai3 * B(3,3);
  168.    }
  169.  
  170. #undef A
  171. #undef B
  172. #undef P
  173. }
  174.  
  175.  
  176.  
  177. /*
  178.  * Compute the inverse of a 4x4 matrix.
  179.  *
  180.  * From an algorithm by V. Strassen, 1969, _Numerishe Mathematik_, vol. 13,
  181.  * pp. 354-356.
  182.  * 60 multiplies, 24 additions, 10 subtractions, 8 negations, 2 divisions,
  183.  * 48 assignments, _0_ branches
  184.  *
  185.  * This implementation by Scott McCaskill
  186.  */
  187.  
  188. typedef GLfloat Mat2[2][2];
  189.  
  190. enum {
  191.     M00 = 0, M01 = 4, M02 = 8, M03 = 12,
  192.     M10 = 1, M11 = 5, M12 = 9, M13 = 13,
  193.     M20 = 2, M21 = 6, M22 = 10,M23 = 14,
  194.     M30 = 3, M31 = 7, M32 = 11,M33 = 15
  195. };
  196.  
  197. static void invert_matrix_general( const GLfloat *m, GLfloat *out )
  198. {
  199.    Mat2 r1, r2, r3, r4, r5, r6, r7;
  200.    const GLfloat * A = m;
  201.    GLfloat *       C = out;
  202.    GLfloat one_over_det;
  203.  
  204.    /*
  205.     * A is the 4x4 source matrix (to be inverted).
  206.     * C is the 4x4 destination matrix
  207.     * a11 is the 2x2 matrix in the upper left quadrant of A
  208.     * a12 is the 2x2 matrix in the upper right quadrant of A
  209.     * a21 is the 2x2 matrix in the lower left quadrant of A
  210.     * a22 is the 2x2 matrix in the lower right quadrant of A
  211.     * similarly, cXX are the 2x2 quadrants of the destination matrix
  212.     */
  213.  
  214.    /* R1 = inverse( a11 ) */
  215.    one_over_det = 1.0f / ( ( A[M00] * A[M11] ) - ( A[M10] * A[M01] ) );
  216.    r1[0][0] = one_over_det * A[M11];
  217.    r1[0][1] = one_over_det * -A[M01];
  218.    r1[1][0] = one_over_det * -A[M10];
  219.    r1[1][1] = one_over_det * A[M00];
  220.  
  221.    /* R2 = a21 x R1 */
  222.    r2[0][0] = A[M20] * r1[0][0] + A[M21] * r1[1][0];
  223.    r2[0][1] = A[M20] * r1[0][1] + A[M21] * r1[1][1];
  224.    r2[1][0] = A[M30] * r1[0][0] + A[M31] * r1[1][0];
  225.    r2[1][1] = A[M30] * r1[0][1] + A[M31] * r1[1][1];
  226.  
  227.    /* R3 = R1 x a12 */
  228.    r3[0][0] = r1[0][0] * A[M02] + r1[0][1] * A[M12];
  229.    r3[0][1] = r1[0][0] * A[M03] + r1[0][1] * A[M13];
  230.    r3[1][0] = r1[1][0] * A[M02] + r1[1][1] * A[M12];
  231.    r3[1][1] = r1[1][0] * A[M03] + r1[1][1] * A[M13];
  232.  
  233.    /* R4 = a21 x R3 */
  234.    r4[0][0] = A[M20] * r3[0][0] + A[M21] * r3[1][0];
  235.    r4[0][1] = A[M20] * r3[0][1] + A[M21] * r3[1][1];
  236.    r4[1][0] = A[M30] * r3[0][0] + A[M31] * r3[1][0];
  237.    r4[1][1] = A[M30] * r3[0][1] + A[M31] * r3[1][1];
  238.  
  239.    /* R5 = R4 - a22 */
  240.    r5[0][0] = r4[0][0] - A[M22];
  241.    r5[0][1] = r4[0][1] - A[M23];
  242.    r5[1][0] = r4[1][0] - A[M32];
  243.    r5[1][1] = r4[1][1] - A[M33];
  244.  
  245.    /* R6 = inverse( R5 ) */
  246.    one_over_det = 1.0f / ( ( r5[0][0] * r5[1][1] ) - ( r5[1][0] * r5[0][1] ) );
  247.    r6[0][0] = one_over_det * r5[1][1];
  248.    r6[0][1] = one_over_det * -r5[0][1];
  249.    r6[1][0] = one_over_det * -r5[1][0];
  250.    r6[1][1] = one_over_det * r5[0][0];
  251.  
  252.    /* c12 = R3 x R6 */
  253.    C[M02] = r3[0][0] * r6[0][0] + r3[0][1] * r6[1][0];
  254.    C[M03] = r3[0][0] * r6[0][1] + r3[0][1] * r6[1][1];
  255.    C[M12] = r3[1][0] * r6[0][0] + r3[1][1] * r6[1][0];
  256.    C[M13] = r3[1][0] * r6[0][1] + r3[1][1] * r6[1][1];
  257.  
  258.    /* c21 = R6 x R2 */
  259.    C[M20] = r6[0][0] * r2[0][0] + r6[0][1] * r2[1][0];
  260.    C[M21] = r6[0][0] * r2[0][1] + r6[0][1] * r2[1][1];
  261.    C[M30] = r6[1][0] * r2[0][0] + r6[1][1] * r2[1][0];
  262.    C[M31] = r6[1][0] * r2[0][1] + r6[1][1] * r2[1][1];
  263.  
  264.    /* R7 = R3 x c21 */
  265.    r7[0][0] = r3[0][0] * C[M20] + r3[0][1] * C[M30];
  266.    r7[0][1] = r3[0][0] * C[M21] + r3[0][1] * C[M31];
  267.    r7[1][0] = r3[1][0] * C[M20] + r3[1][1] * C[M30];
  268.    r7[1][1] = r3[1][0] * C[M21] + r3[1][1] * C[M31];
  269.  
  270.    /* c11 = R1 - R7 */
  271.    C[M00] = r1[0][0] - r7[0][0];
  272.    C[M01] = r1[0][1] - r7[0][1];
  273.    C[M10] = r1[1][0] - r7[1][0];
  274.    C[M11] = r1[1][1] - r7[1][1];
  275.  
  276.    /* c22 = -R6 */
  277.    C[M22] = -r6[0][0];
  278.    C[M23] = -r6[0][1];
  279.    C[M32] = -r6[1][0];
  280.    C[M33] = -r6[1][1];
  281. }
  282.  
  283.  
  284. /*
  285.  * Invert matrix m.  This algorithm contributed by Stephane Rehel
  286.  * <rehel@worldnet.fr>
  287.  */
  288. static void invert_matrix( const GLfloat *m, GLfloat *out )
  289. {
  290. /* NB. OpenGL Matrices are COLUMN major. */
  291. #define MAT(m,r,c) (m)[(c)*4+(r)]
  292.  
  293. /* Here's some shorthand converting standard (row,column) to index. */
  294. #define m11 MAT(m,0,0)
  295. #define m12 MAT(m,0,1)
  296. #define m13 MAT(m,0,2)
  297. #define m14 MAT(m,0,3)
  298. #define m21 MAT(m,1,0)
  299. #define m22 MAT(m,1,1)
  300. #define m23 MAT(m,1,2)
  301. #define m24 MAT(m,1,3)
  302. #define m31 MAT(m,2,0)
  303. #define m32 MAT(m,2,1)
  304. #define m33 MAT(m,2,2)
  305. #define m34 MAT(m,2,3)
  306. #define m41 MAT(m,3,0)
  307. #define m42 MAT(m,3,1)
  308. #define m43 MAT(m,3,2)
  309. #define m44 MAT(m,3,3)
  310.  
  311.    register GLfloat det;
  312.    GLfloat tmp[16]; /* Allow out == in. */
  313.  
  314.    if( m41 != 0. || m42 != 0. || m43 != 0. || m44 != 1. ) {
  315.       invert_matrix_general(m, out);
  316.       return;
  317.    }
  318.  
  319.    /* Inverse = adjoint / det. (See linear algebra texts.)*/
  320.  
  321.    tmp[0]= m22 * m33 - m23 * m32;
  322.    tmp[1]= m23 * m31 - m21 * m33;
  323.    tmp[2]= m21 * m32 - m22 * m31;
  324.  
  325.    /* Compute determinant as early as possible using these cofactors. */
  326.    det= m11 * tmp[0] + m12 * tmp[1] + m13 * tmp[2];
  327.  
  328.    /* Run singularity test. */
  329.    if (det == 0.0F) {
  330.       /* printf("invert_matrix: Warning: Singular matrix.\n"); */
  331.       MEMCPY( out, Identity, 16*sizeof(GLfloat) );
  332.    }
  333.    else {
  334.       GLfloat d12, d13, d23, d24, d34, d41;
  335.       register GLfloat im11, im12, im13, im14;
  336.  
  337.       det= 1. / det;
  338.  
  339.       /* Compute rest of inverse. */
  340.       tmp[0] *= det;
  341.       tmp[1] *= det;
  342.       tmp[2] *= det;
  343.       tmp[3]  = 0.;
  344.  
  345.       im11= m11 * det;
  346.       im12= m12 * det;
  347.       im13= m13 * det;
  348.       im14= m14 * det;
  349.       tmp[4] = im13 * m32 - im12 * m33;
  350.       tmp[5] = im11 * m33 - im13 * m31;
  351.       tmp[6] = im12 * m31 - im11 * m32;
  352.       tmp[7] = 0.;
  353.  
  354.       /* Pre-compute 2x2 dets for first two rows when computing */
  355.       /* cofactors of last two rows. */
  356.       d12 = im11*m22 - m21*im12;
  357.       d13 = im11*m23 - m21*im13;
  358.       d23 = im12*m23 - m22*im13;
  359.       d24 = im12*m24 - m22*im14;
  360.       d34 = im13*m24 - m23*im14;
  361.       d41 = im14*m21 - m24*im11;
  362.  
  363.       tmp[8] =  d23;
  364.       tmp[9] = -d13;
  365.       tmp[10] = d12;
  366.       tmp[11] = 0.;
  367.  
  368.       tmp[12] = -(m32 * d34 - m33 * d24 + m34 * d23);
  369.       tmp[13] =  (m31 * d34 + m33 * d41 + m34 * d13);
  370.       tmp[14] = -(m31 * d24 + m32 * d41 + m34 * d12);
  371.       tmp[15] =  1.;
  372.  
  373.       MEMCPY(out, tmp, 16*sizeof(GLfloat));
  374.   }
  375.  
  376. #undef m11
  377. #undef m12
  378. #undef m13
  379. #undef m14
  380. #undef m21
  381. #undef m22
  382. #undef m23
  383. #undef m24
  384. #undef m31
  385. #undef m32
  386. #undef m33
  387. #undef m34
  388. #undef m41
  389. #undef m42
  390. #undef m43
  391. #undef m44
  392. #undef MAT
  393. }
  394.  
  395.  
  396.  
  397. /*
  398.  * Determine if the given matrix is the identity matrix.
  399.  */
  400. static GLboolean is_identity( const GLfloat m[16] )
  401. {
  402.    if (   m[0]==1.0F && m[4]==0.0F && m[ 8]==0.0F && m[12]==0.0F
  403.        && m[1]==0.0F && m[5]==1.0F && m[ 9]==0.0F && m[13]==0.0F
  404.        && m[2]==0.0F && m[6]==0.0F && m[10]==1.0F && m[14]==0.0F
  405.        && m[3]==0.0F && m[7]==0.0F && m[11]==0.0F && m[15]==1.0F) {
  406.       return GL_TRUE;
  407.    }
  408.    else {
  409.       return GL_FALSE;
  410.    }
  411. }
  412.  
  413.  
  414. /*
  415.  * Examine the current modelview matrix to determine its type.
  416.  * Later we use the matrix type to optimize vertex transformations.
  417.  */
  418. void gl_analyze_modelview_matrix( GLcontext *ctx )
  419. {
  420.    const GLfloat *m = ctx->ModelViewMatrix;
  421.    if (is_identity(m)) {
  422.       ctx->ModelViewMatrixType = MATRIX_IDENTITY;
  423.    }
  424.    else if (                 m[4]==0.0F && m[ 8]==0.0F               
  425.         && m[1]==0.0F               && m[ 9]==0.0F
  426.         && m[2]==0.0F && m[6]==0.0F && m[10]==1.0F && m[14]==0.0F
  427.         && m[3]==0.0F && m[7]==0.0F && m[11]==0.0F && m[15]==1.0F) {
  428.       ctx->ModelViewMatrixType = MATRIX_2D_NO_ROT;
  429.    }
  430.    else if (                               m[ 8]==0.0F               
  431.         &&                             m[ 9]==0.0F
  432.         && m[2]==0.0F && m[6]==0.0F && m[10]==1.0F && m[14]==0.0F
  433.         && m[3]==0.0F && m[7]==0.0F && m[11]==0.0F && m[15]==1.0F) {
  434.       ctx->ModelViewMatrixType = MATRIX_2D;
  435.    }
  436.    else if (m[3]==0.0F && m[7]==0.0F && m[11]==0.0F && m[15]==1.0F) {
  437.       ctx->ModelViewMatrixType = MATRIX_3D;
  438.    }
  439.    else {
  440.       ctx->ModelViewMatrixType = MATRIX_GENERAL;
  441.    }
  442.  
  443.    invert_matrix( ctx->ModelViewMatrix, ctx->ModelViewInv );
  444.    ctx->NewModelViewMatrix = GL_FALSE;
  445. }
  446.  
  447.  
  448.  
  449. /*
  450.  * Examine the current projection matrix to determine its type.
  451.  * Later we use the matrix type to optimize vertex transformations.
  452.  */
  453. void gl_analyze_projection_matrix( GLcontext *ctx )
  454. {
  455.    /* look for common-case ortho and perspective matrices */
  456.    const GLfloat *m = ctx->ProjectionMatrix;
  457.    if (is_identity(m)) {
  458.       ctx->ProjectionMatrixType = MATRIX_IDENTITY;
  459.    }
  460.    else if (                 m[4]==0.0F && m[8] ==0.0F
  461.         && m[1]==0.0F               && m[9] ==0.0F
  462.         && m[2]==0.0F && m[6]==0.0F
  463.         && m[3]==0.0F && m[7]==0.0F && m[11]==0.0F && m[15]==1.0F) {
  464.       ctx->ProjectionMatrixType = MATRIX_ORTHO;
  465.    }
  466.    else if (                 m[4]==0.0F                 && m[12]==0.0F
  467.         && m[1]==0.0F                               && m[13]==0.0F
  468.         && m[2]==0.0F && m[6]==0.0F
  469.         && m[3]==0.0F && m[7]==0.0F && m[11]==-1.0F && m[15]==0.0F) {
  470.       ctx->ProjectionMatrixType = MATRIX_PERSPECTIVE;
  471.    }
  472.    else {
  473.       ctx->ProjectionMatrixType = MATRIX_GENERAL;
  474.    }
  475.  
  476.    ctx->NewProjectionMatrix = GL_FALSE;
  477. }
  478.  
  479.  
  480.  
  481. /*
  482.  * Examine the current texture matrix to determine its type.
  483.  * Later we use the matrix type to optimize texture coordinate transformations.
  484.  */
  485. void gl_analyze_texture_matrix( GLcontext *ctx )
  486. {
  487.    const GLfloat *m = ctx->TextureMatrix;
  488.    if (is_identity(m)) {
  489.       ctx->TextureMatrixType = MATRIX_IDENTITY;
  490.    }
  491.    else if (                               m[ 8]==0.0F               
  492.         &&                             m[ 9]==0.0F
  493.         && m[2]==0.0F && m[6]==0.0F && m[10]==1.0F && m[14]==0.0F
  494.         && m[3]==0.0F && m[7]==0.0F && m[11]==0.0F && m[15]==1.0F) {
  495.       ctx->TextureMatrixType = MATRIX_2D;
  496.    }
  497.    else if (m[3]==0.0F && m[7]==0.0F && m[11]==0.0F && m[15]==1.0F) {
  498.       ctx->TextureMatrixType = MATRIX_3D;
  499.    }
  500.    else {
  501.       ctx->TextureMatrixType = MATRIX_GENERAL;
  502.    }
  503.  
  504.    ctx->NewTextureMatrix = GL_FALSE;
  505. }
  506.  
  507.  
  508.  
  509. void gl_Frustum( GLcontext *ctx,
  510.          GLdouble left, GLdouble right,
  511.          GLdouble bottom, GLdouble top,
  512.          GLdouble nearval, GLdouble farval )
  513. {
  514.    GLfloat x, y, a, b, c, d;
  515.    GLfloat m[16];
  516.  
  517.    if (nearval<=0.0 || farval<=0.0) {
  518.       gl_error( ctx,  GL_INVALID_VALUE, "glFrustum(near or far)" );
  519.    }
  520.  
  521.    x = (2.0*nearval) / (right-left);
  522.    y = (2.0*nearval) / (top-bottom);
  523.    a = (right+left) / (right-left);
  524.    b = (top+bottom) / (top-bottom);
  525.    c = -(farval+nearval) / ( farval-nearval);
  526.    d = -(2.0*farval*nearval) / (farval-nearval);  /* error? */
  527.  
  528. #define M(row,col)  m[col*4+row]
  529.    M(0,0) = x;     M(0,1) = 0.0F;  M(0,2) = a;      M(0,3) = 0.0F;
  530.    M(1,0) = 0.0F;  M(1,1) = y;     M(1,2) = b;      M(1,3) = 0.0F;
  531.    M(2,0) = 0.0F;  M(2,1) = 0.0F;  M(2,2) = c;      M(2,3) = d;
  532.    M(3,0) = 0.0F;  M(3,1) = 0.0F;  M(3,2) = -1.0F;  M(3,3) = 0.0F;
  533. #undef M
  534.  
  535.    gl_MultMatrixf( ctx, m );
  536.  
  537.  
  538.    /* Need to keep a stack of near/far values in case the user push/pops
  539.     * the projection matrix stack so that we can call Driver.NearFar()
  540.     * after a pop.
  541.     */
  542.    ctx->NearFarStack[ctx->ProjectionStackDepth][0] = nearval;
  543.    ctx->NearFarStack[ctx->ProjectionStackDepth][1] = farval;
  544.  
  545.    if (ctx->Driver.NearFar) {
  546.       (*ctx->Driver.NearFar)( ctx, nearval, farval );
  547.    }
  548. }
  549.  
  550.  
  551. void gl_Ortho( GLcontext *ctx,
  552.            GLdouble left, GLdouble right,
  553.            GLdouble bottom, GLdouble top,
  554.            GLdouble nearval, GLdouble farval )
  555. {
  556.    GLfloat x, y, z;
  557.    GLfloat tx, ty, tz;
  558.    GLfloat m[16];
  559.  
  560.    x = 2.0 / (right-left);
  561.    y = 2.0 / (top-bottom);
  562.    z = -2.0 / (farval-nearval);
  563.    tx = -(right+left) / (right-left);
  564.    ty = -(top+bottom) / (top-bottom);
  565.    tz = -(farval+nearval) / (farval-nearval);
  566.  
  567. #define M(row,col)  m[col*4+row]
  568.    M(0,0) = x;     M(0,1) = 0.0F;  M(0,2) = 0.0F;  M(0,3) = tx;
  569.    M(1,0) = 0.0F;  M(1,1) = y;     M(1,2) = 0.0F;  M(1,3) = ty;
  570.    M(2,0) = 0.0F;  M(2,1) = 0.0F;  M(2,2) = z;     M(2,3) = tz;
  571.    M(3,0) = 0.0F;  M(3,1) = 0.0F;  M(3,2) = 0.0F;  M(3,3) = 1.0F;
  572. #undef M
  573.  
  574.    gl_MultMatrixf( ctx, m );
  575.  
  576.    if (ctx->Driver.NearFar) {
  577.       (*ctx->Driver.NearFar)( ctx, nearval, farval );
  578.    }
  579. }
  580.  
  581.  
  582. void gl_MatrixMode( GLcontext *ctx, GLenum mode )
  583. {
  584.    if (INSIDE_BEGIN_END(ctx)) {
  585.       gl_error( ctx,  GL_INVALID_OPERATION, "glMatrixMode" );
  586.       return;
  587.    }
  588.    switch (mode) {
  589.       case GL_MODELVIEW:
  590.       case GL_PROJECTION:
  591.       case GL_TEXTURE:
  592.      ctx->Transform.MatrixMode = mode;
  593.      break;
  594.       default:
  595.      gl_error( ctx,  GL_INVALID_ENUM, "glMatrixMode" );
  596.    }
  597. }
  598.  
  599.  
  600.  
  601. void gl_PushMatrix( GLcontext *ctx )
  602. {
  603.    if (INSIDE_BEGIN_END(ctx)) {
  604.       gl_error( ctx,  GL_INVALID_OPERATION, "glPushMatrix" );
  605.       return;
  606.    }
  607.    switch (ctx->Transform.MatrixMode) {
  608.       case GL_MODELVIEW:
  609.      if (ctx->ModelViewStackDepth>=MAX_MODELVIEW_STACK_DEPTH-1) {
  610.         gl_error( ctx,  GL_STACK_OVERFLOW, "glPushMatrix");
  611.         return;
  612.      }
  613.      MEMCPY( ctx->ModelViewStack[ctx->ModelViewStackDepth],
  614.          ctx->ModelViewMatrix,
  615.          16*sizeof(GLfloat) );
  616.      ctx->ModelViewStackDepth++;
  617.      break;
  618.       case GL_PROJECTION:
  619.      if (ctx->ProjectionStackDepth>=MAX_PROJECTION_STACK_DEPTH) {
  620.         gl_error( ctx,  GL_STACK_OVERFLOW, "glPushMatrix");
  621.         return;
  622.      }
  623.      MEMCPY( ctx->ProjectionStack[ctx->ProjectionStackDepth],
  624.          ctx->ProjectionMatrix,
  625.          16*sizeof(GLfloat) );
  626.      ctx->ProjectionStackDepth++;
  627.  
  628.      /* Save near and far projection values */
  629.      ctx->NearFarStack[ctx->ProjectionStackDepth][0]
  630.         = ctx->NearFarStack[ctx->ProjectionStackDepth-1][0];
  631.      ctx->NearFarStack[ctx->ProjectionStackDepth][1]
  632.         = ctx->NearFarStack[ctx->ProjectionStackDepth-1][1];
  633.      break;
  634.       case GL_TEXTURE:
  635.      if (ctx->TextureStackDepth>=MAX_TEXTURE_STACK_DEPTH) {
  636.         gl_error( ctx,  GL_STACK_OVERFLOW, "glPushMatrix");
  637.         return;
  638.      }
  639.      MEMCPY( ctx->TextureStack[ctx->TextureStackDepth],
  640.          ctx->TextureMatrix,
  641.          16*sizeof(GLfloat) );
  642.      ctx->TextureStackDepth++;
  643.      break;
  644.       default:
  645.      gl_problem(ctx, "Bad matrix mode in gl_PushMatrix");
  646.    }
  647. }
  648.  
  649.  
  650.  
  651. void gl_PopMatrix( GLcontext *ctx )
  652. {
  653.    if (INSIDE_BEGIN_END(ctx)) {
  654.       gl_error( ctx,  GL_INVALID_OPERATION, "glPopMatrix" );
  655.       return;
  656.    }
  657.    switch (ctx->Transform.MatrixMode) {
  658.       case GL_MODELVIEW:
  659.      if (ctx->ModelViewStackDepth==0) {
  660.         gl_error( ctx,  GL_STACK_UNDERFLOW, "glPopMatrix");
  661.         return;
  662.      }
  663.      ctx->ModelViewStackDepth--;
  664.      MEMCPY( ctx->ModelViewMatrix,
  665.          ctx->ModelViewStack[ctx->ModelViewStackDepth],
  666.          16*sizeof(GLfloat) );
  667.      ctx->NewModelViewMatrix = GL_TRUE;
  668.      break;
  669.       case GL_PROJECTION:
  670.      if (ctx->ProjectionStackDepth==0) {
  671.         gl_error( ctx,  GL_STACK_UNDERFLOW, "glPopMatrix");
  672.         return;
  673.      }
  674.      ctx->ProjectionStackDepth--;
  675.      MEMCPY( ctx->ProjectionMatrix,
  676.          ctx->ProjectionStack[ctx->ProjectionStackDepth],
  677.          16*sizeof(GLfloat) );
  678.      ctx->NewProjectionMatrix = GL_TRUE;
  679.  
  680.      /* Device driver near/far values */
  681.      {
  682.         GLfloat nearVal = ctx->NearFarStack[ctx->ProjectionStackDepth][0];
  683.         GLfloat farVal  = ctx->NearFarStack[ctx->ProjectionStackDepth][1];
  684.         if (ctx->Driver.NearFar) {
  685.            (*ctx->Driver.NearFar)( ctx, nearVal, farVal );
  686.         }
  687.      }
  688.      break;
  689.       case GL_TEXTURE:
  690.      if (ctx->TextureStackDepth==0) {
  691.         gl_error( ctx,  GL_STACK_UNDERFLOW, "glPopMatrix");
  692.         return;
  693.      }
  694.      ctx->TextureStackDepth--;
  695.      MEMCPY( ctx->TextureMatrix,
  696.          ctx->TextureStack[ctx->TextureStackDepth],
  697.          16*sizeof(GLfloat) );
  698.      ctx->NewTextureMatrix = GL_TRUE;
  699.      break;
  700.       default:
  701.      gl_problem(ctx, "Bad matrix mode in gl_PopMatrix");
  702.    }
  703. }
  704.  
  705.  
  706.  
  707. void gl_LoadIdentity( GLcontext *ctx )
  708. {
  709.    if (INSIDE_BEGIN_END(ctx)) {
  710.       gl_error( ctx,  GL_INVALID_OPERATION, "glLoadIdentity" );
  711.       return;
  712.    }
  713.    switch (ctx->Transform.MatrixMode) {
  714.       case GL_MODELVIEW:
  715.      MEMCPY( ctx->ModelViewMatrix, Identity, 16*sizeof(GLfloat) );
  716.      MEMCPY( ctx->ModelViewInv, Identity, 16*sizeof(GLfloat) );
  717.      ctx->ModelViewMatrixType = MATRIX_IDENTITY;
  718.      ctx->NewModelViewMatrix = GL_FALSE;
  719.      break;
  720.       case GL_PROJECTION:
  721.      MEMCPY( ctx->ProjectionMatrix, Identity, 16*sizeof(GLfloat) );
  722.      ctx->ProjectionMatrixType = MATRIX_IDENTITY;
  723.      ctx->NewProjectionMatrix = GL_FALSE;
  724.      break;
  725.       case GL_TEXTURE:
  726.      MEMCPY( ctx->TextureMatrix, Identity, 16*sizeof(GLfloat) );
  727.      ctx->TextureMatrixType = MATRIX_IDENTITY;
  728.      ctx->NewTextureMatrix = GL_FALSE;
  729.      break;
  730.       default:
  731.      gl_problem(ctx, "Bad matrix mode in gl_LoadIdentity");
  732.    }
  733. }
  734.  
  735.  
  736. void gl_LoadMatrixf( GLcontext *ctx, const GLfloat *m )
  737. {
  738.    if (INSIDE_BEGIN_END(ctx)) {
  739.       gl_error( ctx,  GL_INVALID_OPERATION, "glLoadMatrix" );
  740.       return;
  741.    }
  742.    switch (ctx->Transform.MatrixMode) {
  743.       case GL_MODELVIEW:
  744.      MEMCPY( ctx->ModelViewMatrix, m, 16*sizeof(GLfloat) );
  745.      ctx->NewModelViewMatrix = GL_TRUE;
  746.      break;
  747.       case GL_PROJECTION:
  748.      MEMCPY( ctx->ProjectionMatrix, m, 16*sizeof(GLfloat) );
  749.      ctx->NewProjectionMatrix = GL_TRUE;
  750.      break;
  751.       case GL_TEXTURE:
  752.      MEMCPY( ctx->TextureMatrix, m, 16*sizeof(GLfloat) );
  753.      ctx->NewTextureMatrix = GL_TRUE;
  754.      break;
  755.       default:
  756.      gl_problem(ctx, "Bad matrix mode in gl_LoadMatrixf");
  757.    }
  758. }
  759.  
  760.  
  761.  
  762. void gl_MultMatrixf( GLcontext *ctx, const GLfloat *m )
  763. {
  764.    if (INSIDE_BEGIN_END(ctx)) {
  765.       gl_error( ctx,  GL_INVALID_OPERATION, "glMultMatrix" );
  766.       return;
  767.    }
  768.    switch (ctx->Transform.MatrixMode) {
  769.       case GL_MODELVIEW:
  770.      matmul( ctx->ModelViewMatrix, ctx->ModelViewMatrix, m );
  771.      ctx->NewModelViewMatrix = GL_TRUE;
  772.      break;
  773.       case GL_PROJECTION:
  774.      matmul( ctx->ProjectionMatrix, ctx->ProjectionMatrix, m );
  775.      ctx->NewProjectionMatrix = GL_TRUE;
  776.      break;
  777.       case GL_TEXTURE:
  778.      matmul( ctx->TextureMatrix, ctx->TextureMatrix, m );
  779.      ctx->NewTextureMatrix = GL_TRUE;
  780.      break;
  781.       default:
  782.      gl_problem(ctx, "Bad matrix mode in gl_MultMatrixf");
  783.    }
  784. }
  785.  
  786.  
  787.  
  788. /*
  789.  * Generate a 4x4 transformation matrix from glRotate parameters.
  790.  */
  791. void gl_rotation_matrix( GLfloat angle, GLfloat x, GLfloat y, GLfloat z,
  792.              GLfloat m[] )
  793. {
  794.    /* This function contributed by Erich Boleyn (erich@uruk.org) */
  795.    GLfloat mag, s, c;
  796.    GLfloat xx, yy, zz, xy, yz, zx, xs, ys, zs, one_c;
  797.  
  798.    s = sin( angle * DEG2RAD );
  799.    c = cos( angle * DEG2RAD );
  800.  
  801.    mag = GL_SQRT( x*x + y*y + z*z );
  802.  
  803.    if (mag == 0.0) {
  804.       /* generate an identity matrix and return */
  805.       MEMCPY(m, Identity, sizeof(GLfloat)*16);
  806.       return;
  807.    }
  808.  
  809.    x /= mag;
  810.    y /= mag;
  811.    z /= mag;
  812.  
  813. #define M(row,col)  m[col*4+row]
  814.  
  815.    /*
  816.     *     Arbitrary axis rotation matrix.
  817.     *
  818.     *  This is composed of 5 matrices, Rz, Ry, T, Ry', Rz', multiplied
  819.     *  like so:  Rz * Ry * T * Ry' * Rz'.  T is the final rotation
  820.     *  (which is about the X-axis), and the two composite transforms
  821.     *  Ry' * Rz' and Rz * Ry are (respectively) the rotations necessary
  822.     *  from the arbitrary axis to the X-axis then back.  They are
  823.     *  all elementary rotations.
  824.     *
  825.     *  Rz' is a rotation about the Z-axis, to bring the axis vector
  826.     *  into the x-z plane.  Then Ry' is applied, rotating about the
  827.     *  Y-axis to bring the axis vector parallel with the X-axis.  The
  828.     *  rotation about the X-axis is then performed.  Ry and Rz are
  829.     *  simply the respective inverse transforms to bring the arbitrary
  830.     *  axis back to it's original orientation.  The first transforms
  831.     *  Rz' and Ry' are considered inverses, since the data from the
  832.     *  arbitrary axis gives you info on how to get to it, not how
  833.     *  to get away from it, and an inverse must be applied.
  834.     *
  835.     *  The basic calculation used is to recognize that the arbitrary
  836.     *  axis vector (x, y, z), since it is of unit length, actually
  837.     *  represents the sines and cosines of the angles to rotate the
  838.     *  X-axis to the same orientation, with theta being the angle about
  839.     *  Z and phi the angle about Y (in the order described above)
  840.     *  as follows:
  841.     *
  842.     *  cos ( theta ) = x / sqrt ( 1 - z^2 )
  843.     *  sin ( theta ) = y / sqrt ( 1 - z^2 )
  844.     *
  845.     *  cos ( phi ) = sqrt ( 1 - z^2 )
  846.     *  sin ( phi ) = z
  847.     *
  848.     *  Note that cos ( phi ) can further be inserted to the above
  849.     *  formulas:
  850.     *
  851.     *  cos ( theta ) = x / cos ( phi )
  852.     *  sin ( theta ) = y / sin ( phi )
  853.     *
  854.     *  ...etc.  Because of those relations and the standard trigonometric
  855.     *  relations, it is pssible to reduce the transforms down to what
  856.     *  is used below.  It may be that any primary axis chosen will give the
  857.     *  same results (modulo a sign convention) using thie method.
  858.     *
  859.     *  Particularly nice is to notice that all divisions that might
  860.     *  have caused trouble when parallel to certain planes or
  861.     *  axis go away with care paid to reducing the expressions.
  862.     *  After checking, it does perform correctly under all cases, since
  863.     *  in all the cases of division where the denominator would have
  864.     *  been zero, the numerator would have been zero as well, giving
  865.     *  the expected result.
  866.     */
  867.  
  868.    xx = x * x;
  869.    yy = y * y;
  870.    zz = z * z;
  871.    xy = x * y;
  872.    yz = y * z;
  873.    zx = z * x;
  874.    xs = x * s;
  875.    ys = y * s;
  876.    zs = z * s;
  877.    one_c = 1.0F - c;
  878.  
  879.    M(0,0) = (one_c * xx) + c;
  880.    M(0,1) = (one_c * xy) - zs;
  881.    M(0,2) = (one_c * zx) + ys;
  882.    M(0,3) = 0.0F;
  883.  
  884.    M(1,0) = (one_c * xy) + zs;
  885.    M(1,1) = (one_c * yy) + c;
  886.    M(1,2) = (one_c * yz) - xs;
  887.    M(1,3) = 0.0F;
  888.  
  889.    M(2,0) = (one_c * zx) - ys;
  890.    M(2,1) = (one_c * yz) + xs;
  891.    M(2,2) = (one_c * zz) + c;
  892.    M(2,3) = 0.0F;
  893.  
  894.    M(3,0) = 0.0F;
  895.    M(3,1) = 0.0F;
  896.    M(3,2) = 0.0F;
  897.    M(3,3) = 1.0F;
  898.  
  899. #undef M
  900. }
  901.  
  902.  
  903.  
  904. void gl_Rotatef( GLcontext *ctx,
  905.          GLfloat angle, GLfloat x, GLfloat y, GLfloat z )
  906. {
  907.    GLfloat m[16];
  908.    gl_rotation_matrix( angle, x, y, z, m );
  909.    gl_MultMatrixf( ctx, m );
  910. }
  911.  
  912.  
  913.  
  914. /*
  915.  * Execute a glScale call
  916.  */
  917. void gl_Scalef( GLcontext *ctx, GLfloat x, GLfloat y, GLfloat z )
  918. {
  919.    GLfloat *m;
  920.  
  921.    if (INSIDE_BEGIN_END(ctx)) {
  922.       gl_error( ctx,  GL_INVALID_OPERATION, "glScale" );
  923.       return;
  924.    }
  925.    switch (ctx->Transform.MatrixMode) {
  926.       case GL_MODELVIEW:
  927.      m = ctx->ModelViewMatrix;
  928.      ctx->NewModelViewMatrix = GL_TRUE;
  929.      break;
  930.       case GL_PROJECTION:
  931.      m = ctx->ProjectionMatrix;
  932.      ctx->NewProjectionMatrix = GL_TRUE;
  933.      break;
  934.       case GL_TEXTURE:
  935.      m = ctx->TextureMatrix;
  936.      ctx->NewTextureMatrix = GL_TRUE;
  937.      break;
  938.       default:
  939.      gl_problem(ctx, "Bad matrix mode in gl_Scalef");
  940.      return;
  941.    }
  942.    m[0] *= x;   m[4] *= y;   m[8]  *= z;
  943.    m[1] *= x;   m[5] *= y;   m[9]  *= z;
  944.    m[2] *= x;   m[6] *= y;   m[10] *= z;
  945.    m[3] *= x;   m[7] *= y;   m[11] *= z;
  946. }
  947.  
  948.  
  949.  
  950. /*
  951.  * Execute a glTranslate call
  952.  */
  953. void gl_Translatef( GLcontext *ctx, GLfloat x, GLfloat y, GLfloat z )
  954. {
  955.    GLfloat *m;
  956.    if (INSIDE_BEGIN_END(ctx)) {
  957.       gl_error( ctx, GL_INVALID_OPERATION, "glTranslate" );
  958.       return;
  959.    }
  960.    switch (ctx->Transform.MatrixMode) {
  961.       case GL_MODELVIEW:
  962.      m = ctx->ModelViewMatrix;
  963.      ctx->NewModelViewMatrix = GL_TRUE;
  964.      break;
  965.       case GL_PROJECTION:
  966.      m = ctx->ProjectionMatrix;
  967.      ctx->NewProjectionMatrix = GL_TRUE;
  968.      break;
  969.       case GL_TEXTURE:
  970.      m = ctx->TextureMatrix;
  971.      ctx->NewTextureMatrix = GL_TRUE;
  972.      break;
  973.       default:
  974.      gl_problem(ctx, "Bad matrix mode in gl_Translatef");
  975.      return;
  976.    }
  977.  
  978.    m[12] = m[0] * x + m[4] * y + m[8]  * z + m[12];
  979.    m[13] = m[1] * x + m[5] * y + m[9]  * z + m[13];
  980.    m[14] = m[2] * x + m[6] * y + m[10] * z + m[14];
  981.    m[15] = m[3] * x + m[7] * y + m[11] * z + m[15];
  982. }
  983.  
  984.  
  985.  
  986.  
  987. /*
  988.  * Define a new viewport and reallocate auxillary buffers if the size of
  989.  * the window (color buffer) has changed.
  990.  */
  991. void gl_Viewport( GLcontext *ctx,
  992.           GLint x, GLint y, GLsizei width, GLsizei height )
  993. {
  994.    if (width<0 || height<0) {
  995.       gl_error( ctx,  GL_INVALID_VALUE, "glViewport" );
  996.       return;
  997.    }
  998.    if (INSIDE_BEGIN_END(ctx)) {
  999.       gl_error( ctx,  GL_INVALID_OPERATION, "glViewport" );
  1000.       return;
  1001.    }
  1002.  
  1003.    /* clamp width, and height to implementation dependent range */
  1004.    width  = CLAMP( width,  1, MAX_WIDTH );
  1005.    height = CLAMP( height, 1, MAX_HEIGHT );
  1006.    /* Save viewport */
  1007.    ctx->Viewport.X = x;
  1008.    ctx->Viewport.Width = width;
  1009.    ctx->Viewport.Y = y;
  1010.    ctx->Viewport.Height = height;
  1011.  
  1012.    /* compute scale and bias values */
  1013.    ctx->Viewport.Sx = (GLfloat) width / 2.0F;
  1014.    ctx->Viewport.Tx = ctx->Viewport.Sx + x;
  1015.    ctx->Viewport.Sy = (GLfloat) height / 2.0F;
  1016.    ctx->Viewport.Ty = ctx->Viewport.Sy + y;
  1017.  
  1018.    ctx->NewState |= NEW_ALL;   /* just to be safe */
  1019.  
  1020.    /* Check if window/buffer has been resized and if so, reallocate the
  1021.     * ancillary buffers.
  1022.     */
  1023.    gl_ResizeBuffersMESA(ctx);
  1024. }
  1025.