ParticleFilterPlacement Class Reference

Inheritance diagram for ParticleFilterPlacement:

BLManager::SteppingAction BLManager::TrackingAction BLCallback

List of all members.


Detailed Description

class ParticleFilterPlacement is one placement of a ParticleFilter element.

It inherits from G4Decay which does the heavy lifting, all we do here is decide whether or not to decay. We also kill tracks in BLManager::SteppingAction, because you cannot do that inside PostStepGetPhysicalInteractionLength().

Because geometry construction time is too early to manipulate physics processes, it registers a callback with BLManager to do it before the reference particle.

Public Member Functions

 ParticleFilterPlacement (BLCMDparticlefilter *_pf, G4VPhysicalVolume *_pv, G4String _name)
void callback (int type)
G4VProcess * clone ()
G4bool IsApplicable (const G4ParticleDefinition &)
G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)
G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
void UserSteppingAction (const G4Step *step)
void PreUserTrackingAction (const G4Track *track)
void PostUserTrackingAction (const G4Track *track)

Private Attributes

BLCMDparticlefilterfilter
G4VPhysicalVolume * vol
unsigned long nDecayed
unsigned long nKilled
int nWait
int referenceWait


Constructor & Destructor Documentation

ParticleFilterPlacement::ParticleFilterPlacement ( BLCMDparticlefilter _pf,
G4VPhysicalVolume *  _pv,
G4String  _name 
)

References filter, BLManager::getObject(), nDecayed, nKilled, BLManager::registerCallback(), and vol.

00437                                                         : 
00438                                         G4Decay(_name), BLCallback()
00439 {
00440         filter=_pf;
00441         vol=_pv; 
00442         nDecayed = 0;
00443         nKilled = 0;
00444 
00445         // now is too early to manage the physics processes
00446         BLManager::getObject()->registerCallback(this,0);
00447 
00448         // also register a callback after tracking, to display statistics
00449         BLManager::getObject()->registerCallback(this,2);
00450 }


Member Function Documentation

void ParticleFilterPlacement::callback ( int  type  )  [virtual]

Reimplemented from BLCallback.

References filter, BLElement::getName(), IsApplicable(), nDecayed, and nKilled.

00453 {
00454         if(type == 0) {
00455             // add ourselves to the process list for each applicable particle
00456             // and remove the Decay process if necessary.
00457             G4ParticleTable::G4PTblDicIterator *theParticleIterator =
00458                         G4ParticleTable::GetParticleTable()->GetIterator();
00459             theParticleIterator->reset();
00460             while((*theParticleIterator)()) {
00461                 G4ParticleDefinition *pd = theParticleIterator->value();
00462                 if(IsApplicable(*pd)) {
00463                         G4ProcessManager *pmgr = pd->GetProcessManager();
00464 again:                  G4ProcessVector *pv = pmgr->GetProcessList();
00465                         int n = pv->size();
00466                         for(int i=0; i<n; ++i) {
00467                                 G4VProcess *process = (*pv)[i];
00468                                 if(process->GetProcessName()=="Decay") {
00469                                         pmgr->RemoveProcess(process);
00470                                         goto again;
00471                                 }
00472                         }
00473                         pmgr->AddProcess(this);
00474                         pmgr->SetProcessOrdering(this,idxPostStep);
00475                         pmgr->SetProcessOrdering(this,idxAtRest);
00476 //for(int i=0; i<pv->size(); ++i) printf("ProcessName='%s'\n",(*pv)[i]->GetProcessName().c_str());
00477                 }
00478             }
00479         } else if(type == 2) {
00480                 printf("particlefilter %s: %lu decays  %lu killed\n",
00481                         filter->getName().c_str(),nDecayed,nKilled);
00482         }
00483 }

G4VProcess* ParticleFilterPlacement::clone (  )  [inline]

00155 { return this; }

G4bool ParticleFilterPlacement::IsApplicable ( const G4ParticleDefinition &  pd  ) 

References BLCMDparticlefilter::decaySet, and filter.

Referenced by callback().

00486 {
00487         if(filter->decaySet.count(pd.GetParticleName()) > 0)
00488                 return G4Decay::IsApplicable(pd);
00489         return false;
00490 }

G4double ParticleFilterPlacement::AtRestGetPhysicalInteractionLength ( const G4Track &  track,
G4ForceCondition *  condition 
)

References filter, BLElement::getName(), nDecayed, BLCMDparticlefilter::steppingVerbose, and vol.

00494 {
00495         // we are only applied to particlies listed in decaySet.
00496         if(track.GetVolume() == vol || track.GetNextVolume() == vol) {
00497                 ++nDecayed;
00498                 *condition = Forced;
00499                 if(filter->steppingVerbose) 
00500                         printf("particlefilter '%s' forced Decay\n",
00501                                         filter->getName().c_str());
00502                 return DBL_MIN;
00503         }
00504 
00505         *condition = NotForced;
00506         return DBL_MAX;
00507 }

G4double ParticleFilterPlacement::PostStepGetPhysicalInteractionLength ( const G4Track &  track,
G4double  previousStepSize,
G4ForceCondition *  condition 
)

