~diresu/blender/blender-command-port

« back to all changes in this revision

Viewing changes to intern/elbeem/intern/ntl_vector3dim.h

  • Committer: theeth
  • Date: 2008-10-14 16:52:04 UTC
  • Revision ID: vcs-imports@canonical.com-20081014165204-r32w2gm6s0osvdhn
copy back trunk

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/******************************************************************************
 
2
 *
 
3
 * El'Beem - Free Surface Fluid Simulation with the Lattice Boltzmann Method
 
4
 * Copyright 2003-2006 Nils Thuerey
 
5
 *
 
6
 * Basic vector class used everywhere, either blitz or inlined GRAPA class
 
7
 *
 
8
 *****************************************************************************/
 
9
#ifndef NTL_VECTOR3DIM_H
 
10
#define NTL_VECTOR3DIM_H
 
11
 
 
12
// this serves as the main include file
 
13
// for all kinds of stuff that might be required
 
14
// under windos there seem to be strange 
 
15
// errors when including the STL header too
 
16
// late...
 
17
#include <iostream>
 
18
#include <map>
 
19
#include <vector>
 
20
#include <string>
 
21
#include <sstream>
 
22
#include <math.h>
 
23
#include <string.h>
 
24
#include <stdio.h>
 
25
#include <stdlib.h>
 
26
 
 
27
// hack for MSVC6.0 compiler
 
28
#ifdef _MSC_VER
 
29
#if _MSC_VER < 1300
 
30
#define for     if(false); else for
 
31
#define map     std::map
 
32
#define vector  std::vector
 
33
#define string  std::string
 
34
// use this define for MSVC6 stuff hereafter
 
35
#define USE_MSVC6FIXES
 
36
#else // _MSC_VER < 1300 , 7.0 or higher
 
37
using std::map;
 
38
using std::vector;
 
39
using std::string;
 
40
#endif
 
41
#else // not MSVC6
 
42
// for proper compilers...
 
43
using std::map;
 
44
using std::vector;
 
45
using std::string;
 
46
#endif // MSVC6
 
47
 
 
48
#ifdef __APPLE_CC__
 
49
// apple
 
50
#else
 
51
#ifdef WIN32
 
52
 
 
53
// windows values missing, see below
 
54
#ifndef snprintf
 
55
#define snprintf _snprintf
 
56
#endif
 
57
#ifndef bool
 
58
#define bool int
 
59
#endif
 
60
#ifndef false
 
61
#define false 0
 
62
#endif
 
63
#ifndef true
 
64
#define true 1
 
65
#endif
 
66
 
 
67
#else // WIN32
 
68
 
 
69
// floating point limits for linux,*bsd etc...
 
70
#include <float.h>
 
71
 
 
72
#endif // WIN32
 
73
#endif // __APPLE_CC__
 
74
 
 
75
// windos, hardcoded limits for now...
 
76
// for e.g. MSVC compiler...
 
77
// some of these defines can be needed
 
78
// for linux systems as well (e.g. FLT_MAX)
 
79
#ifndef __FLT_MAX__
 
80
#       ifdef FLT_MAX  // try to use it instead
 
81
#               define __FLT_MAX__ FLT_MAX
 
82
#       else // FLT_MAX
 
83
#               define __FLT_MAX__ 3.402823466e+38f
 
84
#       endif // FLT_MAX
 
85
#endif // __FLT_MAX__
 
86
#ifndef __DBL_MAX__
 
87
#       ifdef DBL_MAX // try to use it instead
 
88
#               define __DBL_MAX__ DBL_MAX
 
89
#       else // DBL_MAX
 
90
#               define __DBL_MAX__ 1.7976931348623158e+308
 
91
#       endif // DBL_MAX
 
92
#endif // __DBL_MAX__
 
93
 
 
94
#ifndef M_PI
 
95
#define M_PI 3.1415926536
 
96
#define M_E  2.7182818284
 
97
#endif
 
98
 
 
99
// make sure elbeem plugin def is valid
 
100
#if ELBEEM_BLENDER==1
 
101
#ifndef ELBEEM_PLUGIN
 
102
#define ELBEEM_PLUGIN 1
 
103
#endif // !ELBEEM_PLUGIN
 
104
#endif // ELBEEM_BLENDER==1
 
105
 
 
106
// make sure GUI support is disabled for plugin use
 
107
#if ELBEEM_PLUGIN==1
 
108
#ifndef NOGUI
 
109
#define NOGUI 1
 
110
#endif // !NOGUI
 
111
#endif // ELBEEM_PLUGIN==1
 
112
 
 
113
 
 
114
// basic inlined vector class
 
115
template<class Scalar>
 
116
class ntlVector3Dim
 
