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

Revision 2575, 10.0 KB checked in by bittner, 17 years ago (diff)

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