home *** CD-ROM | disk | FTP | other *** search
/ Amiga MA Magazine 1998 #6 / amigamamagazinepolishissue1998.iso / coders / mesa-1.2.8 / src-glu / project.c < prev    next >
C/C++ Source or Header  |  1996-05-27  |  7KB  |  261 lines
  1. /* project.c */
  2.  
  3. /*
  4.  * Mesa 3-D graphics library
  5.  * Version:  1.2
  6.  * Copyright (C) 1995  Brian Paul  (brianp@ssec.wisc.edu)
  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. $Id: project.c,v 1.10 1996/05/15 15:39:16 brianp Exp $
  26.  
  27. $Log: project.c,v $
  28.  * Revision 1.10  1996/05/15  15:39:16  brianp
  29.  * changed some return types from int to GLint
  30.  *
  31.  * Revision 1.9  1995/07/26  15:05:58  brianp
  32.  * replaced memcpy() calls with MEMCPY() macro
  33.  *
  34.  * Revision 1.8  1995/05/22  16:56:20  brianp
  35.  * Release 1.2
  36.  *
  37.  * Revision 1.7  1995/05/16  19:17:21  brianp
  38.  * minor changes to allow compilation with real OpenGL headers
  39.  *
  40.  * Revision 1.6  1995/04/28  14:38:36  brianp
  41.  * added return statement to project/unproject functions
  42.  *
  43.  * Revision 1.5  1995/03/10  20:03:35  brianp
  44.  * fixed -y bug in gluUnProject and gluProject
  45.  *
  46.  * Revision 1.4  1995/03/10  17:01:56  brianp
  47.  * new matmul and invert_matrix function from Thomas Malik
  48.  *
  49.  * Revision 1.3  1995/03/04  19:39:18  brianp
  50.  * version 1.1 beta
  51.  *
  52.  * Revision 1.2  1995/02/24  15:54:46  brianp
  53.  * converted all GLfloats to GLdoubles
  54.  *
  55.  * Revision 1.1  1995/02/24  15:47:09  brianp
  56.  * Initial revision
  57.  *
  58.  */
  59.  
  60.  
  61. #include <stdio.h>
  62. #include <string.h>
  63. #include <math.h>
  64. #include "gluP.h"
  65.  
  66.  
  67. /*
  68.  * This code was contributed by Marc Buffat (buffat@mecaflu.ec-lyon.fr).
  69.  * Thanks Marc!!!
  70.  */
  71.  
  72.  
  73.  
  74. /* implementation de gluProject et gluUnproject */
  75. /* M. Buffat 17/2/95 */
  76.  
  77.  
  78.  
  79. /*
  80.  * Transform a point (column vector) by a 4x4 matrix.  I.e.  out = m * in
  81.  * Input:  m - the 4x4 matrix
  82.  *         in - the 4x1 vector
  83.  * Output:  out - the resulting 4x1 vector.
  84.  */
  85. static void transform_point( GLdouble out[4], const GLdouble m[16],
  86.                  const GLdouble in[4] )
  87. {
  88. #define M(row,col)  m[col*4+row]
  89.    out[0] = M(0,0) * in[0] + M(0,1) * in[1] + M(0,2) * in[2] + M(0,3) * in[3];
  90.    out[1] = M(1,0) * in[0] + M(1,1) * in[1] + M(1,2) * in[2] + M(1,3) * in[3];
  91.    out[2] = M(2,0) * in[0] + M(2,1) * in[1] + M(2,2) * in[2] + M(2,3) * in[3];
  92.    out[3] = M(3,0) * in[0] + M(3,1) * in[1] + M(3,2) * in[2] + M(3,3) * in[3];
  93. #undef M
  94. }
  95.  
  96.  
  97.  
  98.  
  99. /*
  100.  * Perform a 4x4 matrix multiplication  (product = a x b).
  101.  * Input:  a, b - matrices to multiply
  102.  * Output:  product - product of a and b
  103.  */
  104. static void matmul( GLdouble *product, const GLdouble *a, const GLdouble *b )
  105. {
  106.    /* This matmul was contributed by Thomas Malik */
  107.    GLdouble temp[16];
  108.    GLint i;
  109.  
  110. #define A(row,col)  a[(col<<2)+row]
  111. #define B(row,col)  b[(col<<2)+row]
  112. #define T(row,col)  temp[(col<<2)+row]
  113.  
  114.    /* i-te Zeile */
  115.    for (i = 0; i < 4; i++)
  116.      {
  117.     T(i, 0) = A(i, 0) * B(0, 0) + A(i, 1) * B(1, 0) + A(i, 2) * B(2, 0) + A(i, 3) * B(3, 0);
  118.     T(i, 1) = A(i, 0) * B(0, 1) + A(i, 1) * B(1, 1) + A(i, 2) * B(2, 1) + A(i, 3) * B(3, 1);
  119.     T(i, 2) = A(i, 0) * B(0, 2) + A(i, 1) * B(1, 2) + A(i, 2) * B(2, 2) + A(i, 3) * B(3, 2);
  120.     T(i, 3) = A(i, 0) * B(0, 3) + A(i, 1) * B(1, 3) + A(i, 2) * B(2, 3) + A(i, 3) * B(3, 3);
  121.      }
  122.  
  123. #undef A
  124. #undef B
  125. #undef T
  126.    MEMCPY( product, temp, 16*sizeof(GLdouble) );
  127. }
  128.  
  129.  
  130.  
  131.  
  132. /*
  133.  * Find the inverse of the 4 by 4 matrix b using gausian elimination
  134.  * and return it in a.
  135.  *
  136.  * This function was contributed by Thomas Malik (malik@rhrk.uni-kl.de).
  137.  * Thanks Thomas!
  138.  */
  139. static void invert_matrix(const GLdouble *b,GLdouble * a)
  140. {
  141.   static GLdouble identity[16] =
  142.     {
  143.       1.0, 0.0, 0.0, 0.0,
  144.       0.0, 1.0, 0.0, 0.0,
  145.       0.0, 0.0, 1.0, 0.0,
  146.       0.0, 0.0, 0.0, 1.0
  147.     };
  148.  
  149. #define MAT(m,r,c) ((m)[(c)*4+(r)])
  150.  
  151.   GLdouble val, val2;
  152.   GLint   i, j, k, ind;
  153.   GLdouble tmp[16];
  154.  
  155.   MEMCPY(a,identity,sizeof(double)*16);
  156.   MEMCPY(tmp, b,sizeof(double)*16);
  157.  
  158.   for (i = 0; i != 4; i++) {
  159.  
  160.     val = MAT(tmp,i,i);            /* find pivot */
  161.     ind = i;
  162.     for (j = i + 1; j != 4; j++) {
  163.       if (fabs(MAT(tmp,j,i)) > fabs(val)) {
  164.     ind = j;
  165.     val = MAT(tmp,j,i);
  166.       }
  167.     }
  168.  
  169.     if (ind != i) {            /* swap columns */
  170.       for (j = 0; j != 4; j++) {
  171.     val2 = MAT(a,i,j);
  172.     MAT(a,i,j) = MAT(a,ind,j);
  173.     MAT(a,ind,j) = val2;
  174.     val2 = MAT(tmp,i,j);
  175.     MAT(tmp,i,j) = MAT(tmp,ind,j);
  176.     MAT(tmp,ind,j) = val2;
  177.       }
  178.     }
  179.  
  180.     if (val == 0.0F) {    
  181.       fprintf(stderr,"Singular matrix, no inverse!\n");
  182.       MEMCPY( a, identity, 16*sizeof(GLdouble) );  /* return the identity */
  183.       return;
  184.     }
  185.  
  186.     for (j = 0; j != 4; j++) {
  187.       MAT(tmp,i,j) /= val;
  188.       MAT(a,i,j) /= val;
  189.     }
  190.  
  191.     for (j = 0; j != 4; j++) {        /* eliminate column */
  192.       if (j == i)
  193.     continue;
  194.       val = MAT(tmp,j,i);
  195.       for (k = 0; k != 4; k++) {
  196.     MAT(tmp,j,k) -= MAT(tmp,i,k) * val;
  197.     MAT(a,j,k) -= MAT(a,i,k) * val;
  198.       }
  199.     }
  200.   }
  201. #undef MAT
  202. }
  203.  
  204.  
  205.  
  206.  
  207. /* projection du point (objx,objy,obz) sur l'ecran (winx,winy,winz) */
  208. GLint gluProject(GLdouble objx,GLdouble objy,GLdouble objz,
  209.                  const GLdouble model[16],const GLdouble proj[16],
  210.                  const GLint viewport[4],
  211.                  GLdouble *winx,GLdouble *winy,GLdouble *winz)
  212. {
  213.     /* matrice de transformation */
  214.     GLdouble in[4],out[4];
  215.  
  216.     /* initilise la matrice et le vecteur a transformer */
  217.     in[0]=objx; in[1]=objy; in[2]=objz; in[3]=1.0;
  218.     transform_point(out,model,in);
  219.     transform_point(in,proj,out);
  220.  
  221.     /* d'ou le resultat normalise entre -1 et 1*/
  222.     in[0]/=in[3];in[1]/=in[3];in[2]/=in[3];
  223.  
  224.     /* en coordonnees ecran */
  225.     *winx = viewport[0]+(1+in[0])*viewport[2]/2;
  226.     *winy = viewport[1]+(1+in[1])*viewport[3]/2;
  227.     /* entre 0 et 1 suivant z */
  228.     *winz = (1+in[2])/2;
  229.     return GL_TRUE;
  230. }
  231.  
  232.  
  233.  
  234. /* transformation du point ecran (winx,winy,winz) en point objet */
  235. GLint gluUnProject(GLdouble winx,GLdouble winy,GLdouble winz,
  236.                    const GLdouble model[16],const GLdouble proj[16],
  237.                    const GLint viewport[4],
  238.                    GLdouble *objx,GLdouble *objy,GLdouble *objz)
  239. {
  240.     /* matrice de transformation */
  241.     GLdouble m[16], A[16];
  242.     GLdouble in[4],out[4];
  243.  
  244.     /* transformation coordonnees normalisees entre -1 et 1 */
  245.     in[0]=(winx-viewport[0])*2/viewport[2] - 1.0;
  246.     in[1]=(winy-viewport[1])*2/viewport[3] - 1.0;
  247.     in[2]=2*winz - 1.0;
  248.     in[3]=1.0;
  249.  
  250.     /* calcul transformation inverse */
  251.     matmul(A,proj,model); invert_matrix(A,m);
  252.  
  253.     /* d'ou les coordonnees objets */
  254.     transform_point(out,m,in);
  255.     *objx=out[0]/out[3];
  256.     *objy=out[1]/out[3];
  257.     *objz=out[2]/out[3];
  258.     return GL_TRUE;
  259. }
  260.  
  261.