117
{
 
118
public:
 
119
  // Constructor
 
120
  inline ntlVector3Dim(void );
 
121
  // Copy-Constructor
 
122
  inline ntlVector3Dim(const ntlVector3Dim<Scalar> &v );
 
123
  inline ntlVector3Dim(const float *);
 
124
  inline ntlVector3Dim(const double *);
 
125
  // construct a vector from one Scalar
 
126
  inline ntlVector3Dim(Scalar);
 
127
  // construct a vector from three Scalars
 
128
  inline ntlVector3Dim(Scalar, Scalar, Scalar);
 
129
 
 
130
        // get address of array for OpenGL
 
131
        Scalar *getAddress() { return value; }
 
132
 
 
133
  // Assignment operator
 
134
  inline const ntlVector3Dim<Scalar>& operator=  (const ntlVector3Dim<Scalar>& v);
 
135
  // Assignment operator
 
136
  inline const ntlVector3Dim<Scalar>& operator=  (Scalar s);
 
137
  // Assign and add operator
 
138
  inline const ntlVector3Dim<Scalar>& operator+= (const ntlVector3Dim<Scalar>& v);
 
139
  // Assign and add operator
 
140
  inline const ntlVector3Dim<Scalar>& operator+= (Scalar s);
 
141
  // Assign and sub operator
 
142
  inline const ntlVector3Dim<Scalar>& operator-= (const ntlVector3Dim<Scalar>& v);
 
143
  // Assign and sub operator
 
144
  inline const ntlVector3Dim<Scalar>& operator-= (Scalar s);
 
145
  // Assign and mult operator
 
146
  inline const ntlVector3Dim<Scalar>& operator*= (const ntlVector3Dim<Scalar>& v);
 
147
  // Assign and mult operator
 
148
  inline const ntlVector3Dim<Scalar>& operator*= (Scalar s);
 
149
  // Assign and div operator
 
150
  inline const ntlVector3Dim<Scalar>& operator/= (const ntlVector3Dim<Scalar>& v);
 
151
  // Assign and div operator
 
152
  inline const ntlVector3Dim<Scalar>& operator/= (Scalar s);
 
153
 
 
154
 
 
155
  // unary operator
 
156
  inline ntlVector3Dim<Scalar> operator- () const;
 
157
 
 
158
  // binary operator add
 
159
  inline ntlVector3Dim<Scalar> operator+ (const ntlVector3Dim<Scalar>&) const;
 
160
  // binary operator add
 
161
  inline ntlVector3Dim<Scalar> operator+ (Scalar) const;
 
162
  // binary operator sub
 
163
  inline ntlVector3Dim<Scalar> operator- (const ntlVector3Dim<Scalar>&) const;
 
164
  // binary operator sub
 
165
  inline ntlVector3Dim<Scalar> operator- (Scalar) const;
 
166
  // binary operator mult
 
167
  inline ntlVector3Dim<Scalar> operator* (const ntlVector3Dim<Scalar>&) const;
 
168
  // binary operator mult
 
169
  inline ntlVector3Dim<Scalar> operator* (Scalar) const;
 
170
  // binary operator div
 
171
  inline ntlVector3Dim<Scalar> operator/ (const ntlVector3Dim<Scalar>&) const;
 
172
  // binary operator div
 
173
  inline ntlVector3Dim<Scalar> operator/ (Scalar) const;
 
174
 
 
175
  // Projection normal to a vector
 
176
  inline ntlVector3Dim<Scalar>    getOrthogonalntlVector3Dim() const;
 
177
  // Project into a plane
 
178
  inline const ntlVector3Dim<Scalar>& projectNormalTo(const ntlVector3Dim<Scalar> &v);
 
179
  
 
180
  // minimize
 
181
  inline const ntlVector3Dim<Scalar> &minimize(const ntlVector3Dim<Scalar> &);
 
182
  // maximize
 
183
  inline const ntlVector3Dim<Scalar> &maximize(const ntlVector3Dim<Scalar> &);
 
184
  
 
185
  // access operator
 
186
  inline Scalar& operator[](unsigned int i);
 
187
  // access operator
 
188
  inline const Scalar& operator[](unsigned int i) const;
 
189
 
 
190
protected:
 
191
  
 
192
private:
 
193
  Scalar value[3];  //< Storage of vector values
 
194
};
 
195
 
 
196
 
 
197
 
 
198
 
 
199
//------------------------------------------------------------------------------
 
200
// STREAM FUNCTIONS
 
201
//------------------------------------------------------------------------------
 
202
 
 
203
 
 
204
 
 
205
//! global string for formatting vector output in utilities.cpp
 
206
extern const char *globVecFormatStr;
 
207
 
 
208
/*************************************************************************
 
209
  Outputs the object in human readable form using the format
 
210
  [x,y,z]
 
211
  */
 
212
template<class Scalar>
 
213
std::ostream&
 
214
operator<<( std::ostream& os, const ntlVector3Dim<Scalar>& i )
 
215
{
 
216
        char buf[256];
 
217
        snprintf(buf,256,globVecFormatStr,i[0],i[1],i[2]);
 
218
        os << string(buf); 
 
219
  //os << '[' << i[0] << ", " << i[1] << ", " << i[2] << ']';
 
220
  return os;
 
221
}
 
222
 
 
223
 
 
224
 
 
225
/*************************************************************************
 
226
  Reads the contents of the object from a stream using the same format
 
227
  as the output operator.
 
228
  */
 
229
template<class Scalar>
 
230
std::istream&
 
231
operator>>( std::istream& is, ntlVector3Dim<Scalar>& i )
 
232
{
 
233
  char c;
 
234
  char dummy[3];
 
235
  is >> c >> i[0] >> dummy >> i[1] >> dummy >> i[2] >> c;
 
236
  return is;
 
237
}
 
238
 
 
239
 
 
240
//------------------------------------------------------------------------------
 
241
// VECTOR inline FUNCTIONS
 
242
//------------------------------------------------------------------------------
 
243
 
 
244
 
 
245
 
 
246
/*************************************************************************
 
247
  Constructor.
 
248
  */
 
249
template<class Scalar>
 
250
inline ntlVector3Dim<Scalar>::ntlVector3Dim( void )
 
251
{
 
252
  value[0] = value[1] = value[2] = 0;
 
253
}
 
