G4WeightCutOffProcess Class Reference

#include <G4WeightCutOffProcess.hh>

Inheritance diagram for G4WeightCutOffProcess:

G4VProcess

Public Member Functions

 G4WeightCutOffProcess (G4double wsurvival, G4double wlimit, G4double isource, G4VIStore *istore, const G4String &aName="WeightCutOffProcess", G4bool para=false)
virtual ~G4WeightCutOffProcess ()
void SetParallelWorld (G4String parallelWorldName)
void SetParallelWorld (G4VPhysicalVolume *parallelWorld)
void StartTracking (G4Track *)
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
const G4StringGetName () const
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)

Detailed Description

Definition at line 59 of file G4WeightCutOffProcess.hh.


Constructor & Destructor Documentation

G4WeightCutOffProcess::G4WeightCutOffProcess ( G4double  wsurvival,
G4double  wlimit,
G4double  isource,
G4VIStore istore,
const G4String aName = "WeightCutOffProcess",
G4bool  para = false 
)

Definition at line 55 of file G4WeightCutOffProcess.cc.

References FatalException, G4cout, G4endl, G4Exception(), G4PathFinder::GetInstance(), G4Step::GetPostStepPoint(), G4Step::GetPreStepPoint(), G4VProcess::GetProcessName(), G4TransportationManager::GetTransportationManager(), G4VProcess::pParticleChange, and G4VProcess::verboseLevel.

00061   : G4VProcess(aName), 
00062     fParticleChange(new G4ParticleChange),
00063     fWeightSurvival(wsurvival),
00064     fWeightLimit(wlimit),
00065     fSourceImportance(isource),
00066     fIStore(istore),
00067     //    fGCellFinder(aGCellFinder),
00068    fGhostNavigator(0), fNavigatorID(-1), fFieldTrack('0')
00069 {
00070   if (!fParticleChange)
00071   {
00072     G4Exception("G4WeightCutOffProcess::G4WeightCutOffProcess()",
00073                 "FatalError", FatalException,
00074                 "Failed to allocate G4ParticleChange !");
00075   }
00076 
00077   G4VProcess::pParticleChange = fParticleChange;
00078 
00079   fGhostStep = new G4Step();
00080   fGhostPreStepPoint = fGhostStep->GetPreStepPoint();
00081   fGhostPostStepPoint = fGhostStep->GetPostStepPoint();
00082 
00083   fTransportationManager = G4TransportationManager::GetTransportationManager();
00084   fPathFinder = G4PathFinder::GetInstance();
00085 
00086   if (verboseLevel>0)
00087   {
00088     G4cout << GetProcessName() << " is created " << G4endl;
00089   }
00090 
00091   paraflag = para;
00092 
00093 }

G4WeightCutOffProcess::~G4WeightCutOffProcess (  )  [virtual]

Definition at line 95 of file G4WeightCutOffProcess.cc.

00096 {
00097   delete fParticleChange;
00098   delete fGhostStep;
00099 }


Member Function Documentation

G4VParticleChange * G4WeightCutOffProcess::AlongStepDoIt ( const G4Track ,
const G4Step  
) [virtual]

Implements G4VProcess.

Definition at line 387 of file G4WeightCutOffProcess.cc.

References G4VParticleChange::Initialize(), and G4VProcess::pParticleChange.

00388 {
00389   // Dummy ParticleChange ie: does nothing
00390   // Expecting G4Transportation to move the track
00391   pParticleChange->Initialize(track);
00392   return pParticleChange; 
00393 
00394   //  return 0;
00395 }

G4double G4WeightCutOffProcess::AlongStepGetPhysicalInteractionLength ( const G4Track ,
G4double  ,
G4double  ,
G4double ,
G4GPILSelection  
) [virtual]

Implements G4VProcess.

Definition at line 298 of file G4WeightCutOffProcess.cc.

References CandidateForSelection, G4Navigator::ComputeSafety(), G4PathFinder::ComputeStep(), DBL_MAX, G4Track::GetCurrentStepNumber(), G4Track::GetVolume(), kDoNot, kSharedOther, kSharedTransport, kUnique, NotCandidateForSelection, and G4FieldTrackUpdator::Update().

