source: GTP/trunk/Lib/Vis/Preprocessing/src/IntelRayCaster.cpp @ 2729

Revision 2729, 11.7 KB checked in by mattausch, 16 years ago (diff)

worked on gvs
Line 
1#include "IntelRayCaster.h"
2#include "VssRay.h"
3#include "Preprocessor.h"
4#include "SceneGraph.h"
5
6#ifdef GTP_INTERNAL
7#include "ArchModeler2MLRT.hxx"
8
9#define DEBUG_RAYCAST 0
10
11
12namespace GtpVisibilityPreprocessor {
13
14
15FILE *fileOut = 0;
16bool saveRays = false;
17bool saveMutationRays = false;
18const int saveRaysStart = 3000000;
19int cntSavedRaysFLUSH = 0;
20unsigned long int cntSavedRays = 0;
21const int intSavedLIMIT = 1024;
22
23void
24InitSaving()
25{
26  fileOut = fopen("fileRays.txt", "wb");
27  if (!fileOut) {
28    cout << "ERROR: file fileRays.txt cannot be opened .. exiting" << endl;
29    exit(3);
30  }
31}
32void
33FinishSaving()
34{
35  fclose(fileOut);
36}
37
38IntelRayCaster::IntelRayCaster(const Preprocessor &preprocessor,
39                               const string externKdTree):
40RayCaster(preprocessor)
41{
42  if (!InitRayCast(externKdTree))
43    cout << "warning: intel ray tracer could not be initialized!" << endl;
44  if (saveRays || saveMutationRays)
45    InitSaving();
46}
47
48
49IntelRayCaster::~IntelRayCaster()
50{
51  if (saveRays)
52    FinishSaving();
53}
54
55
56bool IntelRayCaster::InitRayCast(const string externKdTree)
57{
58  cout<<"Intel ray cast file: " << externKdTree << endl;
59       
60  return mlrtaLoadAS(externKdTree.c_str());
61}
62
63
64// Using packet of 4 rays supposing that these are coherent
65// We give a box to which each ray is clipped to before the
66// ray shooting is computed !
67
68// Using packet of 4 rays supposing that these are coherent
69// We give a box to which each ray is clipped to before the
70// ray shooting is computed !
71void IntelRayCaster::CastRaysPacket4(const Vector3 &minBox,
72                                                                         const Vector3 &maxBox,
73                                                                         const Vector3 origin4[],
74                                                                         const Vector3 dirs4[],
75                                                                         int     result4[],
76                                                                         float   dist4[])
77{
78  for (int i = 0; i < 4; i++) {
79    mlrtaStoreRayASEye4(&origin4[i].x, &dirs4[i].x, i);
80  }
81 
82  rawCastTimer.Entry();
83  mlrtaTraverseGroupASEye4(&minBox.x, &maxBox.x, result4, dist4);
84  rawCastTimer.Exit();
85
86  if (saveMutationRays) {
87    fprintf(fileOut, "I %4.7f %4.7f %4.7f %4.7f %4.7f %4.7f\n",
88                        minBox.x, minBox.y, minBox.z, maxBox.x, maxBox.y, maxBox.z);
89    for (int i = 0; i < 4; i++) {
90      fprintf(fileOut, "%d %4.7f %4.7f %4.7f %4.7f %4.7f %4.7f %d %4.7f\n",
91                          cntSavedRays,
92                          origin4[i].x,
93                          origin4[i].y,
94                          origin4[i].z,
95                          dirs4[i].x,
96                          dirs4[i].y,
97                          dirs4[i].z,
98                          (result4[i] != -1) ? 1 : 0,
99                          (result4[i] != -1) ? dist4[i] : 0);
100      cntSavedRays ++;
101      cntSavedRaysFLUSH++;
102    }
103    if (cntSavedRaysFLUSH > intSavedLIMIT) {
104      fflush(fileOut);
105      cntSavedRaysFLUSH = 0;
106    }
107  } // if we are saving rays
108
109  return;
110}
111
112
113 
114int IntelRayCaster::CastRay(const SimpleRay &simpleRay,
115                                                        VssRayContainer &vssRays,
116                                                        const AxisAlignedBox3 &box,
117                                                        const bool castDoubleRay,
118                                                        const bool pruneInvalidRays
119                                                        )
120{ 
121  //cout << "intel ray" << endl;
122  VssRay *vssRay  = NULL;
123  int hits = 0;
124  int hittriangle;
125  Intersection hitA(simpleRay.mOrigin), hitB(simpleRay.mOrigin);
126 
127  float dist;
128  float dist1, dist2;
129  double normal[3];
130
131  rawCastTimer.Entry();
132 
133  hittriangle = mlrtaIntersectAS(
134                                 &simpleRay.mOrigin.x,
135                                 &simpleRay.mDirection.x,
136                                 normal,
137                                 dist);
138
139  rawCastTimer.Exit();
140  dist1 = dist;
141 
142  if (hittriangle != -1 ) {
143    Intersectable *intersect = mPreprocessor.GetParentObject(hittriangle);
144   
145    if (intersect)
146    {
147      hitA.mObject = intersect;
148      hitA.mNormal = Vector3(normal[0], normal[1], normal[2]);
149      // Get the normal of that face
150      //                Mesh *mesh = ((MeshInstance *)objectA)->GetMesh();
151      //                normalA = mesh->GetFacePlane(mFaceParents[forward_hit_triangles[i]].mFaceIndex).mNormal;
152      //-rays[index+i].mDirection; // $$ temporary
153      hitA.mPoint = simpleRay.Extrap(dist);
154    }
155  }
156 
157  if (castDoubleRay)
158  {
159    Vector3 dir = -simpleRay.mDirection;
160
161        rawCastTimer.Entry();
162
163        hittriangle = mlrtaIntersectAS(
164                                   &simpleRay.mOrigin.x,
165                                   &dir.x,
166                                   normal,
167                                   dist);
168       
169        rawCastTimer.Exit();
170
171    dist2 = dist;
172
173    Intersectable *intersect = mPreprocessor.GetParentObject(hittriangle);
174               
175    if (intersect)
176    {
177      hitB.mObject = intersect;
178      hitB.mNormal = Vector3(normal[0], normal[1], normal[2]);
179      // Get the normal of that face
180      //                Mesh *mesh = ((MeshInstance *)objectB)->GetMesh();
181      //                normalA = mesh->GetFacePlane(mFaceParents[forward_hit_triangles[i]].mFaceIndex).mNormal;
182      //-rays[index+i].mDirection; // $$ temporary
183      hitB.mPoint = simpleRay.Extrap(dist);
184    }
185  }
186
187  if (saveRays && preprocessor->mTotalRaysCast > saveRaysStart) {   
188    if (castDoubleRay)
189        {
190      fprintf(fileOut, "D %d %4.7f %4.7f %4.7f %4.7f %4.7f %4.7f %d %4.7f %d %4.7f\n",
191              cntSavedRays,
192              simpleRay.mOrigin.x,
193              simpleRay.mOrigin.y,
194              simpleRay.mOrigin.z,
195              simpleRay.mDirection.x,
196              simpleRay.mDirection.y,
197              simpleRay.mDirection.z,
198              hitA.mObject ? 1 : 0,
199              hitA.mObject ? dist1 : 0,
200              hitB.mObject ? 1 : 0,
201              hitB.mObject ? dist2 : 0
202              );
203        }
204    else
205      fprintf(fileOut, "S %d %4.7f %4.7f %4.7f %4.7f %4.7f %4.7f %d %4.7f\n",
206              cntSavedRays,
207              simpleRay.mOrigin.x,
208              simpleRay.mOrigin.y,
209              simpleRay.mOrigin.z,
210              simpleRay.mDirection.x,
211              simpleRay.mDirection.y,
212              simpleRay.mDirection.z,
213              hitA.mObject ? 1 : 0,
214              hitA.mObject ? dist1 : 0
215              );
216    cntSavedRays++;
217    cntSavedRaysFLUSH++;
218    if (cntSavedRaysFLUSH > intSavedLIMIT) {
219      fflush(fileOut);
220      cntSavedRaysFLUSH = 0;
221    }
222  }     
223 
224  return ProcessRay(
225                    simpleRay,
226                    hitA,
227                    hitB,
228                    vssRays,
229                    box,
230                    castDoubleRay,
231                    pruneInvalidRays
232                    );
233}
234 
235
236void IntelRayCaster::CastRays16(
237                                SimpleRayContainer &rays,
238                                VssRayContainer &vssRays,
239                                const AxisAlignedBox3 &sbox,
240                                const bool castDoubleRay,
241                                const bool pruneInvalidRays)
242{
243  CastRays16(rays, 0, vssRays, sbox, castDoubleRay, pruneInvalidRays);
244}
245 
246void IntelRayCaster::CastRays16(
247                                SimpleRayContainer &rays,
248                                const int offset,
249                                VssRayContainer &vssRays,
250                                const AxisAlignedBox3 &sbox,
251                                const bool castDoubleRay,
252                                const bool pruneInvalidRays)
253{ 
254  int i, k;
255  const int num = 16;
256 
257#if DEBUG_RAYCAST
258  Debug<<"C16 "<<flush;
259  static int counter=0;
260  Debug<<counter++<<endl;
261#endif
262 
263  static int forward_hit_triangles[16];
264  static float forward_dist[16];
265 
266  static int backward_hit_triangles[16];
267  static float backward_dist[16];
268 
269 
270  Vector3 min = mPreprocessor.mSceneGraph->GetBox().Min();
271  Vector3 max = mPreprocessor.mSceneGraph->GetBox().Max();
272
273  for (k=offset, i=0; i < num; i++, k++) {
274#if DEBUG_RAYCAST
275    if (counter == 43964) {
276      Debug<<rays[k].mOrigin<<" "<<rays[k].mDirection<<endl;
277    }
278#endif
279
280    mlrtaStoreRayAS16(&rays[k].mOrigin.x,
281                      &rays[k].mDirection.x,
282                      i);
283  }
284 
285#if DEBUG_RAYCAST
286  Debug<<"TA\n"<<flush;
287#endif
288  rawCastTimer.Entry();
289
290  mlrtaTraverseGroupAS16(&min.x,
291                         &max.x,
292                         forward_hit_triangles,
293                         forward_dist);
294 
295  rawCastTimer.Exit();
296 
297#if DEBUG_RAYCAST
298  Debug<<"TAB\n"<<flush;
299#endif
300
301  if (castDoubleRay) {
302    for (k=offset, i=0; i < num; i++, k++)  {
303      Vector3 dir = -rays[k].mDirection;
304      mlrtaStoreRayAS16(&rays[k].mOrigin.x,
305                        &dir.x,
306                        i);
307    }
308   
309#if DEBUG_RAYCAST
310    Debug<<"TB\n"<<flush;
311#endif
312    rawCastTimer.Entry();
313    mlrtaTraverseGroupAS16(&min.x,
314                           &max.x,
315                           backward_hit_triangles,
316                           backward_dist);
317
318        rawCastTimer.Exit();
319  }
320       
321#if DEBUG_RAYCAST
322  Debug<<"BBB\n"<<flush;
323#endif
324
325  if (saveRays && preprocessor->mTotalRaysCast > saveRaysStart) {
326    if (castDoubleRay)
327      fprintf(fileOut, "G\n");
328    else
329      fprintf(fileOut, "H\n");
330  }
331 
332  for (i=0, k=offset; i < num; i++, k++)
333  {
334    Intersection hitA(rays[k].mOrigin), hitB(rays[k].mOrigin);
335
336#if DEBUG_RAYCAST
337    Debug<<"FH\n"<<flush;
338#endif
339
340    Intersectable *intersect =
341      mPreprocessor.GetParentObject(forward_hit_triangles[i]);
342
343    if (intersect)
344    {
345      hitA.mObject = intersect;
346      // Get the normal of that face
347      hitA.mNormal = mPreprocessor.GetParentNormal(forward_hit_triangles[i]);
348     
349      //-rays[index+i].mDirection; // $$ temporary
350      hitA.mPoint = rays[k].Extrap(forward_dist[i]);
351    }
352   
353#if DEBUG_RAYCAST
354    Debug<<"BH\n"<<flush;
355#endif
356
357    if (castDoubleRay)
358    {
359      Intersectable *intersect =
360                  mPreprocessor.GetParentObject(backward_hit_triangles[i]);
361     
362      if (intersect)
363      {
364                  hitB.mObject = intersect;
365                  hitB.mNormal = mPreprocessor.GetParentNormal(backward_hit_triangles[i]);
366                 
367                  // normalB = rays[i].mDirection; // $$ temporary
368                  hitB.mPoint = rays[k].Extrap(-backward_dist[i]);
369      }
370    }
371
372
373    if (saveRays && preprocessor->mTotalRaysCast > saveRaysStart) {
374      if (castDoubleRay)     
375        fprintf(fileOut, "%d %4.7f %4.7f %4.7f %4.7f %4.7f %4.7f %d %4.7f %d %4.7f\n",
376                cntSavedRays,
377                rays[k].mOrigin.x,
378                rays[k].mOrigin.y,
379                rays[k].mOrigin.z,
380                rays[k].mDirection.x,
381                rays[k].mDirection.y,
382                rays[k].mDirection.z,
383                hitA.mObject ? 1 : 0,
384                hitA.mObject ? forward_dist[i] : 0,
385                hitB.mObject ? 1 : 0,
386                hitB.mObject ? -backward_dist[i] : 0
387                );
388      else
389        fprintf(fileOut, "%d %4.7f %4.7f %4.7f %4.7f %4.7f %4.7f %d %4.7f\n",
390                cntSavedRays,
391                rays[k].mOrigin.x,
392                rays[k].mOrigin.y,
393                rays[k].mOrigin.z,
394                rays[k].mDirection.x,
395                rays[k].mDirection.y,
396                rays[k].mDirection.z,
397                hitA.mObject ? 1 : 0,
398                hitA.mObject ? forward_dist[i] : 0
399                );
400      cntSavedRays++;
401      cntSavedRaysFLUSH++;
402    }
403
404   
405#if DEBUG_RAYCAST
406    Debug<<"PR\n"<<flush;
407#endif
408
409#if 1
410    ProcessRay(rays[k],
411               hitA,
412               hitB,
413               vssRays,
414               sbox,
415               castDoubleRay,
416               pruneInvalidRays
417               );
418#endif
419  } // for
420
421  if (saveRays) {
422    if (cntSavedRaysFLUSH > intSavedLIMIT) {
423      fflush(fileOut);
424      cntSavedRaysFLUSH = 0;
425    }
426  }
427
428#if DEBUG_RAYCAST
429  Debug<<"C16F\n"<<flush;
430#endif
431}
432
433
434
435
436
437void
438IntelRayCaster::CastSimpleForwardRays(
439                                      SimpleRayContainer &rays,
440                                      const AxisAlignedBox3 &sbox
441                                      )
442{
443  int hit_triangles[16];
444  float dist[16];
445  Vector3 normals[16];
446  Vector3 min = sbox.Min();
447  Vector3 max = sbox.Max();
448 
449  int packets = (int)rays.size() / 16;
450 
451  int i, j, k = 0;
452  Vector3 dir;
453 
454  for (i=0; i < packets; i++) {
455    for (j=0; j < 16; j++, k++)
456      mlrtaStoreRayAS16(&rays[k].mOrigin.x,
457                        &rays[k].mDirection.x,
458                        j);
459    rawCastTimer.Entry();
460
461    mlrtaTraverseGroupAS16(&min.x,
462                           &max.x,
463                           hit_triangles,
464                           dist);       
465        rawCastTimer.Exit();
466  } // for
467
468
469  for (; k < rays.size(); k++) {
470        double normal[3];
471        rawCastTimer.Entry();
472        hit_triangles[0] = mlrtaIntersectAS(
473                                            &rays[k].mOrigin.x,
474                                            &rays[k].mDirection.x,
475                                            normal,
476                                            dist[0]);
477        rawCastTimer.Exit();
478  }
479 
480}
481
482
483void
484IntelRayCaster::CastRays(SimpleRayContainer &rays,
485                                                 VssRayContainer &vssRays,
486                                                 const AxisAlignedBox3 &sbox,
487                                                 bool castDoubleRay,
488                                                 bool pruneInvalidRays)
489{
490
491  int buckets = (int)rays.size() / 16;
492  int offset = 0;
493
494#if 0
495  int time = GetTime();
496  CastSimpleForwardRays(rays, sbox);
497  cout<<1e-3*2*rays.size()/TimeDiff(time, GetTime())<<" Mrays/sec"<<endl;
498#endif
499 
500  for (int i=0; i < buckets; i++, offset+=16) {
501        CastRays16(rays, offset, vssRays, sbox, castDoubleRay, pruneInvalidRays);
502
503        if ((int)rays.size() > 100000 && i % (100000/16) == 0)
504          cout<<"\r"<<offset<<"/"<<(int)rays.size()<<"\r";
505  }
506
507  for (; offset < (int)rays.size(); offset++)
508        CastRay(rays[offset], vssRays, sbox, castDoubleRay, pruneInvalidRays);
509
510}
511
512}
513
514#endif
Note: See TracBrowser for help on using the repository browser.