254
 
 
255
 
 
256
 
 
257
/*************************************************************************
 
258
  Copy-Constructor.
 
259
  */
 
260
template<class Scalar>
 
261
inline ntlVector3Dim<Scalar>::ntlVector3Dim( const ntlVector3Dim<Scalar> &v )
 
262
{
 
263
  value[0] = v.value[0];
 
264
  value[1] = v.value[1];
 
265
  value[2] = v.value[2];
 
266
}
 
267
template<class Scalar>
 
268
inline ntlVector3Dim<Scalar>::ntlVector3Dim( const float *fvalue)
 
269
{
 
270
  value[0] = (Scalar)fvalue[0];
 
271
  value[1] = (Scalar)fvalue[1];
 
272
  value[2] = (Scalar)fvalue[2];
 
273
}
 
274
template<class Scalar>
 
275
inline ntlVector3Dim<Scalar>::ntlVector3Dim( const double *fvalue)
 
276
{
 
277
  value[0] = (Scalar)fvalue[0];
 
278
  value[1] = (Scalar)fvalue[1];
 
279
  value[2] = (Scalar)fvalue[2];
 
280
}
 
281
 
 
282
 
 
283
 
 
284
/*************************************************************************
 
285
  Constructor for a vector from a single Scalar. All components of
 
286
  the vector get the same value.
 
287
  \param s The value to set
 
288
  \return The new vector
 
289
  */
 
290
template<class Scalar>
 
291
inline ntlVector3Dim<Scalar>::ntlVector3Dim(Scalar s )
 
292
{
 
293
  value[0]= s;
 
294
  value[1]= s;
 
295
  value[2]= s;
 
296
}
 
297
 
 
298
 
 
299
/*************************************************************************
 
300
  Constructor for a vector from three Scalars.
 
301
  \param s1 The value for the first vector component
 
302
  \param s2 The value for the second vector component
 
303
  \param s3 The value for the third vector component
 
304
  \return The new vector
 
305
  */
 
306
template<class Scalar>
 
307
inline ntlVector3Dim<Scalar>::ntlVector3Dim(Scalar s1, Scalar s2, Scalar s3)
 
308
{
 
309
  value[0]= s1;
 
310
  value[1]= s2;
 
311
  value[2]= s3;
 
312
}
 
313
 
 
314
 
 
315
/*************************************************************************
 
316
  Compute the vector product of two 3D vectors
 
317
  \param v Second vector to compute the product with
 
318
  \return A new vector with the product values
 
319
  */
 
320
/*template<class Scalar>
 
321
inline ntlVector3Dim<Scalar> 
 
322
ntlVector3Dim<Scalar>::operator^( const ntlVector3Dim<Scalar> &v ) const
 
323
{
 
324
  return ntlVector3Dim<Scalar>(value[1]*v.value[2] - value[2]*v.value[1],
 
325
                        value[2]*v.value[0] - value[0]*v.value[2],
 
326
                        value[0]*v.value[1] - value[1]*v.value[0]);
 
327
}*/
 
328
 
 
329
 
 
330
/*************************************************************************
 
331
  Copy a ntlVector3Dim componentwise.
 
332
  \param v vector with values to be copied
 
333
  \return Reference to self
 
334
  */
 
335
template<class Scalar>
 
336
inline const ntlVector3Dim<Scalar>&
 
337
ntlVector3Dim<Scalar>::operator=( const ntlVector3Dim<Scalar> &v )
 
338
{
 
339
  value[0] = v.value[0];
 
340
  value[1] = v.value[1];
 
341
  value[2] = v.value[2];  
 
342
  return *this;
 
343
}
 
344
 
 
345
 
 
346
/*************************************************************************
 
347
  Copy a Scalar to each component.
 
348
  \param s The value to copy
 
349
  \return Reference to self
 
350
  */
 
351
template<class Scalar>
 
352
inline const ntlVector3Dim<Scalar>&
 
353
ntlVector3Dim<Scalar>::operator=(Scalar s)
 
354
{
 
355
  value[0] = s;
 
356
  value[1] = s;
 
357
  value[2] = s;  
 
358
  return *this;
 
359
}
 
360
 
 
361
 
 
362
/*************************************************************************
 
363
  Add another ntlVector3Dim componentwise.
 
364
  \param v vector with values to be added
 
365
  \return Reference to self
 
366
  */
 
367
template<class Scalar>
 
368
inline const ntlVector3Dim<Scalar>&
 
369
ntlVector3Dim<Scalar>::operator+=( const ntlVector3Dim<Scalar> &v )
 
370
{
 
371
  value[0] += v.value[0];
 
372
  value[1] += v.value[1];
 
373
  value[2] += v.value[2];  
 
374
  return *this;
 
375
}
 
376
 
 
377
 
 
378
/*************************************************************************
 
379
  Add a Scalar value to each component.
 
380
  \param s Value to add
 
381
  \return Reference to self
 
382
  */
 
383
template<class Scalar>
 
384
inline const ntlVector3Dim<Scalar>&
 
385
ntlVector3Dim<Scalar>::operator+=(Scalar s)
 
386
{
 
387
  value[0] += s;
 
388
  value[1] += s;
 
389
  value[2] += s;  
 
390
  return *this;
 
391
}
 
392
 
 
393
 
 
394
/*************************************************************************
 
395
  Subtract another vector componentwise.
 
396
  \param v vector of values to subtract
 
397
  \return Reference to self
 
398
  */
 
399
template<class Scalar>
 
400
inline const ntlVector3Dim<Scalar>&
 
