source: GTP/trunk/Lib/Vis/Preprocessing/src/havran/raypack.h @ 2592

Revision 2592, 12.9 KB checked in by bittner, 17 years ago (diff)

havran ray caster update

Line 
1// ============================================================================
2// $Id: raypack.h $
3//
4// raypack.h
5//        CRayPacket class - core of the ray-packet traversal routines
6//
7// Class: CRayPacket2x2
8//
9// REPLACEMENT_STRING
10//
11// Initial coding by Vlastimil Havran, 2006. The data design is in fact
12// Jakko Biker layout as propose in the article on Intel Web Site in year 2005
13// http://www.intel.com/cd/ids/developer/asmo-na/eng/245711.htm?page=1
14
15#ifndef __RAYPACK_H__
16#define __RAYPACK_H__
17
18#include <cassert>
19
20namespace GtpVisibilityPreprocessor {
21
22#include "Vector3.h"
23#include "ktbconf.h"
24 
25#ifdef _USE_HAVRAN_SSE
26 
27#ifdef __SSE__
28
29// System headers for SSE
30#ifdef __INTEL_COMPILER
31#include <xmmintrin.h>
32#else
33// We assume GNU GCC compiler 3.4 or higher
34#include <xmmintrin.h>
35#endif
36
37
38// forward declarations
39#define SSE_INTRINSIC
40#ifdef SSE_INTRINSIC
41#define ALIGN16             __declspec(align(16))
42
43#ifdef _MSC_VER
44#define NEW_ALIGN16(type,n)  ((type*)_aligned_malloc((n)*sizeof(type),16))
45#define FREE_ALIGN16(array)  if(array){_aligned_free(array);(array)=0;}
46#else
47#define NEW_ALIGN16(type,n)  ((type*)malloc((n)*sizeof(type)))
48#define FREE_ALIGN16(array)  if (array) { free(array);(array)=0;}
49#endif // _MSC_VC
50
51#define PAD_FOUR(h)           mulFour(h)
52union extract_m128
53{
54  __m128 m;
55  float f[4];
56};
57#else
58#define NEW_ALIGN16(type,n)  ((type*)malloc((n)*sizeof(type)))
59#define ALIGN16             
60#define FREE_ALIGN16(array)  if(array){free(array);(array)=0;}
61#define PAD_FOUR(h)           (h)
62#endif
63 
64// -------------------------------------------------------------------
65// RayPacket2x2 class.  A set of 4 ray is defined by a location and a
66// direction. The direction is always normalized (length == 1) during
67// use.
68// -------------------------------------------------------------------
69class RayPacket2x2
70{
71public:
72  enum {
73    // The number of rays in one packet
74    PACKSIZE = 4
75  };
76 
77  // constructor
78  RayPacket2x2(
79    // origin and the direction of rays
80    const Vector3 orf[],
81    const Vector3 dirf[],
82    // The same type for all four rays
83    const int _type,
84    // All the rays has to start at the same cell
85    const void *_originCell = NULL,
86    // All the rays has to start at the same object
87    const Intersectable *_startObject = NULL,
88    // and also for shadow rays finish at the same object
89    const Intersectable *_stopObject = NULL
90  )
91  {
92    // location
93    ox[0] = orf[0].x; ox[1] = orf[1].x; ox[2] = orf[2].x; ox[3] = orf[3].x;
94    oy[0] = orf[0].y; oy[1] = orf[1].y; oy[2] = orf[2].y; oy[3] = orf[3].y;
95    oz[0] = orf[0].z; oz[1] = orf[1].z; oz[2] = orf[2].z; oz[3] = orf[3].z;
96    // direction, we assume to be normalized
97    dx[0] = dirf[0].x; dx[1] = dirf[1].x; dx[2] = dirf[2].x; dx[3] = dirf[3].x;
98    dy[0] = dirf[0].y; dy[1] = dirf[1].y; dy[2] = dirf[2].y; dy[3] = dirf[3].y;
99    dz[0] = dirf[0].z; dz[1] = dirf[1].z; dz[2] = dirf[2].z; dz[3] = dirf[3].z;
100    // Other components
101    ttype = _type;
102    origin = _originCell;
103    termination = NULL;
104    startObject = _startObject;
105    stopObject = _stopObject;
106    Init();
107  }
108  // dummy constructor
109  RayPacket2x2() {}
110
111  // Inititalize the ray again when already constructed
112  void Init(
113    // origin and the direction of rays
114    const Vector3 orf[],
115    const Vector3 dirf[],
116    // The same type for all four rays
117    const int _type,
118    // All the rays has to start at the same cell
119    const void *_originCell = NULL,
120    // All the rays has to start at the same object
121    const Intersectable *_startObject = NULL,
122    // and also for shadow rays finish at the same object
123    const Intersectable *_stopObject = NULL,
124    // if the direction of rays is normalized or not
125    bool dirNormalized = false)
126  {
127    // location
128    ox[0] = orf[0].x; ox[1] = orf[1].x; ox[2] = orf[2].x; ox[3] = orf[3].x;
129    oy[0] = orf[0].y; oy[1] = orf[1].y; oy[2] = orf[2].y; oy[3] = orf[3].y;
130    oz[0] = orf[0].z; oz[1] = orf[1].z; oz[2] = orf[2].z; oz[3] = orf[3].z;
131    // direction
132    if (dirNormalized) {
133      // already normalized
134      dx[0] = dirf[0].x; dx[1] = dirf[1].x; dx[2] = dirf[2].x; dx[3] = dirf[3].x;
135      dy[0] = dirf[0].y; dy[1] = dirf[1].y; dy[2] = dirf[2].y; dy[3] = dirf[3].y;
136      dz[0] = dirf[0].z; dz[1] = dirf[1].z; dz[2] = dirf[2].z; dz[3] = dirf[3].z;
137    }
138    else {
139      dx[0] = dirf[0].x; dx[1] = dirf[1].x; dx[2] = dirf[2].x; dx[3] = dirf[3].x;
140      dy[0] = dirf[0].y; dy[1] = dirf[1].y; dy[2] = dirf[2].y; dy[3] = dirf[3].y;
141      dz[0] = dirf[0].z; dz[1] = dirf[1].z; dz[2] = dirf[2].z; dz[3] = dirf[3].z;
142      std::cerr << "Normalization not yet implemented" << std::endl;
143      abort();
144    }
145    // other components
146    ttype = _type;
147    origin = _originCell;
148    termination = NULL;
149    startObject = _startObject;
150    stopObject = _stopObject;
151    Init();
152  }
153
154  // Computes the inverted direction of the rays, used optionally by
155  // a ray traversal algorithm. This has to be reconsidered, if it
156  // is really valuable. !!!
157  void ComputeInvertedDir() const;
158
159  // Computes the sign of the rays and returns false if all the directions
160  // for all three axes are the same, but it could be different among axes,
161  // but for one axis all 4 rays must have the same direction
162  bool ComputeDirSign() const;
163 
164  // the cell in the ASDS, where ray starts from
165  void SetOrigin(const void *c) {origin = c;}
166  const void *GetOrigin() const  { return origin; }
167
168  // the cell in the ASDS, where ray finishes the walk
169  void SetTermination(const void *c) {termination = c; }
170  const void* GetTermination() const { return termination;}
171
172  // the object on which the ray starts at
173  const Intersectable* GetStartObject() const { return startObject;}
174  const Intersectable* GetStopObject() const { return stopObject;}
175
176  void SetStartObject(const Intersectable *newStartObject) {
177    startObject = newStartObject;
178  }
179  void SetStopObject(const Intersectable *newStopObject) {
180    stopObject = newStopObject;
181  }
182  int GetType() const { return ttype; }
183
184  // Reading and Setting origin of the ray and direction
185  // Ray origin
186  inline void SetLoc(int i, const Vector3 &l);
187  Vector3  GetLoc(int i) const;
188  // Direction
189  void SetDir(int i, const Vector3 &ndir);
190  Vector3  GetDir(int i) const;
191
192  // Retuns an object that is intersected by i-th ray
193  Intersectable* GetObject(int i) const;
194  void SetObject(int i, Intersectable* obj);
195  // Retuns a signed distance that is intersected by i-th ray
196  float GetT(int i) const;
197  void SetT(int i, float t);
198  // Retuns maximum signed distance that is intersected by i-th ray
199  float GetMaxT(int i) const;
200  void SetMaxT(int i, float t);
201 
202  // make such operation to slightly change the ray direction
203  // in case any component of ray direction is zero. This is
204  // carried out for all the rays in a packet
205  void CorrectZeroComponents();
206
207  // Returns the sign of the direction if this was precomputed
208  const int &GetSign(int axis) const { return sign_dir[axis];}
209
210  // Reset the result of intersection
211  void ResetObjects() {
212    obj[0] = obj[1] = obj[2] = obj[3] = 0;
213  }
214 
215private:
216  // Here it is crucial the layout of the rays in memory
217  // The layout by Jakko Biker is used
218  typedef float real;
219  union
220  {
221    struct
222    {
223      // ox[N],oy[N],oz[N] - origin of the ray N
224      union { real ox[4]; __m128 ox4; };
225      union { real oy[4]; __m128 oy4; };
226      union { real oz[4]; __m128 oz4; };
227    };
228    __m128 orig[3];
229  };
230  union
231  {
232    struct
233    {
234      // dx[N],dy[N],dz[N] - direction of the ray N
235      union { real dx[4]; __m128 dx4; };
236      union { real dy[4]; __m128 dy4; };
237      union { real dz[4]; __m128 dz4; };
238    };
239    __m128 dir[3];
240  };
241
242#define _USE_INVDIR_RP
243#ifdef _USE_INVDIR_RP
244  union
245  {
246    struct
247    {
248      // idx[N],idy[N],idz[N] - direction of the ray N
249      // inverted dir - maybe an overkill for SSE
250      // to be checked !
251      union { real idx[4]; __m128 idx4; };
252      union { real idy[4]; __m128 idy4; };
253      union { real idz[4]; __m128 idz4; };
254    };
255    __m128 idir[3];
256  };
257#endif
258
259  // The auxiliary and result values for traversal
260     
261  // Here is the result - currently computed signed distance
262  union { real t[4]; __m128 t4; };
263  // and so far minimum signed distance computed. This is required
264  // for computing ray object intersections
265  union { real tmax[4]; __m128 tmax4; };
266  // Here are the pointers to the objects that correspond to tmax[4]
267  // above. They can be different for all the rays !
268  union { Intersectable* obj[4]; __m128 obj4; };
269
270  friend class CKTBTraversal;
271 
272  // Type of the ray: primary, shadow, dummy etc., see ERayType above
273  int ttype;
274
275  // The sign of direction to be used in some algorithms. The sign
276  // has to be the same for all the rays in all the components of the
277  // direction vector !!!!
278  mutable int      sign_dir[3];
279 
280  // I should have some abstract cell data type !!! here
281  // corresponds to the spatial elementary cell
282  const void *origin;
283  const void *termination;
284 
285  // the object on which surface a ray starts from
286  const Intersectable *startObject;
287  // the object on which surface a ray ends, for computation
288  // of the visibility queries between two points
289  const Intersectable *stopObject;
290
291  /// Precompute some CRay parameters. Most of them used for ropes traversal.
292  inline void Init();
293
294  // Precompute some values that are necessary.
295  inline void Precompute();
296};
297
298// --------------------------------------------------------------------------
299//  RayPacket2x2::SetLoc()
300// --------------------------------------------------------------------------
301inline void
302RayPacket2x2::SetLoc(int i, const Vector3 &l)
303{
304  assert( (i>=0) && (i<4));
305  ox[i] = l.x;
306  oy[i] = l.y;
307  oz[i] = l.z;
308}
309
310inline void
311RayPacket2x2::SetDir(int i, const Vector3 &ndr)
312{
313  // We assume that the direction is normalized !!!
314  assert( (i>=0) && (i<4));
315  dx[i] = ndr.x;
316  dy[i] = ndr.y;
317  dz[i] = ndr.z;
318}
319
320inline Vector3
321RayPacket2x2::GetLoc(int i) const
322{
323  assert( (i>=0) && (i<4));
324  return Vector3(ox[i],oy[i],oz[i]);
325}
326
327inline Vector3
328RayPacket2x2::GetDir(int i) const
329{
330  assert( (i>=0) && (i<4));
331  return Vector3(dx[i],dy[i],dz[i]);
332}
333
334// --------------------------------------------------------------------------
335//  RayPacket2x2::Precompute()
336// --------------------------------------------------------------------------
337inline void
338RayPacket2x2::Precompute()
339{
340  // initialize inverted dir ?
341#ifdef _USE_INVDIR_RP
342  //
343  const float epsAdd = -1e-25;
344  idx[0] = 1.0f / (dx[0] + epsAdd);
345  idx[1] = 1.0f / (dx[1] + epsAdd);
346  idx[2] = 1.0f / (dx[2] + epsAdd);
347  idx[3] = 1.0f / (dx[3] + epsAdd);
348 
349  idy[0] = 1.0f / (dy[0] + epsAdd);
350  idy[1] = 1.0f / (dy[1] + epsAdd);
351  idy[2] = 1.0f / (dy[2] + epsAdd);
352  idy[3] = 1.0f / (dy[3] + epsAdd);
353 
354  idz[0] = 1.0f / (dz[0] + epsAdd);
355  idz[1] = 1.0f / (dz[1] + epsAdd);
356  idz[2] = 1.0f / (dz[2] + epsAdd);
357  idz[3] = 1.0f / (dz[3] + epsAdd);
358#endif
359}
360
361// --------------------------------------------------------------------------
362//  RayPacket2x2::Init()
363// --------------------------------------------------------------------------
364inline void
365RayPacket2x2::Init()
366{
367  // apply the standard precomputation
368  Precompute();
369}
370
371// Computes the sign of the rays and returns false if all the directions
372// for all three axes are the same, but it could be different among axes,
373// but for one axis all 4 rays must have the same direction
374inline bool
375RayPacket2x2::ComputeDirSign() const
376{
377  // Set the sign of the direction 1 when negative
378  sign_dir[0] = (dx[0] < 0.0f);
379  sign_dir[1] = (dy[0] < 0.0f);
380  sign_dir[2] = (dz[0] < 0.0f);
381  for (int i = 1; i < 4; i++) {
382    if (sign_dir[0] != (dx[i] < 0.0f))
383      return true; // different direction in x
384    if (sign_dir[1] != (dy[i] < 0.0f))
385      return true; // different direction in y
386    if (sign_dir[2] != (dz[i] < 0.0f))
387      return true; // different direction in z
388  }// for
389
390  // Returns false if all 4 rays have the consistent direction
391  return false;
392}
393
394inline Intersectable*
395RayPacket2x2::GetObject(int i) const
396{
397  assert( (i>=0) && (i<4));
398  return obj[i];
399}
400
401inline void
402RayPacket2x2::SetObject(int i, Intersectable* object)
403{
404  assert( (i>=0) && (i<4));
405  obj[i] = object;
406}
407
408inline float
409RayPacket2x2::GetT(int i) const
410{
411  assert( (i>=0) && (i<4));
412  return t[i];
413}
414
415inline void
416RayPacket2x2::SetT(int i, float tnew)
417{
418  assert( (i>=0) && (i<4));
419  t[i] = tnew;
420}
421
422inline float
423RayPacket2x2::GetMaxT(int i) const
424{
425  assert( (i>=0) && (i<4));
426  return tmax[i];
427}
428
429inline void
430RayPacket2x2::SetMaxT(int i, float tmaxnew)
431{
432  assert( (i>=0) && (i<4));
433  tmax[i] = tmaxnew;
434}
435#endif // __SSE__
436
437#endif // _USE_HAVRAN_SSE
438 
439} // namespace
440
441#endif // __RAYPACK_H__
442
Note: See TracBrowser for help on using the repository browser.