References filter, BLElement::getName(), nDecayed, BLCMDparticlefilter::steppingVerbose, and vol.

00512 {
00513         // we are only applied to particlies listed in decaySet.
00514         if(track.GetVolume() == vol || track.GetNextVolume() == vol) {
00515                 ++nDecayed;
00516                 *condition = Forced;
00517                 if(filter->steppingVerbose) 
00518                         printf("particlefilter '%s' forced Decay\n",
00519                                         filter->getName().c_str());
00520                 return DBL_MIN;
00521         }
00522 
00523         *condition = NotForced;
00524         return DBL_MAX;
00525 }

void ParticleFilterPlacement::UserSteppingAction ( const G4Step *  step  )  [virtual]

Implements BLManager::SteppingAction.

References BEAM, BLEvaluator::evaluate(), filter, BLManager::getObject(), BLCMDparticlefilter::keepSet, BLCMDparticlefilter::killSet, nKilled, nWait, referenceWait, BLCMDparticlefilter::require, BLCMDparticlefilter::steppingVerbose, VISUAL, and vol.

00527 {
00528         // handle all states
00529 
00530         G4Track *track = step->GetTrack();
00531         G4VPhysicalVolume *preVol = 
00532                                 step->GetPreStepPoint()->GetPhysicalVolume();
00533         G4VPhysicalVolume *postVol = 
00534                                 step->GetPostStepPoint()->GetPhysicalVolume();
00535 
00536         // only handle the step that enters the volume
00537         if(preVol == postVol || postVol != vol) return;
00538 
00539 //printf("ParticleFilterPlacement nWait=%d referenceWait=%d\n",nWait,referenceWait);
00540         // wait until track has hit us nWait times.
00541         if(BLManager::getObject()->getState() == BEAM) {
00542                 if(--nWait > 0) return;
00543         } else if(BLManager::getObject()->getState() == VISUAL) {
00544                 if(--nWait > 0) return;
00545         } else {
00546                 if(--referenceWait > 0) return;
00547         }
00548 //printf("ParticleFilterPlacement-2 nWait=%d referenceWait=%d\n",nWait,referenceWait);
00549 
00550         // we get called for every particle transiting our volume, so we
00551         // need to check for presence in killSet and keepSet, then require.
00552         G4ParticleDefinition *def = track->GetDefinition();
00553         if(filter->killSet.count(def->GetParticleName()) > 0 ||
00554                         (filter->keepSet.size() > 0 &&
00555                          filter->keepSet.count(def->GetParticleName()) == 0)) {
00556                 track->SetTrackStatus(fStopAndKill);
00557                 ++nKilled;
00558                 if(filter->steppingVerbose) 
00559                         printf("ParticleFilterPlacement killed track\n");
00560         } else if(filter->require != "") {
00561                 static BLEvaluator *e=0;
00562                 if(!e)
00563                         e = new BLEvaluator();
00564                 G4ThreeVector position = track->GetPosition();
00565                 e->setVariable("x",position[0]);
00566                 e->setVariable("y",position[1]);
00567                 e->setVariable("z",position[2]);
00568                 e->setVariable("t",track->GetGlobalTime());
00569                 G4ThreeVector mom = track->GetMomentum();
00570                 e->setVariable("Px",mom[0]);
00571                 e->setVariable("Py",mom[1]);
00572                 e->setVariable("Pz",mom[2]);
00573                 e->setVariable("Ptot",sqrt(mom[0]*mom[0]+mom[1]*mom[1]+
00574                                                         mom[2]*mom[2]));
00575                 e->setVariable("PDGid",def->GetPDGEncoding());
00576                 e->setVariable("TrkId",BLManager::getObject()->
00577                                 getExternalTrackID(track));
00578                 e->setVariable("Parent",BLManager::getObject()->
00579                                 getExternalParentID(track));
00580                 if(e->evaluate(filter->require) == 0.0) {
00581                         track->SetTrackStatus(fStopAndKill);
00582                         ++nKilled;
00583                         if(filter->steppingVerbose) 
00584                                 printf("ParticleFilterPlacement killed track\n");
00585                 }
00586         }
00587 }

void ParticleFilterPlacement::PreUserTrackingAction ( const G4Track *  track  )  [virtual]

Implements BLManager::TrackingAction.

References filter, BLCMDparticlefilter::nWait, nWait, BLCMDparticlefilter::referenceWait, and referenceWait.

00590 {
00591         nWait = filter->nWait;
00592         referenceWait = filter->referenceWait;
00593 //printf("ParticleFilterPlacement-PreTracking nWait=%d referenceWait=%d\n",nWait,referenceWait);
00594 }

void ParticleFilterPlacement::PostUserTrackingAction ( const G4Track *  track  )  [inline, virtual]

Implements BLManager::TrackingAction.

00163 { }


Member Data Documentation

G4VPhysicalVolume* ParticleFilterPlacement::vol [private]

unsigned long ParticleFilterPlacement::nDecayed [private]

unsigned long ParticleFilterPlacement::nKilled [private]


The documentation for this class was generated from the following file:
g4beamline