401
ntlVector3Dim<Scalar>::operator-=( const ntlVector3Dim<Scalar> &v )
 
402
{
 
403
  value[0] -= v.value[0];
 
404
  value[1] -= v.value[1];
 
405
  value[2] -= v.value[2];  
 
406
  return *this;
 
407
}
 
408
 
 
409
 
 
410
/*************************************************************************
 
411
  Subtract a Scalar value from each component.
 
412
  \param s Value to subtract
 
413
  \return Reference to self
 
414
  */
 
415
template<class Scalar>
 
416
inline const ntlVector3Dim<Scalar>&
 
417
ntlVector3Dim<Scalar>::operator-=(Scalar s)
 
418
{
 
419
  value[0]-= s;
 
420
  value[1]-= s;
 
421
  value[2]-= s;  
 
422
  return *this;
 
423
}
 
424
 
 
425
 
 
426
/*************************************************************************
 
427
  Multiply with another vector componentwise.
 
428
  \param v vector of values to multiply with
 
429
  \return Reference to self
 
430
  */
 
431
template<class Scalar>
 
432
inline const ntlVector3Dim<Scalar>&
 
433
ntlVector3Dim<Scalar>::operator*=( const ntlVector3Dim<Scalar> &v )
 
434
{
 
435
  value[0] *= v.value[0];
 
436
  value[1] *= v.value[1];
 
437
  value[2] *= v.value[2];  
 
438
  return *this;
 
439
}
 
440
 
 
441
 
 
442
/*************************************************************************
 
443
  Multiply each component with a Scalar value.
 
444
  \param s Value to multiply with
 
445
  \return Reference to self
 
446
  */
 
447
template<class Scalar>
 
448
inline const ntlVector3Dim<Scalar>&
 
449
ntlVector3Dim<Scalar>::operator*=(Scalar s)
 
450
{
 
451
  value[0] *= s;
 
452
  value[1] *= s;
 
453
  value[2] *= s;  
 
454
  return *this;
 
455
}
 
456
 
 
457
 
 
458
/*************************************************************************
 
459
  Divide by another ntlVector3Dim componentwise.
 
460
  \param v vector of values to divide by
 
461
  \return Reference to self
 
462
  */
 
463
template<class Scalar>
 
464
inline const ntlVector3Dim<Scalar>&
 
465
ntlVector3Dim<Scalar>::operator/=( const ntlVector3Dim<Scalar> &v )
 
466
{
 
467
  value[0] /= v.value[0];
 
468
  value[1] /= v.value[1];
 
469
  value[2] /= v.value[2];  
 
470
  return *this;
 
471
}
 
472
 
 
473
 
 
474
/*************************************************************************
 
475
  Divide each component by a Scalar value.
 
476
  \param s Value to divide by
 
477
  \return Reference to self
 
478
  */
 
479
template<class Scalar>
 
480
inline const ntlVector3Dim<Scalar>&
 
481
ntlVector3Dim<Scalar>::operator/=(Scalar s)
 
482
{
 
483
  value[0] /= s;
 
484
  value[1] /= s;
 
485
  value[2] /= s;
 
486
  return *this;
 
487
}
 
488
 
 
489
 
 
490
//------------------------------------------------------------------------------
 
491
// unary operators
 
492
//------------------------------------------------------------------------------
 
493
 
 
494
 
 
495
/*************************************************************************
 
496
  Build componentwise the negative this vector.
 
497
  \return The new (negative) vector
 
498
  */
 
499
template<class Scalar>
 
500
inline ntlVector3Dim<Scalar>
 
501
ntlVector3Dim<Scalar>::operator-() const
 
502
{
 
503
  return ntlVector3Dim<Scalar>(-value[0], -value[1], -value[2]);
 
504
}
 
505
 
 
506
 
 
507
 
 
508
//------------------------------------------------------------------------------
 
509
// binary operators
 
510
//------------------------------------------------------------------------------
 
511
 
 
512
 
 
513
/*************************************************************************
 
514
  Build a vector with another vector added componentwise.
 
515
  \param v The second vector to add
 
516
  \return The sum vector
 
517
  */
 
518
template<class Scalar>
 
519
inline ntlVector3Dim<Scalar>
 
520
ntlVector3Dim<Scalar>::operator+( const ntlVector3Dim<Scalar> &v ) const
 
521
{
 
522
  return ntlVector3Dim<Scalar>(value[0]+v.value[0],
 
523
                        value[1]+v.value[1],
 
524
                        value[2]+v.value[2]);
 
525
}
 
526
 
 
527
 
 
528
/*************************************************************************
 
529
  Build a vector with a Scalar value added to each component.
 
530
  \param s The Scalar value to add
 
531
  \return The sum vector
 
532
  */
 
533
template<class Scalar>
 
534
inline ntlVector3Dim<Scalar>
 
535
ntlVector3Dim<Scalar>::operator+(Scalar s) const
 
536
{
 
537
  return ntlVector3Dim<Scalar>(value[0]+s,
 
538
                        value[1]+s,
 
539
                        value[2]+s);
 
540
}
 
541
 
 
542
 
 
543
/*************************************************************************
 
544
  Build a vector with another vector subtracted componentwise.
 
545
  \param v The second vector to subtract
 
546
  \return The difference vector
 
547
  */
 
548
template<class Scalar>
 
549
inline ntlVector3Dim<Scalar>
 
550
ntlVector3Dim<Scalar>::operator-( const ntlVector3Dim<Scalar> &v ) const
 
