source: trunk/src/testing/app/volvis/MatrixMath.h @ 4

Revision 4, 17.2 KB checked in by ajaworski, 13 years ago (diff)

Added modified SAGE sources
Line 
1/********************************************************************************
2 * Volatile - Volume Visualization Software for SAGE
3 * Copyright (C) 2004 Electronic Visualization Laboratory,
4 * University of Illinois at Chicago
5 *
6 * All rights reserved.
7 *
8 * Redistribution and use in source and binary forms, with or without
9 * modification, are permitted provided that the following conditions are met:
10 *
11 *  * Redistributions of source code must retain the above copyright
12 *    notice, this list of conditions and the following disclaimer.
13 *  * Redistributions in binary form must reproduce the above
14 *    copyright notice, this list of conditions and the following disclaimer
15 *    in the documentation and/or other materials provided with the distribution.
16 *  * Neither the name of the University of Illinois at Chicago nor
17 *    the names of its contributors may be used to endorse or promote
18 *    products derived from this software without specific prior written permission.
19 *
20 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
23 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
24 * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
25 * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
26 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
27 * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
28 * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
29 * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
30 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31 *
32 * Direct questions, comments etc about Volatile to www.evl.uic.edu/cavern/forum
33 *********************************************************************************/
34
35#ifndef __MatrixMath_defined
36#define __MatrixMath_defined
37#include "VectorMath.h"
38
39//BAD SHALINI...dont assume Opengl - do this by hand
40inline void computeMatrix(float angle, float axis[],GLdouble* mat) {
41        // Have OpenGL compute the new transformation (simple but slow)
42        glPushMatrix();
43        glLoadIdentity();
44        glRotatef(angle, axis[0], axis[1], axis[2]);
45        glMultMatrixd((GLdouble *) mat);
46        glGetDoublev(GL_MODELVIEW_MATRIX, (GLdouble *) mat);
47        glPopMatrix();
48}
49
50inline void computeMatrix(float angle, float axis[],GLfloat* mat) {
51        // Have OpenGL compute the new transformation (simple but slow)
52        glPushMatrix();
53        glLoadIdentity();
54        glRotatef(angle, axis[0], axis[1], axis[2]);
55        glMultMatrixf((GLfloat *) mat);
56        glGetFloatv(GL_MODELVIEW_MATRIX, (GLfloat *) mat);
57        glPopMatrix();
58}
59
60
61// copy two matricies -----------------------------------
62inline void matrixEqual(GLfloat me[16], GLfloat m[16])
63{
64        for(int i=0; i< 16; ++i) me[i] = m[i];
65}
66
67inline void matrixEqual(GLdouble me[16], GLdouble m[16])
68{
69        for(int i=0; i< 16; ++i) me[i] = m[i];
70}
71
72// maxb = ma * mb --------------------------------------
73inline void matrixMult( GLfloat maxb[16], GLfloat ma[16], GLfloat mb[16] )
74{
75        maxb[0] = ma[0]*mb[0] + ma[4]*mb[1] + ma[8]*mb[2] + ma[12]*mb[3];
76        maxb[1] = ma[1]*mb[0] + ma[5]*mb[1] + ma[9]*mb[2] + ma[13]*mb[3];
77        maxb[2] = ma[2]*mb[0] + ma[6]*mb[1] + ma[10]*mb[2] + ma[14]*mb[3];
78        maxb[3] = ma[3]*mb[0] + ma[7]*mb[1] + ma[11]*mb[2] + ma[15]*mb[3];
79
80        maxb[4] = ma[0]*mb[4] + ma[4]*mb[5] + ma[8]*mb[6] + ma[12]*mb[7];
81        maxb[5] = ma[1]*mb[4] + ma[5]*mb[5] + ma[9]*mb[6] + ma[13]*mb[7];
82        maxb[6] = ma[2]*mb[4] + ma[6]*mb[5] + ma[10]*mb[6] + ma[14]*mb[7];
83        maxb[7] = ma[3]*mb[4] + ma[7]*mb[5] + ma[11]*mb[6] + ma[15]*mb[7];
84
85        maxb[8] = ma[0]*mb[8] + ma[4]*mb[9] + ma[8]*mb[10] + ma[12]*mb[11];
86        maxb[9] = ma[1]*mb[8] + ma[5]*mb[9] + ma[9]*mb[10] + ma[13]*mb[11];
87        maxb[10] = ma[2]*mb[8] + ma[6]*mb[9] + ma[10]*mb[10] + ma[14]*mb[11];
88        maxb[11] = ma[3]*mb[8] + ma[7]*mb[9] + ma[11]*mb[10] + ma[15]*mb[11];
89
90        maxb[12] = ma[0]*mb[12] + ma[4]*mb[13] + ma[8]*mb[14] + ma[12]*mb[15];
91        maxb[13] = ma[1]*mb[12] + ma[5]*mb[13] + ma[9]*mb[14] + ma[13]*mb[15];
92        maxb[14] = ma[2]*mb[12] + ma[6]*mb[13] + ma[10]*mb[14] + ma[14]*mb[15];
93        maxb[15] = ma[3]*mb[12] + ma[7]*mb[13] + ma[11]*mb[14] + ma[15]*mb[15];
94}
95
96inline void matrixMult( GLdouble maxb[16], GLdouble ma[16], GLdouble mb[16] )
97{
98        maxb[0] = ma[0]*mb[0] + ma[4]*mb[1] + ma[8]*mb[2] + ma[12]*mb[3];
99        maxb[1] = ma[1]*mb[0] + ma[5]*mb[1] + ma[9]*mb[2] + ma[13]*mb[3];
100        maxb[2] = ma[2]*mb[0] + ma[6]*mb[1] + ma[10]*mb[2] + ma[14]*mb[3];
101        maxb[3] = ma[3]*mb[0] + ma[7]*mb[1] + ma[11]*mb[2] + ma[15]*mb[3];
102
103        maxb[4] = ma[0]*mb[4] + ma[4]*mb[5] + ma[8]*mb[6] + ma[12]*mb[7];
104        maxb[5] = ma[1]*mb[4] + ma[5]*mb[5] + ma[9]*mb[6] + ma[13]*mb[7];
105        maxb[6] = ma[2]*mb[4] + ma[6]*mb[5] + ma[10]*mb[6] + ma[14]*mb[7];
106        maxb[7] = ma[3]*mb[4] + ma[7]*mb[5] + ma[11]*mb[6] + ma[15]*mb[7];
107
108        maxb[8] = ma[0]*mb[8] + ma[4]*mb[9] + ma[8]*mb[10] + ma[12]*mb[11];
109        maxb[9] = ma[1]*mb[8] + ma[5]*mb[9] + ma[9]*mb[10] + ma[13]*mb[11];
110        maxb[10] = ma[2]*mb[8] + ma[6]*mb[9] + ma[10]*mb[10] + ma[14]*mb[11];
111        maxb[11] = ma[3]*mb[8] + ma[7]*mb[9] + ma[11]*mb[10] + ma[15]*mb[11];
112
113        maxb[12] = ma[0]*mb[12] + ma[4]*mb[13] + ma[8]*mb[14] + ma[12]*mb[15];
114        maxb[13] = ma[1]*mb[12] + ma[5]*mb[13] + ma[9]*mb[14] + ma[13]*mb[15];
115        maxb[14] = ma[2]*mb[12] + ma[6]*mb[13] + ma[10]*mb[14] + ma[14]*mb[15];
116        maxb[15] = ma[3]*mb[12] + ma[7]*mb[13] + ma[11]*mb[14] + ma[15]*mb[15];
117}
118
119// compute the inverse of a matrix
120// invm = m^(-1)
121inline void inverseMatrix( GLfloat invm[16], GLfloat m[16] )
122{
123        GLfloat det =
124                m[0]*m[5]*m[10]-
125                m[0]*m[6]*m[9]-
126                m[1]*m[4]*m[10]+
127                m[1]*m[6]*m[8]+
128                m[2]*m[4]*m[9]-
129                m[2]*m[5]*m[8];
130
131        invm[0] = (m[5]*m[10]-m[6]*m[9])/det;
132        invm[1] = (-m[1]*m[10]+m[2]*m[9])/det;
133        invm[2] = (m[1]*m[6]-m[2]*m[5])/det;
134        invm[3] = 0.0;
135
136        invm[4] = (-m[4]*m[10]+m[6]*m[8])/det;
137        invm[5] = (m[0]*m[10]-m[2]*m[8])/det;
138        invm[6] = (-m[0]*m[6]+m[2]*m[4])/det;
139        invm[7] = 0.0;
140
141        invm[8] = (m[4]*m[9]-m[5]*m[8])/det;
142        invm[9] = (-m[0]*m[9]+m[1]*m[8])/det;
143        invm[10] = (m[0]*m[5]-m[1]*m[4])/det;
144        invm[11] = 0.0;
145
146        invm[12] = (-m[4]*m[9]*m[14]+m[4]*m[13]*m[10]+
147                m[5]*m[8]*m[14]-m[5]*m[12]*m[10]-
148                m[6]*m[8]*m[13]+m[6]*m[12]*m[9])/det;
149        invm[13] = (m[0]*m[9]*m[14]-m[0]*m[13]*m[10]-
150                m[1]*m[8]*m[14]+m[1]*m[12]*m[10]+
151                m[2]*m[8]*m[13]-m[2]*m[12]*m[9])/det;
152        invm[14] = (-m[0]*m[5]*m[14]+m[0]*m[13]*m[6]+
153                m[1]*m[4]*m[14]-m[1]*m[12]*m[6]-
154                m[2]*m[4]*m[13]+m[2]*m[12]*m[5])/det;
155        invm[15] = 1.0;
156}
157
158inline void inverseMatrix( GLdouble invm[16], GLdouble m[16] )
159{
160        GLfloat det = (float)(
161                m[0]*m[5]*m[10]-
162                m[0]*m[6]*m[9]-
163                m[1]*m[4]*m[10]+
164                m[1]*m[6]*m[8]+
165                m[2]*m[4]*m[9]-
166                m[2]*m[5]*m[8]);
167
168        invm[0] = (m[5]*m[10]-m[6]*m[9])/det;
169        invm[1] = (-m[1]*m[10]+m[2]*m[9])/det;
170        invm[2] = (m[1]*m[6]-m[2]*m[5])/det;
171        invm[3] = 0.0;
172
173        invm[4] = (-m[4]*m[10]+m[6]*m[8])/det;
174        invm[5] = (m[0]*m[10]-m[2]*m[8])/det;
175        invm[6] = (-m[0]*m[6]+m[2]*m[4])/det;
176        invm[7] = 0.0;
177
178        invm[8] = (m[4]*m[9]-m[5]*m[8])/det;
179        invm[9] = (-m[0]*m[9]+m[1]*m[8])/det;
180        invm[10] = (m[0]*m[5]-m[1]*m[4])/det;
181        invm[11] = 0.0;
182
183        invm[12] = (-m[4]*m[9]*m[14]+m[4]*m[13]*m[10]+
184                m[5]*m[8]*m[14]-m[5]*m[12]*m[10]-
185                m[6]*m[8]*m[13]+m[6]*m[12]*m[9])/det;
186        invm[13] = (m[0]*m[9]*m[14]-m[0]*m[13]*m[10]-
187                m[1]*m[8]*m[14]+m[1]*m[12]*m[10]+
188                m[2]*m[8]*m[13]-m[2]*m[12]*m[9])/det;
189        invm[14] = (-m[0]*m[5]*m[14]+m[0]*m[13]*m[6]+
190                m[1]*m[4]*m[14]-m[1]*m[12]*m[6]-
191                m[2]*m[4]*m[13]+m[2]*m[12]*m[5])/det;
192        invm[15] = 1.0;
193}
194
195//transpose a matrix
196inline void transposeMatrix(GLfloat m[16])
197{
198        GLfloat tmp;
199        tmp = m[1];
200        m[1] = m[4];
201        m[4] = tmp;
202        tmp = m[2];
203        m[2] = m[8];
204        m[8] = tmp;
205        tmp = m[3];
206        m[3] = m[12];
207        m[12] = tmp;
208        tmp = m[6];
209        m[6] = m[9];
210        m[9] = tmp;
211        tmp = m[7];
212        m[7] = m[13];
213        m[13] = tmp;
214        tmp = m[11];
215        m[11] = m[14];
216        m[14] = tmp;
217
218}
219
220inline void transposeMatrix(GLdouble m[16])
221{
222        GLdouble tmp;
223        tmp = m[1];
224        m[1] = m[4];
225        m[4] = tmp;
226        tmp = m[2];
227        m[2] = m[8];
228        m[8] = tmp;
229        tmp = m[3];
230        m[3] = m[12];
231        m[12] = tmp;
232        tmp = m[6];
233        m[6] = m[9];
234        m[9] = tmp;
235        tmp = m[7];
236        m[7] = m[13];
237        m[13] = tmp;
238        tmp = m[11];
239        m[11] = m[14];
240        m[14] = tmp;
241
242}
243
244//angle and axis to a rotation matrix ----------------------------
245inline void axis2Rot( GLfloat m[16], GLfloat k[3], GLfloat theta )
246{
247        float c = (float)cos(theta);
248        float s = (float)sin(theta);
249        float v = 1 - c;
250
251        m[0] = k[0]*k[0]*v + c;
252        m[4] = k[0]*k[1]*v - k[2]*s;
253        m[8] = k[0]*k[2]*v + k[1]*s;
254
255        m[1] = k[0]*k[1]*v + k[2]*s;
256        m[5] = k[1]*k[1]*v + c;
257        m[9] = k[1]*k[2]*v - k[0]*s;
258
259        m[2] = k[0]*k[2]*v - k[1]*s;
260        m[6] = k[1]*k[2]*v + k[0]*s;
261        m[10] = k[2]*k[2]*v + c;
262}
263
264//angle and axis to a rotation matrix ----------------------------
265inline void axis2Rot( GLdouble m[16], GLfloat k[3], GLfloat theta )
266{
267        float c = (float)cos(theta);
268        float s = (float)sin(theta);
269        float v = 1 - c;
270
271        m[0] = k[0]*k[0]*v + c;
272        m[4] = k[0]*k[1]*v - k[2]*s;
273        m[8] = k[0]*k[2]*v + k[1]*s;
274
275        m[1] = k[0]*k[1]*v + k[2]*s;
276        m[5] = k[1]*k[1]*v + c;
277        m[9] = k[1]*k[2]*v - k[0]*s;
278
279        m[2] = k[0]*k[2]*v - k[1]*s;
280        m[6] = k[1]*k[2]*v + k[0]*s;
281        m[10] = k[2]*k[2]*v + c;
282}
283
284
285inline void axis2Rot( GLdouble m[16], GLdouble k[3], GLdouble theta )
286{
287        float c = (float)cos(theta);
288        float s = (float)sin(theta);
289        float v = 1 - c;
290
291        m[0] = k[0]*k[0]*v + c;
292        m[4] = k[0]*k[1]*v - k[2]*s;
293        m[8] = k[0]*k[2]*v + k[1]*s;
294
295        m[1] = k[0]*k[1]*v + k[2]*s;
296        m[5] = k[1]*k[1]*v + c;
297        m[9] = k[1]*k[2]*v - k[0]*s;
298
299        m[2] = k[0]*k[2]*v - k[1]*s;
300        m[6] = k[1]*k[2]*v + k[0]*s;
301        m[10] = k[2]*k[2]*v + c;
302}
303
304// Inverse angle-axis formula.-----------------------------------
305
306inline void rot2Axis( GLfloat k[3], GLfloat *theta, GLfloat m[16] )
307{
308        GLfloat c = (float)(0.5 * (m[0] + m[5] + m[10] - 1.0));
309        GLfloat r1 = m[6] - m[9];
310        GLfloat r2 = m[8] - m[2];
311        GLfloat r3 = m[1] - m[4];
312        GLfloat s = (float)(0.5 * sqrt(r1*r1+r2*r2+r3*r3));
313
314        *theta = (float)atan2(s, c);
315
316        if( fabs(s) > EPSILON )
317    {
318                c = (float)(2.0*s);
319
320                k[0] = r1 / c;
321                k[1] = r2 / c;
322                k[2] = r3 / c;
323    }
324        else
325    {
326                if( c > 0 ) // theta = 0
327                {
328                        k[0] = 0;
329                        k[1] = 0;
330                        k[2] = 1;
331                }
332                else // theta = +-pi: Shepperd procedure
333                {
334                        GLfloat k0_2 = (m[0] + 1)/2;
335                        GLfloat k1_2 = (m[5] + 1)/2;
336                        GLfloat k2_2 = (m[10] + 1)/2;
337
338                        if( k0_2 > k1_2 )
339                        {
340                                if( k0_2 > k2_2 ) // k0 biggest
341                                {
342                                        k[0] = (float)sqrt(k1_2);
343                                        k[1] = (m[1] + m[4])/(4*k[0]);
344                                        k[2] = (m[2] + m[8])/(4*k[0]);
345                                }
346                                else // k2 biggest
347                                {
348                                        k[2] = (float)sqrt(k2_2);
349                                        k[0] = (m[2] + m[8])/(4*k[2]);
350                                        k[1] = (m[6] + m[9])/(4*k[2]);
351                                }
352                        }
353                        else
354                        {
355                                if( k1_2 > k2_2 ) // k1 biggest
356                                {
357                                        k[1] = (float)sqrt(k1_2);
358                                        k[0] = (m[1] + m[4])/(4*k[1]);
359                                        k[2] = (m[6] + m[9])/(4*k[1]);
360                                }
361                                else // k2 biggest
362                                {
363                                        k[2] = (float)sqrt(k2_2);
364                                        k[0] = (m[2] + m[8])/(4*k[2]);
365                                        k[1] = (m[6] + m[9])/(4*k[2]);
366                                }
367                        }
368                }
369    }
370}
371
372//quaternion to rotation matrix
373
374inline void q2R( GLfloat m[16], GLfloat q[4] )
375{
376        m[0] = q[0]*q[0] + q[1]*q[1] - q[2]*q[2] - q[3]*q[3];
377        m[1] = 2*q[1]*q[2] + 2*q[0]*q[3];
378        m[2] = 2*q[1]*q[3] - 2*q[0]*q[2];
379
380        m[4] = 2*q[1]*q[2] - 2*q[0]*q[3];
381        m[5] = q[0]*q[0] + q[2]*q[2] - q[1]*q[1] - q[3]*q[3];
382        m[6] = 2*q[2]*q[3] + 2*q[0]*q[1];
383
384        m[8] = 2*q[1]*q[3] + 2*q[0]*q[2];
385        m[9] = 2*q[2]*q[3] - 2*q[0]*q[1];
386        m[10]= q[0]*q[0] + q[3]*q[3] - q[1]*q[1] - q[2]*q[2];
387}
388
389inline void axisToQuat(float a[3], float phi, float q[4])
390{
391  vNormal(a);
392  vCopy(a,q);
393  vScale(q,(float)(sin(phi/2.0)));
394  q[3] = (float)(cos(phi/2.0));
395}
396
397/////////////////////////////////////////////////////////////////////////////
398//
399//  Quaternions always obey:  a^2 + b^2 + c^2 + d^2 = 1.0
400//  If they don't add up to 1.0, dividing by their magnitued will
401//  renormalize them.
402//
403//  Note: See the following for more information on quaternions:
404//
405//  - Shoemake, K., Animating rotation with quaternion curves, Computer
406//    Graphics 19, No 3 (Proc. SIGGRAPH'85), 245-254, 1985.
407//  - Pletinckx, D., Quaternion calculus as a basic tool in computer
408//    graphics, The Visual Computer 5, 2-13, 1989.
409//
410inline void normalizeQuat(float q[4])
411{
412    float  mag = (q[0]*q[0] + q[1]*q[1] + q[2]*q[2] + q[3]*q[3]);
413    for (int i = 0; i < 4; i++) q[i] /= mag;
414}
415
416/////////////////////////////////////////////////////////////////////////////
417//
418//  Given two rotations, e1 and e2, expressed as quaternion rotations,
419//  figure out the equivalent single rotation and stuff it into dest.
420//
421//  This routine also normalizes the result every RENORMCOUNT times it is
422//  called, to keep error from creeping in.
423//
424//  NOTE: This routine is written so that q1 or q2 may be the same
425//  as dest (or each other).
426//
427inline void addQuats(float q1[4], float q2[4], float dest[4])
428{
429  static int count=0;
430  float t1[4], t2[4], t3[4];
431  float tf[4];
432
433  vCopy(q1,t1);
434  vScale(t1,q2[3]);
435
436  vCopy(q2,t2);
437  vScale(t2,q1[3]);
438
439  vCross(q2,q1,t3);
440  vAdd(t1,t2,tf);
441  vAdd(t3,tf,tf);
442  tf[3] = q1[3] * q2[3] - vDot(q1,q2);
443
444  dest[0] = tf[0];
445  dest[1] = tf[1];
446  dest[2] = tf[2];
447  dest[3] = tf[3];
448
449  if (++count > 97) {
450    count = 0;
451    normalizeQuat(dest);
452  }
453}
454
455
456
457/*
458* Prints matrix to stderr.
459*/
460
461inline void printMatrix( GLfloat m[16] )
462{
463        int i, j;
464
465        for( i=0; i<4; i++ )
466    {
467                for( j=0; j<4; j++ )
468                        fprintf(stderr, "%f ", m[i+j*4]);
469                fprintf(stderr, "\n");
470    }
471        fprintf(stderr, "\n");
472}
473
474inline void printMatrix( GLdouble m[16] )
475{
476        int i, j;
477
478        for( i=0; i<4; i++ )
479    {
480                for( j=0; j<4; j++ )
481                        fprintf(stderr, "%f ", m[i+j*4]);
482                fprintf(stderr, "\n");
483    }
484        fprintf(stderr, "\n");
485}
486
487inline void buildLookAt(GLfloat m[16],
488                                                GLfloat eye[3], GLfloat at[3], GLfloat up[3])
489{
490        //I am not 100% certain that my cross products are in the right order
491
492        GLfloat tmp1[3], tmp2[3], tmp3[3];
493
494        subV3(tmp1, eye, at);
495        GLfloat norm;
496        norm = normV3(tmp1);
497        tmp1[0] /=norm;
498        tmp1[1] /=norm;
499        tmp1[2] /=norm;
500
501        m[2]     =tmp1[0];
502        m[6]     =tmp1[1];
503        m[10]    =tmp1[2];
504
505        crossV3(tmp2, up, tmp1);
506        norm = normV3(tmp2);
507        tmp2[0] /=norm;
508        tmp2[1] /=norm;
509        tmp2[2] /=norm;
510
511        m[0]     =tmp2[0];
512        m[4]     =tmp2[1];
513        m[8]     =tmp2[2];
514
515        crossV3(tmp3, tmp1, tmp2);
516        norm = normV3(tmp3);
517        tmp3[0] /=norm;
518        tmp3[1] /=norm;
519        tmp3[2] /=norm;
520
521        m[1]     =tmp3[0];
522        m[5]     =tmp3[1];
523        m[9]     =tmp3[2];
524
525        m[12]= -eye[0];
526        m[13]= -eye[1];
527        m[14]= -eye[2];
528        m[15]= 1;
529
530        m[3]  = 0;
531        m[7]  = 0;
532        m[11] = 0;
533}
534
535//give the original 8 vertices and the model view matrix
536//This function returns true if the given point is within the bounds of these vertices
537inline int withinBounds3D(GLdouble mv[16], float bbox[8][3], float point[3]) {
538        float min[3], max[3];
539        min[0] = min[1] = min[2] = 10;
540        max[0] = max[1] = max[2] = -10;
541        float rotatedbbox[8][3];
542        for(int i=0; i<8; ++i){
543                translateV3W(rotatedbbox[i], mv, bbox[i]); //get the rotated vol coords
544                //now get the max and min z in view space
545                if(max[2]< MAX(max[2], rotatedbbox[i][2])){
546                        max[2] = MAX(max[2], rotatedbbox[i][2]);
547                }
548                if(min[2] > MIN(min[2], rotatedbbox[i][2])){
549                        min[2] = MIN(min[2], rotatedbbox[i][2]);
550                }
551                //now get the max and min x in view space
552                if(max[0]< MAX(max[0], rotatedbbox[i][0])){
553                        max[0] = MAX(max[0], rotatedbbox[i][0]);
554                }
555                if(min[0] > MIN(min[0], rotatedbbox[i][0])){
556                        min[0] = MIN(min[0], rotatedbbox[i][0]);
557                }
558                //now get the max and min y in view space
559                if(max[1]< MAX(max[1], rotatedbbox[i][1])){
560                        max[1] = MAX(max[1], rotatedbbox[i][1]);
561                }
562                if(min[1] > MIN(min[1], rotatedbbox[i][1])){
563                        min[1] = MIN(min[1], rotatedbbox[i][1]);
564                }
565        }
566        if (point[0]>min[0] && point[0]<max[0]
567                && point[1]>min[1] && point[1]<max[1]) {
568         //     && point[2]>min[2] && point[2]<max[2]) {
569               //pointer is within the bounds
570                return true;
571        }
572        else
573                return false;
574}
575
576//give the original 4 vertices and the model view matrix
577//This function returns true if the given point is within the bounds of these vertices
578inline int withinBounds2D(GLdouble mv[16], float bbox[4][3], float point[3]) {
579        float min[3], max[3];
580        min[0] = min[1] = min[2] = 10;
581        max[0] = max[1] = max[2] = -10;
582        float rotatedbbox[4][3];
583        for(int i=0; i<4; ++i){
584                translateV3W(rotatedbbox[i], mv, bbox[i]); //get the rotated vol coords
585                //now get the max and min z in view space
586                if(max[2]< MAX(max[2], rotatedbbox[i][2])){
587                        max[2] = MAX(max[2], rotatedbbox[i][2]);
588                }
589                if(min[2] > MIN(min[2], rotatedbbox[i][2])){
590                        min[2] = MIN(min[2], rotatedbbox[i][2]);
591                }
592                //now get the max and min x in view space
593                if(max[0]< MAX(max[0], rotatedbbox[i][0])){
594                        max[0] = MAX(max[0], rotatedbbox[i][0]);
595                }
596                if(min[0] > MIN(min[0], rotatedbbox[i][0])){
597                        min[0] = MIN(min[0], rotatedbbox[i][0]);
598                }
599                //now get the max and min y in view space
600                if(max[1]< MAX(max[1], rotatedbbox[i][1])){
601                        max[1] = MAX(max[1], rotatedbbox[i][1]);
602                }
603                if(min[1] > MIN(min[1], rotatedbbox[i][1])){
604                        min[1] = MIN(min[1], rotatedbbox[i][1]);
605                }
606        }
607        if (point[0]>min[0] && point[0]<max[0]
608                && point[1]>min[1] && point[1]<max[1]) {
609//              && point[2]>min[2] && point[2]<max[2]) {
610               //pointer is within the bounds
611                return true;
612        }
613        else
614                return false;
615}
616
617
618#endif
Note: See TracBrowser for help on using the repository browser.