source: GTP/trunk/Lib/Vis/Preprocessing/src/Mutation.cpp @ 2592

Revision 2592, 21.7 KB checked in by bittner, 16 years ago (diff)

havran ray caster update

RevLine 
[1991]1#include "SamplingStrategy.h"
2#include "Ray.h"
3#include "Intersectable.h"
4#include "Preprocessor.h"
5#include "ViewCellsManager.h"
6#include "AxisAlignedBox3.h"
7#include "RssTree.h"
8#include "Vector2.h"
9#include "RndGauss.h"
10#include "Mutation.h"
[2001]11#include "Exporter.h"
[2060]12#include "Environment.h"
[2076]13#include "RayCaster.h"
[1991]14
15#ifdef GTP_INTERNAL
16#include "ArchModeler2MLRT.hxx"
17#endif
18
19namespace GtpVisibilityPreprocessor {
20
21#define MUTATION_USE_CDF 0
[1997]22
[2060]23  //#define SIL_TERMINATION_MUTATION_PROB 0.8f
[1997]24
[1991]25#define EVALUATE_MUTATION_STATS 1
26
[2590]27//#define Q_SEARCH_STEPS 3
[2001]28
[2590]29// VH - commened out !!!!! 8/1/2008
[2583]30#define SORT_RAY_ENTRIES 1
[2011]31
32// use avg ray contribution as importance
33// if 0 the importance is evaluated from the succ of mutations
[2035]34#define USE_AVG_CONTRIBUTION 1
[2011]35
36MutationBasedDistribution::RayEntry &
37MutationBasedDistribution::GetEntry(const int index)
38{
39#if SORT_RAY_ENTRIES
40  return mRays[index];
41#else
42  return mRays[(mBufferStart+index)%mRays.size()];
43#endif
44}
45
[1991]46void
47MutationBasedDistribution::Update(VssRayContainer &vssRays)
48{
49  //  for (int i=0; i < mRays.size(); i++)
50  //    cout<<mRays[i].mMutations<<" ";
51  //  cout<<endl;
[2353]52  cerr<<"Mutation update..."<<endl;
[1991]53  cerr<<"rays = "<<mRays.size()<<endl;
[2123]54  if ((int)mRays.size()) {
[1991]55        cerr<<"Oversampling factors = "<<
56          GetEntry(0).mMutations<<" "<<
57          GetEntry(1).mMutations<<" "<<
58          GetEntry(2).mMutations<<" "<<
59          GetEntry(3).mMutations<<" "<<
60          GetEntry(4).mMutations<<" "<<
61          GetEntry(5).mMutations<<" ... "<<
[2123]62          GetEntry((int)mRays.size()-6).mMutations<<" "<<
63          GetEntry((int)mRays.size()-5).mMutations<<" "<<
64          GetEntry((int)mRays.size()-4).mMutations<<" "<<
65          GetEntry((int)mRays.size()-3).mMutations<<" "<<
66          GetEntry((int)mRays.size()-2).mMutations<<" "<<
67          GetEntry((int)mRays.size()-1).mMutations<<endl;
[1991]68  }
69  int contributingRays = 0;
70
71  int mutationRays = 0;
72  int dummyNcMutations = 0;
73  int dummyCMutations = 0;
[2001]74
75  int reverseCandidates = 0;
[2011]76
77#if 0
78  sort(mRays.begin(), mRays.end());
79  // reset the start of the buffer
80  mBufferStart = 0;
81#endif
82
[1991]83  for (int i=0; i < vssRays.size(); i++) {
84        if (vssRays[i]->mPvsContribution) {
85          // reset the counter of unsuccsseful mutation for a generating ray (if it exists)
86          if (vssRays[i]->mDistribution == MUTATION_BASED_DISTRIBUTION &&
87                  vssRays[i]->mGeneratorId != -1
88                  ) {
89                mRays[vssRays[i]->mGeneratorId].mUnsuccessfulMutations = 0;
90#if EVALUATE_MUTATION_STATS
91                mutationRays++;
92               
[2021]93                Intersectable *newObject = vssRays[i]->mTerminationObject;
[1991]94
[2210]95                Intersectable *oldObject = mRays[vssRays[i]->mGeneratorId].mRay->mTerminationObject;
[2592]96                // the ray generated a contribution although it hits the same object
97                // mark this using a counter
[1991]98                if (oldObject == newObject)
99                  dummyCMutations++;
100#endif
101          }
102          contributingRays++;
103          if (mRays.size() < mMaxRays) {
104                VssRay *newRay = new VssRay(*vssRays[i]);
105                // add this ray
106                newRay->Ref();
107                mRays.push_back(RayEntry(newRay));
108          } else {
[2592]109                // unref the old ray and add the new ray into the mutation buffer
[1991]110                *mRays[mBufferStart].mRay = *vssRays[i];
111                mRays[mBufferStart].mMutations = 0;
[2011]112                mRays[mBufferStart].mUnsuccessfulMutations = 0;
[2001]113                mRays[mBufferStart].ResetReverseMutation();
[1991]114                //              mRays[mBufferStart] = RayEntry(newRay);
115                mBufferStart++;
116                if (mBufferStart >= mMaxRays)
117                  mBufferStart = 0;
118          }
119        } else {
[2592]120          // the ray did not generate any contribution
[1991]121          if (vssRays[i]->mDistribution == MUTATION_BASED_DISTRIBUTION &&
122                  vssRays[i]->mGeneratorId != -1
123                  ) {
[2001]124                // check whether not to store a new backward mutation candidate
[2592]125                // this happens if the new ray is occluded significantly closer than
126                // the generator ray
[2001]127                VssRay *oldRay = mRays[vssRays[i]->mGeneratorId].mRay;
128                VssRay *newRay = vssRays[i];
[1991]129
[2021]130                Intersectable *oldObject = oldRay->mTerminationObject;
[1991]131
[2592]132                // only allow one reverse mutation per generator ray
[2001]133                if (!mRays[newRay->mGeneratorId].HasReverseMutation()) {
134                  if (DotProd(oldRay->GetDir(), newRay->GetDir()) > 0.0f) {
135                  float oldDist = Magnitude(oldRay->mTermination - newRay->mOrigin);
136                  float newDist = Magnitude(newRay->mTermination - newRay->mOrigin);
137                 
[2060]138                  if (newDist < oldDist - oldObject->GetBox().Radius()*mReverseSamplesDistance) {
[2001]139                        Vector3 origin, termination;
140                        if (ComputeReverseMutation(*oldRay, *newRay, origin, termination)) {
141                          mRays[newRay->mGeneratorId].SetReverseMutation(origin, termination);
142                        }
143                       
144                        reverseCandidates++;
145                        //mReverseCandidates
146                  }
147                  }
148                }
149#if EVALUATE_MUTATION_STATS
150                mutationRays++;
151               
[2021]152                Intersectable *newObject = vssRays[i]->mTerminationObject;
[1991]153
154                if (oldObject == newObject)
155                  dummyNcMutations++;
[2001]156#endif
[1991]157          }
158        }
159  }
[2001]160 
[1991]161  if (mutationRays) {
162        cout<<"Mutated rays:"<<mutationRays<<endl;
163        cout<<"Dummy mutations ratio:"<<100.0f*(dummyCMutations + dummyNcMutations)/
164          (float)mutationRays<<"%"<<endl;
165        cout<<"Dummy NC mutations ratio:"<<100.0f*dummyNcMutations/(float)mutationRays<<"%"<<endl;
166        cout<<"Dummy C mutations ratio:"<<100.0f*dummyCMutations/(float)mutationRays<<"%"<<endl;
[2592]167        cout<<"Reverse candidates:"<<100.0f*reverseCandidates/(float)mutationRays<<endl;
[1991]168  }
[2011]169 
[1991]170  float pContributingRays = contributingRays/(float)vssRays.size();
[2011]171 
[2592]172  cout<<"Ratio of contributing rays:"<<pContributingRays<<endl;
[2011]173 
[2060]174  if (mUseUnsuccCountImportance) {
175        // use unsucc mutation samples as feedback on importance
176        for (int i=0; i < mRays.size(); i++) {
177          const float  minImportance = 0.1f;
178          const int minImportanceSamples = 20;
179          mRays[i].mImportance = minImportance +
180                (1-minImportance)*exp(-3.0f*mRays[i].mUnsuccessfulMutations/minImportanceSamples);
181        }
182  } else {
183        float importance = 1.0f;
184        if (mUsePassImportance)
185          importance = 1.0f/(pContributingRays + 1e-5);
186        // set this values for last contributingRays
187        int index = mBufferStart - 1;
[2011]188       
[2060]189        for (int i=0; i < contributingRays; i++, index--) {
190          if (index < 0)
[2124]191                index = (int)mRays.size()-1;
[2060]192          mRays[index].mImportance = importance;
193        }
[2011]194  }
195 
196#if SORT_RAY_ENTRIES
197  long t1 = GetTime();
198  sort(mRays.begin(), mRays.end());
199  // reset the start of the buffer
200  mBufferStart = 0;
[2124]201  mLastIndex = (int)mRays.size();
[2011]202  cout<<"Mutation candidates sorted in "<<TimeDiff(t1, GetTime())<<" ms."<<endl;
203#endif
[2590]204                                                                                                                                                                 
[1991]205#if MUTATION_USE_CDF
206  // compute cdf
207  mRays[0].mCdf = mRays[0].mImportance/(mRays[0].mMutations+1);
208  for (int i=1; i < mRays.size(); i++)
209        mRays[i].mCdf = mRays[i-1].mCdf + mRays[i].mImportance/(mRays[i].mMutations+1);
210 
211  float scale = 1.0f/mRays[i-1].mCdf;
212  for (i=0; i < mRays.size(); i++) {
213        mRays[i].mCdf *= scale;
214  }
215#endif
216 
217  cout<<"Importance = "<<
218        GetEntry(0).mImportance<<" "<<
[2124]219        GetEntry((int)mRays.size()-1).mImportance<<endl;
[2011]220
221  cout<<"Sampling factor = "<<
222        GetEntry(0).GetSamplingFactor()<<" "<<
[2124]223        GetEntry((int)mRays.size()-1).GetSamplingFactor()<<endl;
[2011]224
[1991]225  cerr<<"Mutation update done."<<endl;
226}
227
228
229Vector3
230MutationBasedDistribution::ComputeOriginMutation(const VssRay &ray,
231                                                                                                 const Vector3 &U,
232                                                                                                 const Vector3 &V,
233                                                                                                 const Vector2 vr2,
234                                                                                                 const float radius
235                                                                                                 )
236{
237#if 0
238  return (U*(r[0] - 0.5f) + V*(r[1] - 0.5f))*(2*radius);
239#endif
240
241
242  // Output random variable
243  Vector2 gaussvec2;
244 
245  // Here we apply transform to gaussian, so 2D bivariate
246  // normal distribution
247  //  float sigma = ComputeSigmaFromRadius(radius);
248  float sigma = radius;
249  GaussianOn2D(vr2,
250                           sigma, // input
251                           gaussvec2);  // output
252
253       
254    // Here we tranform the point correctly to 3D space using base
255    // vectors of the 3D space defined by the direction
256    Vector3 shift = gaussvec2.xx * U + gaussvec2.yy * V;
257       
258        //      cout<<shift<<endl;
259        return shift;
260}
261
262Vector3
263MutationBasedDistribution::ComputeTerminationMutation(const VssRay &ray,
264                                                                                                          const Vector3 &U,
265                                                                                                          const Vector3 &V,
266                                                                                                          const Vector2 vr2,
267                                                                                                          const float radius
268                                                                                                          )
269{
270#if 0
271  return (U*(vr2.xx - 0.5f) + V*(vr2.yy - 0.5f))*(4.0f*radius);
272#endif
273  Vector2 gaussvec2;
274#if 1
275  float sigma = radius;
276  GaussianOn2D(vr2,
277                           sigma, // input
278                           gaussvec2);  // output
279  Vector3 shift = gaussvec2.xx * U + gaussvec2.yy * V;
280  //    cout<<shift<<endl;
281  return shift;
282#endif
283#if 0
284  // Here we estimate standard deviation (sigma) from radius
285  float sigma = 1.1f*ComputeSigmaFromRadius(radius);
286  Vector3 vr3(vr2.xx, vr2.yy, RandomValue(0,1));
287  PolarGaussianOnDisk(vr3,
288                                          sigma,
289                                          radius, // input
290                                          gaussvec2); // output
291 
292  // Here we tranform the point correctly to 3D space using base
293  // vectors of the 3D space defined by the direction
294  Vector3 shift = gaussvec2.xx * U + gaussvec2.yy * V;
295 
296  //    cout<<shift<<endl;
297  return shift;
298#endif
299}
300
[2001]301bool
302MutationBasedDistribution::ComputeReverseMutation(
303                                                                                                  const VssRay &oldRay,
304                                                                                                  const VssRay &newRay,
305                                                                                                  Vector3 &origin,
306                                                                                                  Vector3 &termination
307                                                                                                  )
308{
[2063]309#ifdef GTP_INTERNAL
310
[2001]311  // first reconstruct the termination point
312  Vector3 oldDir = Normalize(oldRay.GetDir());
313  Plane3 oldPlane(oldDir, oldRay.mTermination);
314 
315  termination = oldPlane.FindIntersection(newRay.mOrigin,
316                                                                                  newRay.mTermination);
317 
318  // now find the new origin of the ray by casting ray backward from the termination and termining
319  // silhouette point with respect to the occluding object (object containing the newRay termination)
[1991]320
[2001]321  Plane3 newPlane(oldDir, newRay.mTermination);
322
323  Vector3 oldPivot = newPlane.FindIntersection(oldRay.mOrigin,
324                                                                                           oldRay.mTermination);
325
326  Vector3 newPivot = newRay.mTermination;
327  Vector3 line = 2.0f*(oldPivot - newPivot);
328
[2021]329  Intersectable *occluder = newRay.mTerminationObject;
[2001]330 
331  AxisAlignedBox3 box = occluder->GetBox();
[2592]332  // consider slightly larger neighborhood of the occluder
333  // in search for unocclude reverse ray
[2001]334  box.Scale(2.0f);
335 
336  const int packetSize = 4;
337  static int hit_triangles[packetSize];
338  static float dist[packetSize];
339  static Vector3 dirs[packetSize];
340  static Vector3 shifts[packetSize];
[2582]341  static Vector3 origs[packetSize];
[2001]342  // now find the silhouette along the line
343  int i;
344  float left = 0.0f;
345  float right = 1.0f;
[2063]346
[2001]347  // cast rays to find silhouette ray
[2060]348  for (int j=0; j < mSilhouetteSearchSteps; j++) {
[2001]349        for (i=0; i < packetSize; i++) {
350          float r = left + (i+1)*(right-left)/(packetSize+1);
351          shifts[i] = r*line;
352          dirs[i] = Normalize(newPivot + shifts[i] - termination );
[2582]353          origs[i] = termination;
354          // mlrtaStoreRayASEye4(&termination.x, &dirs[i].x, i);
[2001]355        }
[2582]356
357        // mlrtaTraverseGroupASEye4(&box.Min().x, &box.Max().x, hit_triangles, dist);
358        assert(preprocessor->mRayCaster);
359        preprocessor->mRayCaster->CastRaysPacket4(box.Min(),
360                                                                                          box.Max(),
[2583]361                                                                                          origs,
362                                                                                          dirs,
[2582]363                                                                                          hit_triangles,
364                                                                                          dist);       
[2001]365       
366        for (i=0; i < packetSize; i++) {
[2583]367          //cout<<hit_triangles[i]<<endl;
[2001]368          if (hit_triangles[i] == -1) {
[2583]369                // break on first passing ray
[2582]370            break;
[2001]371          }
372        }
373        float rr = left + (i+1)*(right-left)/(packetSize+1);
374        float rl = left + i*(right-left)/(packetSize+1);
375        left = rl;
376        right = rr;
377  }
378
379  float t = right;
380  if (right==1.0f)
381        return false;
382 
383  if (i == packetSize)
384        origin = newPivot + right*line;
385  else
386        origin = newPivot + shifts[i];
387
388  if (0) {
389       
390        static VssRayContainer rRays;
391        static int counter = 0;
392        char filename[256];
393
394        if (counter < 50) {
395          sprintf(filename, "reverse_rays_%03d.x3d", counter++);
396         
397          VssRay tRay(origin, termination, NULL, NULL);
398          rRays.push_back((VssRay *)&oldRay);
399          rRays.push_back((VssRay *)&newRay);
400          rRays.push_back(&tRay);
401         
402          Exporter *exporter = NULL;
403          exporter = Exporter::GetExporter(filename);
404         
405          exporter->SetFilled();
406         
[2021]407          Intersectable *occludee =
408                oldRay.mTerminationObject;
[2001]409         
410          exporter->SetForcedMaterial(RgbColor(0,0,1));
411          exporter->ExportIntersectable(occluder);
412          exporter->SetForcedMaterial(RgbColor(0,1,0));
413          exporter->ExportIntersectable(occludee);
414          exporter->ResetForcedMaterial();
415         
416          exporter->SetWireframe();
417         
418         
419          exporter->ExportRays(rRays, RgbColor(1, 0, 0));
420          delete exporter;
421          rRays.clear();
422        }
423  }
424
[2592]425  return true;
[2063]426#else
[2574]427  cerr << "warning: reverse mutation not supported!" << endl;
[2592]428  return false;
[2063]429#endif
[2001]430 
431 
432  // now the origin and termination is swapped compred to the generator ray
433  // swap(origin, termination);???
434  // -> perhaps not neccessary as the reverse mutation wil only be used once!
435}
436
[1997]437Vector3
438MutationBasedDistribution::ComputeSilhouetteTerminationMutation(const VssRay &ray,
439                                                                                                                                const Vector3 &origin,
440                                                                                                                                const AxisAlignedBox3 &box,
441                                                                                                                                const Vector3 &U,
442                                                                                                                                const Vector3 &V,
443                                                                                                                                const float radius
444                                                                                                                                )
445{
[2063]446#ifdef GTP_INTERNAL
447
[2001]448  const int packetSize = 4;
[1997]449  static int hit_triangles[packetSize];
450  static float dist[packetSize];
451  static Vector3 dirs[packetSize];
[2001]452  static Vector3 shifts[packetSize];
[2582]453  static Vector3 origs[packetSize];
[1997]454  // mutate the
[2071]455  float alpha = RandomValue(0.0f, Real(2.0*M_PI));
[2001]456  //float alpha = vr2.x*2.0f*M_PI;
[1997]457 
458  // direction along which we will mutate the ray
459  Vector3 line = sin(alpha)*U + cos(alpha)*V;
460 
461  //  cout<<line<<endl;
462  // create 16 rays along the selected dir
[2001]463  int i;
464  float left = 0.0f;
465  float right = radius;
[1997]466  // cast rays to find silhouette ray
[2060]467  for (int j=0; j < mSilhouetteSearchSteps; j++) {
[2001]468        for (i=0; i < packetSize; i++) {
469          float r = left + (i+1)*(right-left)/(packetSize+1);
470          shifts[i] = r*line;
471          dirs[i] = Normalize(ray.mTermination + shifts[i] - origin );
[2582]472          origs[i] = origin;
473          //mlrtaStoreRayASEye4(&origin.x, &dirs[i].x, i);
[1997]474        }
[2001]475       
[2582]476        // mlrtaTraverseGroupASEye4(&box.Min().x, &box.Max().x, hit_triangles, dist);
477        assert(preprocessor->mRayCaster);
478        preprocessor->mRayCaster->CastRaysPacket4(box.Min(), box.Max(),
479                                                  origs, dirs, hit_triangles, dist);       
[2001]480       
481        for (i=0; i < packetSize; i++) {
482          if (hit_triangles[i] == -1) {
483                //        if (hit_triangles[i] == -1 || !box.IsInside(origin + dist[i]*dirs[i])) {
484                // break on first passing ray
485                break;
486          }
487        }
488        float rr = left + (i+1)*(right-left)/(packetSize+1);
489        float rl = left + i*(right-left)/(packetSize+1);
490        left = rl;
491        right = rr;
[1997]492  }
493 
494  if (i == packetSize) {
495        //      cerr<<"Warning: hit the same box here should never happen!"<<endl;
496        // shift the ray even a bit more
[2001]497        //cout<<"W"<<i<<endl;
498        //      return (RandomValue(1.0f, 1.5f)*radius)*line;
499        return right*line;
[1997]500  }
[2001]501 
[1997]502  //  cout<<i<<endl;
[2001]503  return shifts[i];
[2063]504
505#else
[2592]506  //cerr << "warning: shiluette mutation not supported" << endl;
507  return Vector3(0, 0, 0);
[2063]508#endif
509 
[1997]510}
511
512
[1991]513bool
514MutationBasedDistribution::GenerateSample(SimpleRay &sray)
515{
516
517  if (mRays.size() == 0) {
518        float rr[5];
[2592]519        // use direction based distribution until we have some mutation candidates
[1991]520        Vector3 origin, direction;
[2076]521
522        for (int i=0; i < 5; i++)
523          rr[i] = RandomValue(0,1);
[1991]524       
525        mPreprocessor.mViewCellsManager->GetViewPoint(origin,
526                                                                                                  Vector3(rr[0], rr[1], rr[2]));
527       
[2592]528       
[1991]529        direction = UniformRandomVector(rr[3], rr[4]);
530       
531        const float pdf = 1.0f;
532        sray = SimpleRay(origin, direction, MUTATION_BASED_DISTRIBUTION, pdf);
533        sray.mGeneratorId = -1;
[2592]534       
[1991]535        return true;
[2592]536  }
[1991]537
538  int index;
539
540#if !MUTATION_USE_CDF
[2011]541#if SORT_RAY_ENTRIES
542  index = mLastIndex - 1;
543  if (index < 0 || index >= mRays.size()-1) {
[2124]544        index = (int)mRays.size() - 1;
[2011]545  } else
546        if (
547                mRays[index].GetSamplingFactor() >= mRays[mLastIndex].GetSamplingFactor()) {
[2592]548          // make another round to equalize the oversampling factor
[2011]549
550          //      cout<<"R2"<<endl;
551          //      cout<<mLastIndex<<endl;
552          //      cout<<index<<endl;
[2124]553          index = (int)mRays.size() - 1;
[2011]554        }
555#else
[1991]556  // get tail of the buffer
557  index = (mLastIndex+1)%mRays.size();
558  if (mRays[index].GetSamplingFactor() >
559          mRays[mLastIndex].GetSamplingFactor()) {
560        // search back for index where this is valid
561        index = (mLastIndex - 1 + mRays.size())%mRays.size();
562        for (int i=0; i < mRays.size(); i++) {
563         
564          //      if (mRays[index].mMutations > mRays[mLastIndex].mMutations)
565          //            break;
566          if (mRays[index].GetSamplingFactor() >
567                  mRays[mLastIndex].GetSamplingFactor() )
568                break;
569          index = (index - 1 + mRays.size())%mRays.size();
570        }
571        // go one step back
572        index = (index+1)%mRays.size();
573  }
[2011]574#endif
[1991]575#else
576  static HaltonSequence iHalton;
577  iHalton.GetNext(1, rr);
578  //rr[0] = RandomValue(0,1);
579  // use binary search to find index with this cdf
[2210]580  int l=0, r=(int)mRays.size()-1;
[1991]581  while(l<r) {
582        int i = (l+r)/2;
583        if (rr[0] < mRays[i].mCdf )
584          r = i;
585        else
586          l = i+1;
587  }
588  index = l;
589  //  if (rr[0] >= mRays[r].mCdf)
590  //    index = r;
591  //  else
592  //    index = l;
593       
594       
595#endif
596  //  cout<<index<<" "<<rr[0]<<" "<<mRays[index].mCdf<<" "<<mRays[(index+1)%mRays.size()].mCdf<<endl;
597
598  mLastIndex = index;
[2011]599  //  Debug<<index<<" "<<mRays[index].GetSamplingFactor()<<endl;
600 
[2001]601  if (mRays[index].HasReverseMutation()) {
602        //cout<<"R "<<mRays[index].mutatedOrigin<<" "<<mRays[index].mutatedTermination<<endl;
603        sray = SimpleRay(mRays[index].mutatedOrigin,
604                                         Normalize(mRays[index].mutatedTermination - mRays[index].mutatedOrigin),
[2105]605                                         MUTATION_BASED_DISTRIBUTION,
606                                         1.0f);
607       
608       
[2001]609        sray.mGeneratorId = index;
610        mRays[index].ResetReverseMutation();
611        mRays[index].mMutations++;
[2002]612        mRays[index].mUnsuccessfulMutations++;
[2574]613       
[2001]614        return true;
615  }
[2574]616
[1991]617  return GenerateMutation(index, sray);
618}
619
620
621
622
623
624bool
625MutationBasedDistribution::GenerateMutationCandidate(const int index,
626                                                                                                         SimpleRay &sray,
627                                                                                                         Intersectable *object,
628                                                                                                         const AxisAlignedBox3 &box
629                                                                                                         )
630{
631  float rr[4];
632
633  VssRay *ray = mRays[index].mRay;
634
[2001]635  //  rr[0] = RandomValue(0.0f,0.99999f);
636  //  rr[1] = RandomValue(0.0f,0.99999f);
637  //  rr[2] = RandomValue(0.0f,0.99999f);
638  //  rr[3] = RandomValue(0.0f,0.99999f);
[1991]639 
640  // mutate the origin
641  Vector3 d = ray->GetDir();
642 
[2022]643  float objectRadius = box.Radius();
[1991]644  //  cout<<objectRadius<<endl;
645  if (objectRadius < Limits::Small)
646        return false;
647 
648  // Compute right handed coordinate system from direction
649  Vector3 U, V;
650  Vector3 nd = Normalize(d);
651  nd.RightHandedBase(U, V);
652 
653  Vector3 origin = ray->mOrigin;
654  Vector3 termination = ray->mTermination; //box.Center(); //ray->mTermination; //box.Center();
[1995]655
[1991]656 
[2022]657  // use probabilitistic approach to decide for the type of mutation
658  float a = RandomValue(0.0f,1.0f);
[2105]659  bool bidirectional = true;
660 
[2060]661  if (mUseSilhouetteSamples && a < mSilhouetteProb) {
[2022]662        termination += ComputeSilhouetteTerminationMutation(*ray,
663                                                                                                                origin,
664                                                                                                                box,
665                                                                                                                U, V,
666                                                                                                                2.0f*objectRadius);
[2105]667        bidirectional = false;
[2022]668  } else {
669        mRays[index].mHalton.GetNext(4, rr);
[2076]670
[2022]671        // fuzzy random mutation
672        origin += ComputeOriginMutation(*ray, U, V,
673                                                                        Vector2(rr[0], rr[1]),
[2060]674                                                                        mMutationRadiusOrigin*objectRadius);
[2022]675       
676        termination += ComputeTerminationMutation(*ray, U, V,
677                                                                                          Vector2(rr[2], rr[3]),
[2060]678                                                                                          mMutationRadiusTermination*objectRadius);
[2022]679  }
680
[1991]681  Vector3 direction = termination - origin;
682 
683  if (Magnitude(direction) < Limits::Small)
684        return false;
685 
686  // shift the origin a little bit
[2592]687#if 1
[1991]688  origin += direction*0.5f;
[2592]689#else
690  // $$JB - 14.1. 2008 increase shift to test HavranRayCaster
691  origin = 0.5f*(origin + termination);
692#endif
[1991]693 
694  direction.Normalize();
695 
696  // $$ jb the pdf is yet not correct for all sampling methods!
697  const float pdf = 1.0f;
698 
699  sray = SimpleRay(origin, direction, MUTATION_BASED_DISTRIBUTION, pdf);
700  sray.mGeneratorId = index;
[2105]701  sray.SetBidirectional(bidirectional);
702 
[2048]703  return true;
[1991]704}
705
706bool
707MutationBasedDistribution::GenerateMutation(const int index, SimpleRay &sray)
708{
709  VssRay *ray = mRays[index].mRay;
710
[2021]711  Intersectable *object = ray->mTerminationObject;
[1991]712 
713  AxisAlignedBox3 box = object->GetBox();
714 
715  if (GenerateMutationCandidate(index, sray, object, box)) {
716    mRays[index].mMutations++;
[2002]717        mRays[index].mUnsuccessfulMutations++;
[1991]718        return true;
719  }
720  return false;
721}
[2060]722 
[1991]723
[2060]724
725MutationBasedDistribution::MutationBasedDistribution(Preprocessor &preprocessor
726                                                                                                         ) :
727  SamplingStrategy(preprocessor)
[1991]728{
[2060]729  mType = MUTATION_BASED_DISTRIBUTION;
730  mBufferStart = 0;
731  mLastIndex = 0;
732
733  Environment::GetSingleton()->GetIntValue("Mutation.bufferSize",
734                                                                                   mMaxRays);
[1991]735 
[2060]736  mRays.reserve(mMaxRays);
[1991]737 
[2060]738  Environment::GetSingleton()->GetFloatValue("Mutation.radiusOrigin",
739                                                                                         mMutationRadiusOrigin);
[1991]740
[2060]741  Environment::GetSingleton()->GetFloatValue("Mutation.radiusTermination",
742                                                                                         mMutationRadiusTermination);
[1991]743 
[2076]744  bool useEnhancedFeatures;
745
746#ifdef GTP_INTERNAL
747  int rayCaster;
[2172]748  Environment::GetSingleton()->GetIntValue("Preprocessor.rayCastMethod",
[2076]749                                                                                   rayCaster);
[2583]750  useEnhancedFeatures = (rayCaster !=
751                                                 RayCaster::INTERNAL_RAYCASTER);
[2076]752#else
753  useEnhancedFeatures = false;
754#endif
[1991]755 
[2076]756  if (useEnhancedFeatures) {
757        Environment::GetSingleton()->GetBoolValue("Mutation.useReverseSamples",
758                                                                                          mUseReverseSamples);
759       
760        Environment::GetSingleton()->GetFloatValue("Mutation.reverseSamplesDistance",
761
762                                                                                           mReverseSamplesDistance);
763
764        Environment::GetSingleton()->GetFloatValue("Mutation.silhouetteProb",
765                                                                                           mSilhouetteProb);
766
767  } else {
768        mUseReverseSamples = false;
769        mReverseSamplesDistance = 1e20f;
770        mSilhouetteProb = 0.0f;
771  }
[2060]772 
773  Environment::GetSingleton()->GetBoolValue("Mutation.useSilhouetteSamples",
774                                                                                        mUseSilhouetteSamples);
[1991]775
[2060]776  Environment::GetSingleton()->GetIntValue("Mutation.silhouetteSearchSteps",
777                                                                                   mSilhouetteSearchSteps);
[1991]778 
[2076]779 
[2060]780  Environment::GetSingleton()->GetBoolValue("Mutation.usePassImportance",
781                                                                                        mUsePassImportance);
[1991]782 
783
[2060]784  Environment::GetSingleton()->GetBoolValue("Mutation.useUnsuccCountImportance",
785                                                                                        mUseUnsuccCountImportance);
[1991]786
787 
[2060]788 
[1991]789}
790
791
[2011]792
[1991]793}
Note: See TracBrowser for help on using the repository browser.