551
{
 
552
  return ntlVector3Dim<Scalar>(value[0]-v.value[0],
 
553
                        value[1]-v.value[1],
 
554
                        value[2]-v.value[2]);
 
555
}
 
556
 
 
557
 
 
558
/*************************************************************************
 
559
  Build a vector with a Scalar value subtracted componentwise.
 
560
  \param s The Scalar value to subtract
 
561
  \return The difference vector
 
562
  */
 
563
template<class Scalar>
 
564
inline ntlVector3Dim<Scalar>
 
565
ntlVector3Dim<Scalar>::operator-(Scalar s ) const
 
566
{
 
567
  return ntlVector3Dim<Scalar>(value[0]-s,
 
568
                        value[1]-s,
 
569
                        value[2]-s);
 
570
}
 
571
 
 
572
 
 
573
 
 
574
/*************************************************************************
 
575
  Build a vector with another vector multiplied by componentwise.
 
576
  \param v The second vector to muliply with
 
577
  \return The product vector
 
578
  */
 
579
template<class Scalar>
 
580
inline ntlVector3Dim<Scalar>
 
581
ntlVector3Dim<Scalar>::operator*( const ntlVector3Dim<Scalar>& v) const
 
582
{
 
583
  return ntlVector3Dim<Scalar>(value[0]*v.value[0],
 
584
                        value[1]*v.value[1],
 
585
                        value[2]*v.value[2]);
 
586
}
 
587
 
 
588
 
 
589
/*************************************************************************
 
590
  Build a ntlVector3Dim with a Scalar value multiplied to each component.
 
591
  \param s The Scalar value to multiply with
 
592
  \return The product vector
 
593
  */
 
594
template<class Scalar>
 
595
inline ntlVector3Dim<Scalar>
 
596
ntlVector3Dim<Scalar>::operator*(Scalar s) const
 
597
{
 
598
  return ntlVector3Dim<Scalar>(value[0]*s, value[1]*s, value[2]*s);
 
599
}
 
600
 
 
601
 
 
602
/*************************************************************************
 
603
  Build a vector divided componentwise by another vector.
 
604
  \param v The second vector to divide by
 
605
  \return The ratio vector
 
606
  */
 
607
template<class Scalar>
 
608
inline ntlVector3Dim<Scalar>
 
609
ntlVector3Dim<Scalar>::operator/(const ntlVector3Dim<Scalar>& v) const
 
610
{
 
611
  return ntlVector3Dim<Scalar>(value[0]/v.value[0],
 
612
                        value[1]/v.value[1],
 
613
                        value[2]/v.value[2]);
 
614
}
 
615
 
 
616
 
 
617
 
 
618
/*************************************************************************
 
619
  Build a vector divided componentwise by a Scalar value.
 
620
  \param s The Scalar value to divide by
 
621
  \return The ratio vector
 
622
  */
 
623
template<class Scalar>
 
624
inline ntlVector3Dim<Scalar>
 
625
ntlVector3Dim<Scalar>::operator/(Scalar s) const
 
626
{
 
627
  return ntlVector3Dim<Scalar>(value[0]/s,
 
628
                        value[1]/s,
 
629
                        value[2]/s);
 
630
}
 
631
 
 
632
 
 
633
 
 
634
 
 
635
 
 
636
/*************************************************************************
 
637
  Get a particular component of the vector.
 
638
  \param i Number of Scalar to get
 
639
  \return Reference to the component
 
640
  */
 
641
template<class Scalar>
 
642
inline Scalar&
 
643
ntlVector3Dim<Scalar>::operator[]( unsigned int i )
 
644
{
 
645
  return value[i];
 
646
}
 
647
 
 
648
 
 
649
/*************************************************************************
 
650
  Get a particular component of a constant vector.
 
651
  \param i Number of Scalar to get
 
652
  \return Reference to the component
 
653
  */
 
654
template<class Scalar>
 
655
inline const Scalar&
 
656
ntlVector3Dim<Scalar>::operator[]( unsigned int i ) const
 
657
{
 
658
  return value[i];
 
659
}
 
660
 
 
661
 
 
662
 
 
663
//------------------------------------------------------------------------------
 
664
// BLITZ compatibility functions
 
665
//------------------------------------------------------------------------------
 
666
 
 
667
 
 
668
 
 
669
/*************************************************************************
 
670
  Compute the scalar product with another vector.
 
671
  \param v The second vector to work with
 
672
  \return The value of the scalar product
 
673
  */
 
674
template<class Scalar>
 
675
inline Scalar dot(const ntlVector3Dim<Scalar> &t, const ntlVector3Dim<Scalar> &v )
 
676
{
 
677
  //return t.value[0]*v.value[0] + t.value[1]*v.value[1] + t.value[2]*v.value[2];
 
678
  return ((t[0]*v[0]) + (t[1]*v[1]) + (t[2]*v[2]));
 
679
}
 
680
 
 
681
 
 
682
/*************************************************************************
 
683
  Calculate the cross product of this and another vector
 
684
 */
 
685
template<class Scalar>
 
686
inline ntlVector3Dim<Scalar> cross(const ntlVector3Dim<Scalar> &t, const ntlVector3Dim<Scalar> &v)
 
687
{
 
688
  ntlVector3Dim<Scalar> cp( 
 
689
                        ((t[1]*v[2]) - (t[2]*v[1])),
 
690
                  ((t[2]*v[0]) - (t[0]*v[2])),
 
691
                  ((t[0]*v[1]) - (t[1]*v[0])) );
 
692
  return cp;
 
693
}
 