00301 {
00302   if(paraflag) {
00303     static G4FieldTrack endTrack('0');
00304     static ELimited eLimited;
00305     
00306     *selection = NotCandidateForSelection;
00307     G4double returnedStep = DBL_MAX;
00308     
00309     if (previousStepSize > 0.)
00310       { fGhostSafety -= previousStepSize; }
00311     //  else
00312     //  { fGhostSafety = -1.; }
00313     if (fGhostSafety < 0.) fGhostSafety = 0.0;
00314     
00315     // ------------------------------------------
00316     // Determination of the proposed STEP LENGTH:
00317     // ------------------------------------------
00318     if (currentMinimumStep <= fGhostSafety && currentMinimumStep > 0.)
00319       {
00320         // I have no chance to limit
00321         returnedStep = currentMinimumStep;
00322         fOnBoundary = false;
00323         proposedSafety = fGhostSafety - currentMinimumStep;
00324       }
00325     else // (currentMinimumStep > fGhostSafety: I may limit the Step)
00326       {
00327         G4FieldTrackUpdator::Update(&fFieldTrack,&track);
00328         //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00329         // ComputeStep
00330         //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00331         returnedStep
00332           = fPathFinder->ComputeStep(fFieldTrack,currentMinimumStep,fNavigatorID,
00333                                      track.GetCurrentStepNumber(),fGhostSafety,eLimited,
00334                                      endTrack,track.GetVolume());
00335         //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00336         if(eLimited == kDoNot)
00337           {
00338             // Track is not on the boundary
00339             fOnBoundary = false;
00340             fGhostSafety = fGhostNavigator->ComputeSafety(endTrack.GetPosition());
00341           }
00342         else
00343           {
00344             // Track is on the boundary
00345             fOnBoundary = true;
00346             proposedSafety = fGhostSafety;
00347           }
00348         //xbug? proposedSafety = fGhostSafety;
00349         if(eLimited == kUnique || eLimited == kSharedOther) {
00350           *selection = CandidateForSelection;
00351         }else if (eLimited == kSharedTransport) { 
00352           returnedStep *= (1.0 + 1.0e-9);  
00353           // Expand to disable its selection in Step Manager comparison
00354         }
00355         
00356       }
00357 
00358   // ----------------------------------------------
00359   // Returns the fGhostSafety as the proposedSafety
00360   // The SteppingManager will take care of keeping
00361   // the smallest one.
00362   // ----------------------------------------------
00363     return returnedStep;
00364 
00365   } else {
00366     return DBL_MAX;
00367     //not sensible!    return -1.0;
00368   }
00369 
00370 }

G4VParticleChange * G4WeightCutOffProcess::AtRestDoIt ( const G4Track ,
const G4Step  
) [virtual]

Implements G4VProcess.

Definition at line 381 of file G4WeightCutOffProcess.cc.

00382 {
00383   return 0;
00384 }

G4double G4WeightCutOffProcess::AtRestGetPhysicalInteractionLength ( const G4Track ,
G4ForceCondition  
) [virtual]

Implements G4VProcess.

Definition at line 374 of file G4WeightCutOffProcess.cc.

00376 {
00377   return -1.0;
00378 }

const G4String & G4WeightCutOffProcess::GetName (  )  const

Definition at line 292 of file G4WeightCutOffProcess.cc.

References G4VProcess::theProcessName.

00293 {
00294   return theProcessName;
00295 }

G4VParticleChange * G4WeightCutOffProcess::PostStepDoIt ( const G4Track ,
const G4Step  
) [virtual]

Implements G4VProcess.

Definition at line 189 of file G4WeightCutOffProcess.cc.

References G4PathFinder::CreateTouchableHandle(), fStopAndKill, G4UniformRand, G4VIStore::GetImportance(), G4StepPoint::GetPhysicalVolume(), G4Step::GetPostStepPoint(), G4VTouchable::GetReplicaNumber(), G4StepPoint::GetTouchable(), G4StepPoint::GetTouchableHandle(), G4Track::GetWeight(), G4ParticleChange::Initialize(), G4VParticleChange::ProposeTrackStatus(), G4VParticleChange::ProposeWeight(), and G4StepPoint::SetTouchableHandle().

00191 {
00192   fParticleChange->Initialize(aTrack);
00193 
00194   if(paraflag) {
00195     fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle();
00196     //xbug?    fOnBoundary = false;
00197     CopyStep(aStep);
00198 
00199     if(fOnBoundary)
00200       {
00201 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00202 // Locate the point and get new touchable
00203 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00204   //??  fPathFinder->Locate(step.GetPostStepPoint()->GetPosition(),
00205   //??                      step.GetPostStepPoint()->GetMomentumDirection());
00206         fNewGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID);
00207 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00208       }
00209     else
00210       {
00211         // Do I need this ??????????????????????????????????????????????????????????
00212         // fGhostNavigator->LocateGlobalPointWithinVolume(track.GetPosition());
00213         // ?????????????????????????????????????????????????????????????????????????
00214         
00215         // fPathFinder->ReLocate(track.GetPosition());
00216         
00217         // reuse the touchable
00218         fNewGhostTouchable = fOldGhostTouchable;
00219       }
00220     
00221     fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
00222     fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
00223 
00224   }
00225 
00226   if(paraflag) {
00227     G4GeometryCell postCell(*(fGhostPostStepPoint->GetPhysicalVolume()), 
00228                             fGhostPostStepPoint->GetTouchable()->GetReplicaNumber());
00229 
00230 
00231     //  G4GeometryCell postCell = fGCellFinder.GetPostGeometryCell(aStep);
00232     //  G4GeometryCell postCell = fGCellFinder.GetPostGeometryCell(fGhostStep);
00233     G4double R = fSourceImportance;
00234     if (fIStore)
00235       {
00236         G4double i = fIStore->GetImportance(postCell);
00237         if (i>0)
00238           {
00239             R/=i;
00240           }
00241       }
00242     G4double w = aTrack.GetWeight();
00243     if (w<R*fWeightLimit)
00244       {
00245         G4double ws = fWeightSurvival*R;
00246         G4double p = w/(ws);
00247         if (G4UniformRand()<p)
00248           {
00249             fParticleChange->ProposeTrackStatus(fStopAndKill);
00250           }
00251         else
00252           {
00253             fParticleChange->ProposeWeight(ws);
00254           }                  
00255       }
00256   } else {
00257 
00258     G4GeometryCell postCell(*(aStep.GetPostStepPoint()->GetPhysicalVolume()), 
00259                             aStep.GetPostStepPoint()->GetTouchable()->GetReplicaNumber());
00260 
00261     //  G4GeometryCell postCell = fGCellFinder.GetPostGeometryCell(aStep);
00262     //  G4GeometryCell postCell = fGCellFinder.GetPostGeometryCell(fGhostStep);
00263     G4double R = fSourceImportance;
00264     if (fIStore)
00265       {
00266         G4double i = fIStore->GetImportance(postCell);
00267         if (i>0)
00268           {
00269             R/=i;
00270           }
00271       }
00272     G4double w = aTrack.GetWeight();
00273     if (w<R*fWeightLimit)
00274       {
00275         G4double ws = fWeightSurvival*R;
00276         G4double p = w/(ws);
00277         if (G4UniformRand()<p)
00278           {
00279             fParticleChange->ProposeTrackStatus(fStopAndKill);
00280           }
00281         else
00282           {
00283             fParticleChange->ProposeWeight(ws);
00284           }                  
00285       }
00286   }
00287 
00288   return fParticleChange;
00289 
00290 }

