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