694
 
 
695
 
 
696
 
 
697
 
 
698
/*************************************************************************
 
699
  Compute a vector that is orthonormal to self. Nothing else can be assumed
 
700
  for the direction of the new vector.
 
701
  \return The orthonormal vector
 
702
  */
 
703
template<class Scalar>
 
704
ntlVector3Dim<Scalar>
 
705
ntlVector3Dim<Scalar>::getOrthogonalntlVector3Dim() const
 
706
{
 
707
  // Determine the  component with max. absolute value
 
708
  int max= (fabs(value[0]) > fabs(value[1])) ? 0 : 1;
 
709
  max= (fabs(value[max]) > fabs(value[2])) ? max : 2;
 
710
 
 
711
  /*************************************************************************
 
712
    Choose another axis than the one with max. component and project
 
713
    orthogonal to self
 
714
    */
 
715
  ntlVector3Dim<Scalar> vec(0.0);
 
716
  vec[(max+1)%3]= 1;
 
717
  vec.normalize();
 
718
  vec.projectNormalTo(this->getNormalized());
 
719
  return vec;
 
720
}
 
721
 
 
722
 
 
723
/*************************************************************************
 
724
  Projects the vector into a plane normal to the given vector, which must
 
725
  have unit length. Self is modified.
 
726
  \param v The plane normal
 
727
  \return The projected vector
 
728
  */
 
729
template<class Scalar>
 
730
inline const ntlVector3Dim<Scalar>&
 
731
ntlVector3Dim<Scalar>::projectNormalTo(const ntlVector3Dim<Scalar> &v)
 
732
{
 
733
  Scalar sprod = dot(*this,v);
 
734
  value[0]= value[0] - v.value[0] * sprod;
 
735
  value[1]= value[1] - v.value[1] * sprod;
 
736
  value[2]= value[2] - v.value[2] * sprod;  
 
737
  return *this;
 
738
}
 
739
 
 
740
 
 
741
 
 
742
//------------------------------------------------------------------------------
 
743
// Other helper functions
 
744
//------------------------------------------------------------------------------
 
745
 
 
746
 
 
747
 
 
748
/*************************************************************************
 
749
  Minimize the vector, i.e. set each entry of the vector to the minimum
 
750
  of both values.
 
751
  \param pnt The second vector to compare with
 
752
  \return Reference to the modified self
 
753
  */
 
754
template<class Scalar>
 
755
inline const ntlVector3Dim<Scalar> &
 
756
ntlVector3Dim<Scalar>::minimize(const ntlVector3Dim<Scalar> &pnt)
 
757
{
 
758
  for (unsigned int i = 0; i < 3; i++)
 
759
    value[i] = MIN(value[i],pnt[i]);
 
760
  return *this;
 
761
}
 
762
 
 
763
 
 
764
 
 
765
/*************************************************************************
 
766
  Maximize the vector, i.e. set each entry of the vector to the maximum
 
767
  of both values.
 
768
  \param pnt The second vector to compare with
 
769
  \return Reference to the modified self
 
770
  */
 
771
template<class Scalar>
 
772
inline const ntlVector3Dim<Scalar> &
 
773
ntlVector3Dim<Scalar>::maximize(const ntlVector3Dim<Scalar> &pnt)
 
774
{
 
775
  for (unsigned int i = 0; i < 3; i++)
 
776
    value[i] = MAX(value[i],pnt[i]);
 
777
  return *this;
 
778
}
 
779
 
 
780
 
 
781
 
 
782
 
 
783
// ----
 
784
 
 
785
// a 3D vector with double precision
 
786
typedef ntlVector3Dim<double>  ntlVec3d; 
 
787
 
 
788
// a 3D vector with single precision
 
789
typedef ntlVector3Dim<float>   ntlVec3f; 
 
790
 
 
791
// a 3D integer vector
 
792
typedef ntlVector3Dim<int>     ntlVec3i; 
 
793
 
 
794
// Color uses single precision fp values
 
795
typedef ntlVec3f ntlColor;
 
796
 
 
797
/* convert a float to double vector */
 
798
template<class T> inline ntlVec3d vec2D(T v) { return ntlVec3d(v[0],v[1],v[2]); }
 
799
template<class T> inline ntlVec3f vec2F(T v) { return ntlVec3f(v[0],v[1],v[2]); }
 
800
template<class T> inline ntlColor vec2Col(T v) { return ntlColor(v[0],v[1],v[2]); }
 
801
 
 
802
 
 
803
 
 
804
/************************************************************************/
 
805
// graphics vector typing
 
806
 
 
807
 
 
808
// use which fp-precision for raytracing? 1=float, 2=double
 
809
 
 
810
/* VECTOR_EPSILON is the minimal vector length
 
811
   In order to be able to discriminate floating point values near zero, and
 
812
   to be sure not to fail a comparison because of roundoff errors, use this
 
813
   value as a threshold.  */
 
814
 
 
815
// use which fp-precision for graphics? 1=float, 2=double
 
816
#ifdef PRECISION_GFX_SINGLE
 
817
#define GFX_PRECISION 1
 
818
#else
 
819
#ifdef PRECISION_GFX_DOUBLE
 
820
#define GFX_PRECISION 2
 
821
#else
 
822
// standard precision for graphics 
 
823
#ifndef GFX_PRECISION
 
824
#define GFX_PRECISION 1
 
825
#endif
 
826
#endif
 
827
#endif
 
828
                
 
829
#if GFX_PRECISION==1
 
830
typedef float gfxReal;
 
831
#define GFX_REAL_MAX __FLT_MAX__
 