G4double G4WeightCutOffProcess::PostStepGetPhysicalInteractionLength ( const G4Track aTrack,
G4double  previousStepSize,
G4ForceCondition condition 
) [virtual]

Implements G4VProcess.

Definition at line 177 of file G4WeightCutOffProcess.cc.

References DBL_MAX, and Forced.

00179 {
00180 //   *condition = Forced;
00181 //   return kInfinity;
00182 
00183 //  *condition = StronglyForced;
00184   *condition = Forced;
00185   return DBL_MAX;
00186 }

void G4WeightCutOffProcess::SetParallelWorld ( G4VPhysicalVolume parallelWorld  ) 

Definition at line 120 of file G4WeightCutOffProcess.cc.

References G4VPhysicalVolume::GetName(), and G4TransportationManager::GetNavigator().

00121 {
00122 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00123 // Get pointer of navigator
00124 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00125   fGhostWorldName = parallelWorld->GetName();
00126   fGhostWorld = parallelWorld;
00127   fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld);
00128 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00129 }

void G4WeightCutOffProcess::SetParallelWorld ( G4String  parallelWorldName  ) 

Definition at line 108 of file G4WeightCutOffProcess.cc.

References G4TransportationManager::GetNavigator(), and G4TransportationManager::GetParallelWorld().

Referenced by G4WeightCutOffConfigurator::Configure().

00109 {
00110 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00111 // Get pointers of the parallel world and its navigator
00112 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00113   fGhostWorldName = parallelWorldName;
00114   fGhostWorld = fTransportationManager->GetParallelWorld(fGhostWorldName);
00115   fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld);
00116 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00117 }

void G4WeightCutOffProcess::StartTracking ( G4Track  )  [virtual]

Reimplemented from G4VProcess.

Definition at line 136 of file G4WeightCutOffProcess.cc.

References G4TransportationManager::ActivateNavigator(), G4PathFinder::CreateTouchableHandle(), FatalException, G4Exception(), G4Track::GetMomentumDirection(), G4Track::GetPosition(), G4PathFinder::PrepareNewTrack(), and G4StepPoint::SetTouchableHandle().

00137 {
00138 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00139 // Activate navigator and get the navigator ID
00140 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00141 // G4cout << " G4ParallelWorldScoringProcess::StartTracking" << G4endl;
00142 
00143   if(paraflag) {
00144     if(fGhostNavigator)
00145       { fNavigatorID = fTransportationManager->ActivateNavigator(fGhostNavigator); }
00146     else
00147       {
00148         G4Exception("G4WeightCutOffProcess::StartTracking",
00149                     "ProcParaWorld000",FatalException,
00150                     "G4WeightCutOffProcess is used for tracking without having a parallel world assigned");
00151       }
00152 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00153 
00154 // G4cout << "G4ParallelWorldScoringProcess::StartTracking <<<<<<<<<<<<<<<<<< " << G4endl;
00155 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00156 // Let PathFinder initialize
00157 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00158     fPathFinder->PrepareNewTrack(trk->GetPosition(),trk->GetMomentumDirection());
00159 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00160 
00161 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00162 // Setup initial touchables for the first step
00163 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00164     fOldGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID);
00165     fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
00166     fNewGhostTouchable = fOldGhostTouchable;
00167     fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
00168     
00169     // Initialize 
00170     fGhostSafety = -1.;
00171     fOnBoundary = false;
00172   }
00173 }


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:53:58 2013 for Geant4 by  doxygen 1.4.7