eic-opticks/u4/U4RandomTools.hh

0001 #pragma once
0002
0003 #include <CLHEP/Units/PhysicalConstants.h>
0004
0005 #include "globals.hh"
0006 #include "Randomize.hh"
0007 #include "G4TwoVector.hh"
0008 #include "G4ThreeVector.hh"
0009 #include "U4RandomDirection.hh"
0010
0011 #ifdef DEBUG_TAG
0012 #include "U4Stack.h"
0013 #include "SEvt.hh"
0014 #endif
0015
0016
0017
0018 // ---------------------------------------------------------------------------
0019 // Returns a random lambertian unit vector (rejection sampling)
0020 //
0021 inline G4ThreeVector U4LambertianRand(const G4ThreeVector& normal)
0022 {
0023   G4ThreeVector vect;
0024   G4double ndotv;
0025   G4int count=0;
0026   const G4int max_trials = 1024;
0027
0028   G4double u_exitloop ;
0029
0030   do
0031   {
0032     ++count;
0033     vect = U4RandomDirection();
0034     ndotv = normal * vect;
0035
0036     if (ndotv < 0.0)
0037     {
0038       vect = -vect;
0039       ndotv = -ndotv;
0040     }
0041
0042     u_exitloop = G4UniformRand() ;
0043 #ifndef PRODUCTION
0044 #ifdef DEBUG_TAG
0045     SEvt::AddTag(1, U4Stack_LambertianRand, u_exitloop );
0046 #endif
0047 #endif
0048
0049   } while (!(u_exitloop < ndotv) && (count < max_trials));
0050
0051   return vect;
0052 }
0053
0054