832
//#define vecF2Gfx(x) (x)
 
833
//#define vecGfx2F(x) (x)
 
834
//#define vecD2Gfx(x) vecD2F(x)
 
835
//#define vecGfx2D(x) vecF2D(x)
 
836
#define VECTOR_EPSILON (1e-5f)
 
837
#else
 
838
typedef double gfxReal;
 
839
#define GFX_REAL_MAX __DBL_MAX__
 
840
//#define vecF2Gfx(x) vecD2F(x)
 
841
//#define vecGfx2F(x) vecF2D(x)
 
842
//#define vecD2Gfx(x) (x)
 
843
//#define vecGfx2D(x) (x)
 
844
#define VECTOR_EPSILON (1e-10)
 
845
#endif
 
846
 
 
847
/* fixed double prec. type, for epxlicitly double values */
 
848
typedef double gfxDouble;
 
849
 
 
850
// a 3D vector for graphics output, typically float?
 
851
typedef ntlVector3Dim<gfxReal>  ntlVec3Gfx; 
 
852
 
 
853
template<class T> inline ntlVec3Gfx vec2G(T v) { return ntlVec3Gfx(v[0],v[1],v[2]); }
 
854
 
 
855
/* get minimal vector length value that can be discriminated.  */
 
856
//inline double getVecEpsilon() { return (double)VECTOR_EPSILON; }
 
857
inline gfxReal getVecEpsilon() { return (gfxReal)VECTOR_EPSILON; }
 
858
 
 
859
#define HAVE_GFXTYPES
 
860
 
 
861
 
 
862
 
 
863
 
 
864
/************************************************************************/
 
865
// HELPER FUNCTIONS, independent of implementation
 
866
/************************************************************************/
 
867
 
 
868
#define VECTOR_TYPE ntlVector3Dim<Scalar>
 
869
 
 
870
 
 
871
/*************************************************************************
 
872
  Compute the length (norm) of the vector.
 
873
  \return The value of the norm
 
874
  */
 
875
template<class Scalar>
 
876
inline Scalar norm( const VECTOR_TYPE &v)
 
877
{
 
878
  Scalar l = v[0]*v[0] + v[1]*v[1] + v[2]*v[2];
 
879
  return (fabs(l-1.) < VECTOR_EPSILON*VECTOR_EPSILON) ? 1. : sqrt(l);
 
880
}
 
881
 
 
882
 
 
883
/*************************************************************************
 
884
  Same as getNorm but doesnt sqrt  
 
885
  */
 
886
template<class Scalar>
 
887
inline Scalar normNoSqrt( const VECTOR_TYPE &v)
 
888
{
 
889
  return v[0]*v[0] + v[1]*v[1] + v[2]*v[2];
 
890
}
 
891
 
 
892
 
 
893
/*************************************************************************
 
894
  Compute a normalized vector based on this vector.
 
895
  \return The new normalized vector
 
896
  */
 
897
template<class Scalar>
 
898
inline VECTOR_TYPE getNormalized( const VECTOR_TYPE &v)
 
899
{
 
900
  Scalar l = v[0]*v[0] + v[1]*v[1] + v[2]*v[2];
 
901
  if (fabs(l-1.) < VECTOR_EPSILON*VECTOR_EPSILON)
 
902
    return v; /* normalized "enough"... */
 
903
  else if (l > VECTOR_EPSILON*VECTOR_EPSILON)
 
904
  {
 
905
    Scalar fac = 1./sqrt(l);
 
906
    return VECTOR_TYPE(v[0]*fac, v[1]*fac, v[2]*fac);
 
907
  }
 
908
  else
 
909
    return VECTOR_TYPE((Scalar)0);
 
910
}
 
911
 
 
912
 
 
913
/*************************************************************************
 
914
  Compute the norm of the vector and normalize it.
 
915
  \return The value of the norm
 
916
  */
 
917
template<class Scalar>
 
918
inline Scalar normalize( VECTOR_TYPE &v) 
 
919
{
 
920
  Scalar norm;
 
921
  Scalar l = v[0]*v[0] + v[1]*v[1] + v[2]*v[2];  
 
922
  if (fabs(l-1.) < VECTOR_EPSILON*VECTOR_EPSILON) {
 
923
    norm = 1.;
 
924
        } else if (l > VECTOR_EPSILON*VECTOR_EPSILON) {
 
925
    norm = sqrt(l);
 
926
    Scalar fac = 1./norm;
 
927
    v[0] *= fac;
 
928
    v[1] *= fac;
 
929
    v[2] *= fac; 
 
930
        } else {
 
931
    v[0]= v[1]= v[2]= 0;
 
932
    norm = 0.;
 
933
  }
 
934
  return (Scalar)norm;
 
935
}
 
936
 
 
937
 
 
938
/*************************************************************************
 
939
  Compute a vector, that is self (as an incoming
 
940
  vector) reflected at a surface with a distinct normal vector. Note
 
941
  that the normal is reversed, if the scalar product with it is positive.
 
942
  \param n The surface normal
 
943
  \return The new reflected vector
 
944
  */
 
945
template<class Scalar>
 
946
inline VECTOR_TYPE reflectVector(const VECTOR_TYPE &t, const VECTOR_TYPE &n) 
 
947
{
 
948
  VECTOR_TYPE nn= (dot(t, n) > 0.0) ? (n*-1.0) : n;
 
949
  return ( t - nn * (2.0 * dot(nn, t)) );
 
950
}
 
951
 
 
952
 
 
953
 
 
954
/*************************************************************************
 
955
 * My own refraction calculation
 
956
 * Taken from Glassner's book, section 5.2 (Heckberts method)
 
957
 */
 
958
template<class Scalar>
 
959
inline VECTOR_TYPE refractVector(const VECTOR_TYPE &t, const VECTOR_TYPE &normal, Scalar nt, Scalar nair, int &refRefl) 
 
960
{
 
961
        Scalar eta = nair / nt;
 
962
        Scalar n = -dot(t, normal);
 
963
        Scalar tt = 1.0 + eta*eta* (n*n-1.0);
 
964
        if(tt<0.0) {
 
965
                // we have total reflection!
 
966
                refRefl = 1;
 
967
        } else {
 
968
                // normal reflection
 
969
                tt = eta*n - sqrt(tt);
 
970
                return( t*eta + normal*tt );
 
971
        }
 
972
        return t;
 
973
}
 
974
        /*double eta = nair / nt;
 
975
        double n = -((*this) | normal);
 
976
        double t = 1.0 + eta*eta* (n*n-1.0);
 
977
        if(t<0.0) {
 
978
                // we have total reflection!
 
979
                refRefl = 1;
 
980
        } else {
 
981
                // normal reflection
 
982
                t = eta*n - sqrt(t);
 
983
                return( (*this)*eta + normal*t );
 
984
        }
 
985
        return (*this);*/
 
986
 
 
987
 
 
988
/*************************************************************************
 
989
  Test two ntlVector3Dims for equality based on the equality of their
 
990
  values within a small threshold.
 
991
  \param c The second vector to compare
 
992
  \return TRUE if both are equal
 
993
  \sa getEpsilon()
 
994
  */
 
995
template<class Scalar>
 
996
inline bool equal(const VECTOR_TYPE &v, const VECTOR_TYPE &c)
 
997
{
 
998
  return (ABS(v[0]-c[0]) + 
 
999
          ABS(v[1]-c[1]) + 
 
1000
          ABS(v[2]-c[2]) < VECTOR_EPSILON);
 
1001
}
 
1002
 
 
1003
 
 
1004
/*************************************************************************
 
1005
 * Assume this vector is an RGB color, and convert it to HSV
 
1006
 */
 
1007
template<class Scalar>
 
1008
inline void rgbToHsv( VECTOR_TYPE &V )
 
1009
{
 
1010
        Scalar h=0,s=0,v=0;
 
1011
        Scalar maxrgb, minrgb, delta;
 
1012
        // convert to hsv...
 
1013
        maxrgb = V[0];
 
1014
        int maxindex = 1;
 
1015
        if(V[2] > maxrgb){ maxrgb = V[2]; maxindex = 2; }
 
1016
        if(V[1] > maxrgb){ maxrgb = V[1]; maxindex = 3; }
 
1017
        minrgb = V[0];
 
1018
        if(V[2] < minrgb) minrgb = V[2];
 
1019
        if(V[1] < minrgb) minrgb = V[1];
 
1020
 
 
1021
        v = maxrgb;
 
1022
        delta = maxrgb-minrgb;
 
1023
 
 
1024
        if(maxrgb > 0) s = delta/maxrgb;
 
1025
        else s = 0;
 
1026
 
 
1027
        h = 0;
 
1028
        if(s > 0) {
 
1029
                if(maxindex == 1) {
 
1030
                        h = ((V[1]-V[2])/delta)  + 0.0; }
 
1031
                if(maxindex == 2) {
 
1032
                        h = ((V[2]-V[0])/delta)  + 2.0; }
 
1033
                if(maxindex == 3) {
 
1034
                        h = ((V[0]-V[1])/delta)  + 4.0; }
 
1035
                h *= 60.0;
 
1036
                if(h < 0.0) h += 360.0;
 
1037
        }
 
1038
 
 
1039
        V[0] = h;
 
1040
        V[1] = s;
 
1041
        V[2] = v;
 
1042
}
 
1043
 
 
1044
/*************************************************************************
 
1045
 * Assume this vector is HSV and convert to RGB
 
1046
 */
 
1047
template<class Scalar>
 
1048
inline void hsvToRgb( VECTOR_TYPE &V )
 
1049
{
 
1050
        Scalar h = V[0], s = V[1], v = V[2];
 
1051
        Scalar r=0,g=0,b=0;
 
1052
        Scalar p,q,t, fracth;
 
1053
        int floorh;
 
1054
        // ...and back to rgb
 
1055
        if(s == 0) {
 
1056
                r = g = b = v; }
 
1057
        else {
 
1058
                h /= 60.0;
 
1059
                floorh = (int)h;
 
1060
                fracth = h - floorh;
 
1061
                p = v * (1.0 - s);
 
1062
                q = v * (1.0 - (s * fracth));
 
1063
                t = v * (1.0 - (s * (1.0 - fracth)));
 
1064
                switch (floorh) {
 
1065
                case 0: r = v; g = t; b = p; break;
 
1066
                case 1: r = q; g = v; b = p; break;
 
1067
                case 2: r = p; g = v; b = t; break;
 
1068
                case 3: r = p; g = q; b = v; break;
 
1069
                case 4: r = t; g = p; b = v; break;
 
1070
                case 5: r = v; g = p; b = q; break;
 
1071
                }
 
1072
        }
 
1073
 
 
1074
        V[0] = r;
 
1075
        V[1] = g;
 
1076
        V[2] = b;
 
1077
}
 
1078
 
 
1079
 
 
1080
 
 
1081
 
 
1082
#endif /* NTL_VECTOR3DIM_HH */