BLRunManager Class Reference

#include <BLRunManager.hh>

Inheritance diagram for BLRunManager:

BLCallback

List of all members.


Detailed Description

class BLRunManager -- enhanced version of G4RunManager.

There are three basic modes of BLRunManager: normal, collective, and enhanced.

Normal mode: This class inherents from G4RunManager, and can thus be used just like it, for standard control of the Geant4 simulation (via function BeamOn()).

Collective mode:

In collective mode a vector<G4TrackData*> is generated from the input file's beam commands, and then all tracks are stepped for one time step after which the BLCollectiveComputation::collectiveStep() is called; this is continued until all tracks are killed. Any secondaries generated can either be kept or discarded (see setCollectiveMode() below). The TrackID-s of primary tracks are preserved, and TrackIDs for secondary tracks are assigned sequentially from 1001, ignoring EventID.

Collective mode is selected by calling setCollectiveMode(), and then behaving externally just like normal mode (i.e. just call BeamOn()).

Enhanced mode: This class also implements functions that permit enhanced control of the simulation -- in particular there is a processOneTrack() function that simply processes a single track.

NOTE: the word "track" is ambiguous, and is commonly used as both a noun and a verb. Here it is a noun only (referring to a G4Track), and the verb is "process". Processing a track consists of stepping it through the simulation until it is killed or suspended, while holding all secondaries it generates. The secondaries are themselves tracks, and must be dealt with before another call to processOneTrack() (see below).

Enhanced run control: All of the usual G4beamline initialization must occur, so the following code is best put into a BLCallback function registered to replace the main loop -- that means it gets called after the post-Reference callbacks are called, and that after all replace-main-loop callbacks are called the program exits (thus avoiding the main loop).

BLRunManager *runmgr = BLRunManager::getObject(); runmgr->beginRun(int runid); --- loop over events --- runmgr->beginEvent(int evid); --- loop over tracks --- --- get a track runmgr->processOneTrack(track); --- handle secondaries (see below) runmgr->endEvent(); runmgr->endRun(); // write out NTuples, etc. (if appropriate)

There can be a loop over runs, if necessary.

Note that these enhanced run control routines do not create or handle trajectories, and ignore them if present. Hence no visualization of tracks is possible.

Note that secondaries are tracks, but they have not yet been assigned either a trackid or a parentid. Those should be assigned before the track is deferred, as all deferred tracks should have those id-s.

Handling secondaries -- do one of: runmgr->discardAllSecondaries(); or runmgr->processAllSecondariesAndDeferredTracks(); or runmgr->deferAllSecondaries(); // followed sometime later by: runmgr->processAllDeferredTracksAndTheirSecondaries(); or while((track=runmgr->popOneSecondary()) != 0) { // NOTE: you cannot call processOneTrack() here! // -- that would discard any remaining secondaries. if(--- condition ---) runmgr->deferOneTrack(track); else delete track; } // followed sometime later by: runmgr->processAllDeferredTracksAndTheirSecondaries();

Public Member Functions

 BLRunManager ()
 Constructor.
virtual ~BLRunManager ()
 Destructor.
virtual void BeamOn (int nEvents, const char *macroFile=0, G4int n_select=-1)
 BeamOn() will track beam events, in normal or collective mode. See setCollectiveMode() below.
virtual void callback (int type)
 callback() from BLCallback.
void beamOnCollective (int nEvents)
 beamOnCollective() processes an entire run in collective mode. In collective mode, a vector of BLTrackData is constructed, and all tracks in it are stepped one time step at a time. This is essentially an inversion of the usual loop over steps inside loops over tracks and events. All secondaries are immediately added to the end of trackVector, if keepSecondaries is true. The Tune and Reference particles are always tracked in normal mode.
void registerCollectiveComputation (BLCollectiveComputation *compute)
 registerCollectiveComputation() will register a computation to be called after every time step in collective mode. Multiple computations are called in the order they are registered. Must be called before beamOnCollective() or steps may be missed (usually called in the command() function of a BLCommand, so the order of commands in the input file is preserved).
std::vector< BLTrackData > & getTrackVector ()
 getTrackVector() returns the trackVector. Note: the track vector and its tracks can be modified inside the collective computation (takes effect immediately). As BLRunManager is a singleton, the trackVector never changes. No track is ever removed from the track vector, but they can have their status set to fStopAndKill which effectively does that.
int getCurrentTrackIndex ()
 getCurrentTrackIndex() returns the index in the track vector of the track currently being tracked. Returns -1 if not in collective mode, or if tracking is not in progress.
bool getCollectiveMode ()
 getCollectiveMode() returns true if in collective mode.
void setCollectiveMode (bool flag=true, bool _keepSecondaries=true, G4double _deltaT=1.0 *ns)
 setCollectiveMode() arranges for calls to BeamOn() to call beamOnCollective() rather than the usual G4RunManager function. This is for BEAM only; TUNE and REFERENCE particles are always tracked normally via G4RunManager::BeamOn(). _keepSecondaries applies ONLY in collective mode (in normal mode secondaries are always kept, in enhanced mode that is up to the user).
G4double getStepTime ()
 getStepTime() will get the ste time. NOTE: inside collectiveStep() this function returns the time of the current step; outside collectiveStep() it returns the next step time (i.e. the time at which tracks should be suspended).
void setStepTime (G4double t)
 setStepTime() will set the time for the next step in time. If called from within collectiveStep(), it may confuse other classes collectiveStep()-s. If called from beginCollectiveTracking() it sets the time of the first step.
G4double getDeltaT ()
 getDeltaT() will get the time interval between time steps
void setDeltaT (G4double dt)
 setDeltaT() will set the time interval between time steps. All tracks will be stepped to the same value of global time. Initiall value is 1.0*ns. deltaT <= 0.0 is an error. Can be called at any time; if called within collectiveStep() it will affect the immediately following time step.
void rejectCollectiveStep ()
 rejectCollectiveStep() will cause the run manager to abandon tracking (if in progress) and to discard all of the current tracks, replacing the vector with a copy saved immediately after the previous call to collectiveStep() returned. Thus the current collective step will be repeated as if no tracking had been performed since the previous return from collectiveStep(). For this to make sense, setDeltaT() should have been called with a smaller value -- this is intended for a collective algorithm to dynamically adapt the value of deltaT. Note: the first step is to bring all beam tracks to a common value of global time, and rejecting this step is ignored.
void beginRun (int runid=0)
 beginRun() begins a run in enhanced run-control mode.
void endRun ()
 endRun() ends a run in enhanced run-control mode.
void beginEvent (int evid=0)
 beginEvent() begins an event in enhanced run-control mode.
void endEvent ()
 endEvent() ends an event in enhanced run-control mode.
bool getNextBeamEventAndTrack (G4Event **pevent, G4Track **ptrack)
 getNextBeamEventAndTrack() will get the next beam event and track. returns true if one is returned, false if none are left. Neithe event actions nor track actions are performed. Both *pevent and *ptrack must be deleted.
void processOneTrack (G4Track *track)
 processOneTrack() tracks a single track until it is suspended or killed. Note that before tracking begins it deletes any held secondaries, so if secondaries are not to be discarded they must be deferred or processed before the next call. You may want to call BLManager::setState(SPECIAL) to prevent this track from being entered into NTuples for virtualdetectors (etc.) encountered during processing. This function MUST NOT BE CALLED if any track is being processed via any means (e.g. the normal BeamOn() is executing).
void discardAllSecondaries ()
 discardAllSecondaries() will discard all secondaries.
int deferAllSecondaries (int secondaryid=10000, int parentid=-1)
 deferAllSecondaries() will defer all secondaries until deferred tracks are processed. Returns # tracks deferred.
void deferOneTrack (G4Track *track)
 deferOneTrack() will defer a single track until deferred tracks are processed.
int processAllSecondariesAndDeferredTracks (int secondaryid=10000, int parentid=-1)
 processAllSecondariesAndDeferredTracks() will defer all secondaries and then process all deferred tracks and their secondaries. returns # tracks processed.
int processAllDeferredTracksAndTheirSecondaries (int trackid=10000)
 processAllDeferredTracksAndTheirSecondaries() will process all deferred tracks, including processing their secondaries. Returns # tracks processed.
G4Track * popOneSecondary ()
 popOneSecondary() will return a pointer to one secondary, removing it from the list of secondaries; returns NULL if no more. Intended to be called in a loop immediately after a call to processOneTrack(). The returned pointer should eventually be deleted. Order is LIFO.
G4Track * popOneDeferredTrack ()
 popOneDeferredTrack() will return a pointer to one deferred track, removing it from the list of deferred tracks; returns NULL if no more. The returned pointer should eventually be deleted. Note that a loop calling this function can defer other tracks or secondaries, which simply extends the loop. Order is LIFO.
void RunInitialization ()
void DoEventLoop (G4int n_event, const char *macroFile, G4int n_select)
void RunTermination ()
G4EventManager * getEventManager ()
G4TrackingManager * getTrackingManager ()
G4StackManager * getStackManager ()
void setCurrentEvent (G4Event *ev)
G4bool getRunAborted ()

Static Public Member Functions

static BLRunManagergetObject ()
 returns the (singleton) BLRunManager. Note that G4RunManager::GetRunManager() returns a pointr to the base class, not this class.

Private Member Functions

void appendTrack (G4Track *track)

Private Attributes

G4TrackingManager * trackManager
G4StackManager * stackManager
std::vector< BLTrackDatatrackVector
std::vector
< BLCollectiveComputation * > 
computeVector
bool collectiveMode
bool keepSecondaries
G4double stepTime
G4double deltaT
G4int currentTrackIndex
bool rejected
int nextSecondaryTrackID

Static Private Attributes

static BLRunManagersingleton = 0


Constructor & Destructor Documentation

BLRunManager::BLRunManager (  ) 

Constructor.

References BLAssert, collectiveMode, currentTrackIndex, deltaT, keepSecondaries, nextSecondaryTrackID, rejected, singleton, stackManager, and stepTime.

Referenced by getObject().

00112                            : G4RunManager(), trackVector(), computeVector()
00113 {
00114         BLAssert(singleton==0);
00115         singleton = this;
00116 
00117         stackManager = 0;
00118         collectiveMode = false;
00119         keepSecondaries = true;
00120         stepTime = 0.0;
00121         deltaT = 1.0*ns;
00122         currentTrackIndex = -1;
00123         rejected = false;
00124         nextSecondaryTrackID = -9999;
00125 }

BLRunManager::~BLRunManager (  )  [virtual]

Destructor.

References computeVector, and trackVector.

00128 {
00129         trackVector.clear();
00130         computeVector.clear();
00131 }


Member Function Documentation

void BLRunManager::appendTrack ( G4Track *  track  )  [private]

References nextSecondaryTrackID, and trackVector.

Referenced by beamOnCollective().

00322 {
00323         G4int evId = currentEvent->GetEventID();
00324         if(track->GetTrackID() <= 0) {
00325                 track->SetTrackID(nextSecondaryTrackID++);
00326         }
00327 
00328         G4Event *tmp = currentEvent;
00329         currentEvent = new G4Event(*currentEvent);
00330 
00331         trackVector.push_back(BLTrackData(currentEvent,track));
00332 
00333         currentEvent = tmp;
00334 }

BLRunManager * BLRunManager::getObject (  )  [static]

returns the (singleton) BLRunManager. Note that G4RunManager::GetRunManager() returns a pointr to the base class, not this class.

References BLRunManager(), and singleton.

Referenced by BLCMDtrackermode::callback(), BLCMDtrace::command(), BLCMDcollective::command(), BLCMDtracker::fitTrack(), BLManager::handleCallbacks(), BLCMDtracker::handlePreviousTracks(), BLCMDtracker::operator()(), BLCMDtracker::PostUserTrackingAction(), TraceNTuple::TraceNTuple(), and TStepLimiter::TStepLimiter().

00101 {
00102         if(!singleton) {
00103                 new BLRunManager();
00104                 // startup sequencing prevents registering callback, so it
00105                 // is hard-coded in BLManager::handleCallbacks(). This is due
00106                 // to the BLManager constructor calling this function.
00107         }
00108 
00109         return singleton; 
00110 }

void BLRunManager::BeamOn ( int  nEvents,
const char *  macroFile = 0,
G4int  n_select = -1 
) [virtual]

BeamOn() will track beam events, in normal or collective mode. See setCollectiveMode() below.

References BEAM, beamOnCollective(), collectiveMode, and BLManager::getObject().

Referenced by BLManager::trackBeam(), and BLManager::trackTuneAndReferenceParticles().

00133                                                                             {
00134         if(collectiveMode && BLManager::getObject()->getState() == BEAM)
00135                 beamOnCollective(nEvents);
00136         else
00137                 G4RunManager::BeamOn(nEvents,macroFile,n_select);
00138 }

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

callback() from BLCallback.

Reimplemented from BLCallback.

References collectiveMode.

Referenced by BLManager::handleCallbacks().

00141 {
00142         // Only perform actions between reference and Beam
00143         if(type != 1) return;
00144 
00145         // in collectiveMode, add a TStepLimiter to every particle
00146         if(collectiveMode) {
00147                 printf("BLRunManager: adding TStepLimiter processes to all particles\n");
00148                 G4ParticleTable::G4PTblDicIterator *theParticleIterator
00149                         = G4ParticleTable::GetParticleTable()->GetIterator();
00150                 theParticleIterator->reset();
00151                 while((*theParticleIterator)()) {
00152                         G4ParticleDefinition *pd = 
00153                                         theParticleIterator->value();
00154                         if(pd->IsShortLived()) continue;
00155                         G4ProcessManager *pmgr = pd->GetProcessManager();
00156                         if(!pmgr) continue;
00157                         TStepLimiter *tsl = new TStepLimiter();
00158                         pmgr->AddProcess(tsl,-1,-1,3);
00159                         pmgr->SetProcessOrdering(tsl,idxPostStep);
00160                         pmgr->SetProcessOrdering(tsl,idxAtRest);
00161                 }
00162         }
00163 }

void BLRunManager::beamOnCollective ( int  nEvents  ) 

beamOnCollective() processes an entire run in collective mode. In collective mode, a vector of BLTrackData is constructed, and all tracks in it are stepped one time step at a time. This is essentially an inversion of the usual loop over steps inside loops over tracks and events. All secondaries are immediately added to the end of trackVector, if keepSecondaries is true. The Tune and Reference particles are always tracked in normal mode.

References appendTrack(), BEAM, BLManager::BeginOfEventAction(), BLAssert, computeVector, currentTrackIndex, deltaT, BLTrackData::event, BLManager::getNextSecondaryTrackID(), BLManager::getObject(), TStepLimiter::getPrevStatus(), BLManager::getPrimaryTrackID(), IDLE, keepSecondaries, nextSecondaryTrackID, processOneTrack(), rejected, RunInitialization(), RunTermination(), BLManager::setEventID(), BLManager::setEventTimeLimit(), BLManager::setState(), stackManager, stepTime, BLTrackData::track, trackManager, and trackVector.

Referenced by BeamOn().

00166 {
00167         BLManager *manager = BLManager::getObject();
00168         G4TransportationManager *transportManager =
00169                         G4TransportationManager::GetTransportationManager();
00170         G4StateManager* stateManager = G4StateManager::GetStateManager();
00171         numberOfEventToBeProcessed = nEvents;
00172 
00173         if(!ConfirmBeamOnCondition()) {
00174                 G4Exception("BLRunManager","Cannot run beam",FatalException,"");
00175         }
00176         printf("================== In Collective Mode ==================\n");
00177         ConstructScoringWorlds();
00178         RunInitialization();
00179         G4ApplicationState currentState = stateManager->GetCurrentState();
00180         if(currentState != G4State_GeomClosed) {
00181                 G4Exception("BLRunManager","Geometry not closed",FatalException,
00182                                                                         "");
00183         }
00184 
00185         manager->setState(BEAM);
00186         stateManager->SetNewState(G4State_EventProc);
00187         trackManager = 0; // unused in collective mode
00188         stackManager = eventManager->GetStackManager();
00189         runAborted = false;
00190 
00191         // the per-event time limit makes no sense in collective mode
00192         manager->setEventTimeLimit(-1); // infinite
00193 
00194         // event loop to create trackVector
00195         G4PrimaryTransformer *transformer = new G4PrimaryTransformer();
00196         G4Navigator* navigator = transportManager->GetNavigatorForTracking();
00197         stepTime = -DBL_MAX;
00198         int nev;
00199         for(nev=0; nev<nEvents; ++nev) {
00200                 currentEvent = GenerateEvent(nev);
00201                 if(!currentEvent || runAborted) break;
00202                 stackManager->PrepareNewEvent();
00203                 manager->BeginOfEventAction(currentEvent);
00204                 G4TrackVector *tv=transformer->GimmePrimaries(currentEvent,0);
00205                 int evId=currentEvent->GetEventID();
00206                 if(evId <= 0) {
00207                         evId = nev;
00208                         currentEvent->SetEventID(evId);
00209                 }
00210                 nextSecondaryTrackID = manager->getNextSecondaryTrackID();
00211                 for(unsigned j=0; j<tv->size(); ++j) {
00212                         G4Track *track = (*tv)[j];
00213                         BLAssert(j==0); // cannot handle multiple primaries 
00214                         track->SetTrackID(manager->getPrimaryTrackID());
00215                         appendTrack(track);
00216                         G4double t = track->GetGlobalTime();
00217                         if(t > stepTime) stepTime = t;
00218                 }
00219                 if(tv->size() == 0) delete currentEvent;
00220                 tv->clear();
00221                 if(runAborted) break;
00222                 transportManager->SetNavigatorForTracking(navigator);
00223                 manager->setEventID(evId+1);
00224         }
00225 
00226         printf("=========== Collective: %d Events, %ld Tracks ==============\n",
00227                                 nev,(long)trackVector.size());
00228 
00229         for(unsigned i=0; i<computeVector.size(); ++i)
00230                 computeVector[i]->beginCollectiveTracking(trackVector);
00231 
00232         // Step loop
00233         runAborted = false;
00234         int nActive;
00235         bool first=true;
00236         G4double saveStepTime = stepTime;
00237         std::vector<G4Track*> saveVector;
00238         do {
00239                 if(runAborted) goto RunAborted;
00240                 nActive = 0;
00241                 rejected = false;
00242                 // note: secondaries can be appended to trackVector during loop
00243                 for(unsigned i=0; i<trackVector.size(); ++i) {
00244                         if(rejected && !first) break;
00245                         BLTrackData *td = &trackVector[i];
00246                         G4Track *track = td->track;
00247                         G4TrackStatus trackStatus = track->GetTrackStatus();
00248 //printf("Collective Track loop: ev=%d trk=%d status=%d fAlive=%d fSuspend=%d\n",td->event->GetEventID(),track->GetTrackID(),(int)trackStatus, fAlive,fSuspend);
00249                         if(trackStatus != fAlive && 
00250                                                 trackStatus != fStopButAlive)
00251                                 continue;
00252                         ++nActive;
00253                         currentEvent = td->event;
00254                         manager->setEventID(td->event->GetEventID());
00255                         currentTrackIndex = i;
00256                         for(unsigned j=0; j<computeVector.size(); ++j)
00257                                 computeVector[j]->beginTrack(trackVector,i);
00258                         if(runAborted) goto RunAborted;
00259                         processOneTrack(track);
00260                         currentTrackIndex = -1;
00261                         if(td->event->IsAborted())
00262                                 track->SetTrackStatus(fKillTrackAndSecondaries);
00263                         if(runAborted) goto RunAborted;
00264                         trackStatus = track->GetTrackStatus();
00265                         if(trackStatus == fSuspend)
00266                                 track->SetTrackStatus(TStepLimiter::getPrevStatus());
00267                         G4TrackVector *tv = trackManager->GimmeSecondaries();
00268                         if(keepSecondaries && trackStatus != fKillTrackAndSecondaries) {
00269                                 for(unsigned j=0; j<tv->size(); ++j)
00270                                         appendTrack((*tv)[j]);
00271                                 tv->clear();
00272                         } else {
00273                                 for(unsigned j=0; j<tv->size(); ++j)
00274                                         delete (*tv)[j];
00275                                 tv->clear();
00276                         }
00277                 }
00278                 if(runAborted) goto RunAborted;
00279                 // don't call collectiveStep() if rejected during tracking
00280                 if(!rejected || first) {
00281                         for(unsigned i=0; i<computeVector.size(); ++i)
00282                                 computeVector[i]->collectiveStep(trackVector);
00283                 }
00284                 if(rejected && !first) {
00285                         // discard current tracks, and restore saved ones
00286                         unsigned n=saveVector.size();
00287                         // first, discard any added secondaries
00288                         while(trackVector.size() > n) {
00289                                 delete trackVector.back().track;
00290                                 delete trackVector.back().event;
00291                                 trackVector.pop_back();
00292                         }
00293                         BLAssert(trackVector.size() == saveVector.size());
00294                         for(unsigned i=0; i<n; ++i) {
00295                                 delete trackVector[i].track;
00296                                 trackVector[i].track = saveVector[i];
00297                         }
00298                         stepTime = saveStepTime;
00299                 } else {
00300                         // delete saved tracks, and save current ones
00301                         for(unsigned i=0; i<saveVector.size(); ++i)
00302                                 delete saveVector[i];
00303                         saveVector.clear();
00304                         for(unsigned i=0; i<trackVector.size(); ++i)
00305                                 saveVector.push_back(
00306                                         new G4Track(*trackVector[i].track));
00307                         saveStepTime = stepTime;
00308                         stepTime += deltaT;
00309                 }
00310                 first = false;
00311         } while(nActive > 0);
00312 
00313 RunAborted:
00314         for(unsigned i=0; i<computeVector.size(); ++i)
00315                 computeVector[i]->endCollectiveTracking(trackVector);
00316 
00317         RunTermination();
00318         manager->setState(IDLE);
00319 }

void BLRunManager::registerCollectiveComputation ( BLCollectiveComputation compute  )  [inline]

registerCollectiveComputation() will register a computation to be called after every time step in collective mode. Multiple computations are called in the order they are registered. Must be called before beamOnCollective() or steps may be missed (usually called in the command() function of a BLCommand, so the order of commands in the input file is preserved).

References computeVector.

Referenced by BLCMDtrace::command(), and BLCMDcollective::command().

00179                 { computeVector.push_back(compute); }

std::vector<BLTrackData>& BLRunManager::getTrackVector (  )  [inline]

getTrackVector() returns the trackVector. Note: the track vector and its tracks can be modified inside the collective computation (takes effect immediately). As BLRunManager is a singleton, the trackVector never changes. No track is ever removed from the track vector, but they can have their status set to fStopAndKill which effectively does that.

References trackVector.

00187 { return trackVector; }

int BLRunManager::getCurrentTrackIndex (  )  [inline]

getCurrentTrackIndex() returns the index in the track vector of the track currently being tracked. Returns -1 if not in collective mode, or if tracking is not in progress.

References currentTrackIndex.

00192 { return currentTrackIndex; }

bool BLRunManager::getCollectiveMode (  )  [inline]

void BLRunManager::setCollectiveMode ( bool  flag = true,
bool  _keepSecondaries = true,
G4double  _deltaT = 1.0*ns 
) [inline]

setCollectiveMode() arranges for calls to BeamOn() to call beamOnCollective() rather than the usual G4RunManager function. This is for BEAM only; TUNE and REFERENCE particles are always tracked normally via G4RunManager::BeamOn(). _keepSecondaries applies ONLY in collective mode (in normal mode secondaries are always kept, in enhanced mode that is up to the user).

References collectiveMode, keepSecondaries, and setDeltaT().

Referenced by BLCMDcollective::command(), and BLManager::trackTuneAndReferenceParticles().

00206                 { collectiveMode=flag;
00207                   keepSecondaries=_keepSecondaries;
00208                   setDeltaT(_deltaT);
00209                 }

G4double BLRunManager::getStepTime (  )  [inline]

getStepTime() will get the ste time. NOTE: inside collectiveStep() this function returns the time of the current step; outside collectiveStep() it returns the next step time (i.e. the time at which tracks should be suspended).

References stepTime.

Referenced by TStepLimiter::AtRestGetPhysicalInteractionLength(), BLCMDcollective::collectiveStep(), and TStepLimiter::PostStepGetPhysicalInteractionLength().

00215 { return stepTime; }

void BLRunManager::setStepTime ( G4double  t  )  [inline]

setStepTime() will set the time for the next step in time. If called from within collectiveStep(), it may confuse other classes collectiveStep()-s. If called from beginCollectiveTracking() it sets the time of the first step.

References stepTime.

00221 { stepTime = t; }

G4double BLRunManager::getDeltaT (  )  [inline]

getDeltaT() will get the time interval between time steps

References deltaT.

00224 { return deltaT; }

void BLRunManager::setDeltaT ( G4double  dt  )  [inline]

setDeltaT() will set the time interval between time steps. All tracks will be stepped to the same value of global time. Initiall value is 1.0*ns. deltaT <= 0.0 is an error. Can be called at any time; if called within collectiveStep() it will affect the immediately following time step.

References BLAssert, and deltaT.

Referenced by BLCMDcollective::beginCollectiveTracking(), BLCMDcollective::command(), and setCollectiveMode().

00231 { BLAssert(dt>0.0); deltaT = dt; }

void BLRunManager::rejectCollectiveStep (  )  [inline]

rejectCollectiveStep() will cause the run manager to abandon tracking (if in progress) and to discard all of the current tracks, replacing the vector with a copy saved immediately after the previous call to collectiveStep() returned. Thus the current collective step will be repeated as if no tracking had been performed since the previous return from collectiveStep(). For this to make sense, setDeltaT() should have been called with a smaller value -- this is intended for a collective algorithm to dynamically adapt the value of deltaT. Note: the first step is to bring all beam tracks to a common value of global time, and rejecting this step is ignored.

References rejected.

00245 { rejected = true; }

void BLRunManager::beginRun ( int  runid = 0  ) 

beginRun() begins a run in enhanced run-control mode.

References BEAM, BLManager::getObject(), RunInitialization(), and BLManager::setState().

Referenced by BLCMDtrackermode::callback(), and BLCMDtracker::handlePreviousTracks().

00337 {
00338         BLManager::getObject()->setState(BEAM);
00339         ConstructScoringWorlds();
00340         RunInitialization();
00341         SetRunIDCounter(runid);
00342         runAborted = false;
00343 }

void BLRunManager::endRun (  ) 

endRun() ends a run in enhanced run-control mode.

References BLManager::getObject(), IDLE, RunTermination(), and BLManager::setState().

Referenced by BLCMDtrackermode::callback(), and BLCMDtracker::handlePreviousTracks().

00346 {
00347         RunTermination();
00348         BLManager::getObject()->setState(IDLE);
00349 }

void BLRunManager::beginEvent ( int  evid = 0  ) 

beginEvent() begins an event in enhanced run-control mode.

References BLManager::BeginOfEventAction(), BLManager::getObject(), stackManager, and trackManager.

Referenced by BLCMDtrackermode::callback(), and BLCMDtracker::handlePreviousTracks().

00352 {
00353         if(currentEvent) delete currentEvent;
00354         currentEvent = new G4Event(evid);
00355         trackManager = eventManager->GetTrackingManager();
00356         stackManager = eventManager->GetStackManager();
00357         BLManager::getObject()->BeginOfEventAction(currentEvent);
00358 }

void BLRunManager::endEvent (  ) 

endEvent() ends an event in enhanced run-control mode.

References BLManager::EndOfEventAction(), BLManager::getObject(), stackManager, and trackManager.

Referenced by BLCMDtrackermode::callback(), and BLCMDtracker::handlePreviousTracks().

00361 {
00362         BLManager::getObject()->EndOfEventAction(currentEvent);
00363         if(currentEvent) delete currentEvent;
00364         currentEvent = 0;
00365         trackManager = 0;
00366         stackManager = 0;
00367 }

bool BLRunManager::getNextBeamEventAndTrack ( G4Event **  pevent,
G4Track **  ptrack 
)

getNextBeamEventAndTrack() will get the next beam event and track. returns true if one is returned, false if none are left. Neithe event actions nor track actions are performed. Both *pevent and *ptrack must be deleted.

References BEAM, BLAssert, BLManager::getObject(), BLManager::getPrimaryParentID(), BLManager::getPrimaryTrackID(), BLManager::setState(), and SPECIAL.

00370 {
00371         static bool more=true;
00372         static int nev=0;
00373         static BLManager *manager=0;
00374         static G4PrimaryTransformer *transformer=0;
00375         if(!manager) {
00376                 manager = BLManager::getObject();
00377                 transformer = new G4PrimaryTransformer();
00378         }
00379 
00380         if(!more) return false;
00381 
00382         manager->setState(BEAM);
00383         currentEvent = GenerateEvent(nev);
00384         if(!currentEvent || runAborted) {
00385                 more = false;
00386                 return false;
00387         }
00388         *pevent = currentEvent;
00389         G4TrackVector *tv=transformer->GimmePrimaries(currentEvent,0);
00390         int evId=currentEvent->GetEventID();
00391         if(evId <= 0) {
00392                 evId = nev;
00393                 currentEvent->SetEventID(evId);
00394         }
00395         for(unsigned j=0; j<tv->size(); ++j) {
00396                 G4Track *track = (*tv)[j];
00397                 BLAssert(j==0); // cannot handle multiple primaries 
00398                 track->SetTrackID(manager->getPrimaryTrackID());
00399                 track->SetParentID(manager->getPrimaryParentID());
00400                 *ptrack = track;
00401         }
00402         tv->clear();
00403         if(runAborted) {
00404                 more = false;
00405                 return false;
00406         }
00407         ++nev;
00408         manager->setState(SPECIAL);
00409         return true;
00410 }

void BLRunManager::processOneTrack ( G4Track *  track  ) 

processOneTrack() tracks a single track until it is suspended or killed. Note that before tracking begins it deletes any held secondaries, so if secondaries are not to be discarded they must be deferred or processed before the next call. You may want to call BLManager::setState(SPECIAL) to prevent this track from being entered into NTuples for virtualdetectors (etc.) encountered during processing. This function MUST NOT BE CALLED if any track is being processed via any means (e.g. the normal BeamOn() is executing).

References BLAssert, stackManager, and trackManager.

Referenced by beamOnCollective(), BLCMDtracker::fitTrack(), BLCMDtracker::operator()(), and processAllDeferredTracksAndTheirSecondaries().

00413 {
00414         if(!trackManager || !stackManager) {
00415                 trackManager = eventManager->GetTrackingManager();
00416                 stackManager = eventManager->GetStackManager();
00417         }
00418         BLAssert(track != 0);
00419         trackManager->ProcessOneTrack(track);
00420 }

void BLRunManager::discardAllSecondaries (  ) 

discardAllSecondaries() will discard all secondaries.

References trackManager.

00423 {
00424         G4TrackVector *secondaries = trackManager->GimmeSecondaries();
00425         for(unsigned j=0; j<secondaries->size(); ++j)
00426                 delete (*secondaries)[j];
00427         secondaries->clear();
00428 }

int BLRunManager::deferAllSecondaries ( int  secondaryid = 10000,
int  parentid = -1 
)

deferAllSecondaries() will defer all secondaries until deferred tracks are processed. Returns # tracks deferred.

References deferOneTrack(), and trackManager.

Referenced by processAllDeferredTracksAndTheirSecondaries(), and processAllSecondariesAndDeferredTracks().

00431 {
00432         int n=0;
00433         G4TrackVector *secondaries = trackManager->GimmeSecondaries();
00434         for(unsigned j=0; j<secondaries->size(); ++j) {
00435                 G4Track *track = (*secondaries)[j];
00436                 track->SetTrackID(secondaryid++);
00437                 track->SetParentID(parentid); 
00438                 deferOneTrack(track);
00439                 ++n;
00440         }
00441         secondaries->clear();
00442 
00443         return n;
00444 }

void BLRunManager::deferOneTrack ( G4Track *  track  ) 

deferOneTrack() will defer a single track until deferred tracks are processed.

References stackManager.

Referenced by deferAllSecondaries().

00447 {
00448         stackManager->PushOneTrack(track);
00449 }

int BLRunManager::processAllSecondariesAndDeferredTracks ( int  secondaryid = 10000,
int  parentid = -1 
) [inline]

processAllSecondariesAndDeferredTracks() will defer all secondaries and then process all deferred tracks and their secondaries. returns # tracks processed.

References deferAllSecondaries(), and processAllDeferredTracksAndTheirSecondaries().

00301                                                                          {
00302                 int first = secondaryid;
00303                 secondaryid += deferAllSecondaries(secondaryid,parentid);
00304                 secondaryid +=
00305                        processAllDeferredTracksAndTheirSecondaries(secondaryid);
00306                 return secondaryid-first;
00307         }

int BLRunManager::processAllDeferredTracksAndTheirSecondaries ( int  trackid = 10000  ) 

processAllDeferredTracksAndTheirSecondaries() will process all deferred tracks, including processing their secondaries. Returns # tracks processed.

References deferAllSecondaries(), processOneTrack(), and stackManager.

Referenced by processAllSecondariesAndDeferredTracks().

00452 {
00453         int first=secondaryid;
00454         for(;;) {
00455                 G4VTrajectory *trajectory=0;
00456                 G4Track *track = stackManager->PopNextTrack(&trajectory);
00457                 if(!track) break;
00458                 processOneTrack(track);
00459                 int parentid = track->GetTrackID();
00460                 delete track;
00461                 secondaryid += deferAllSecondaries(secondaryid,parentid);
00462         }
00463 
00464         return secondaryid-first;
00465 }

G4Track * BLRunManager::popOneSecondary (  ) 

popOneSecondary() will return a pointer to one secondary, removing it from the list of secondaries; returns NULL if no more. Intended to be called in a loop immediately after a call to processOneTrack(). The returned pointer should eventually be deleted. Order is LIFO.

References trackManager.

00468 {
00469         G4TrackVector *secondaries = trackManager->GimmeSecondaries();
00470         int i = secondaries->size();
00471         if(i == 0) return 0;
00472         G4Track *track = secondaries->back();
00473         secondaries->pop_back();
00474         return track;
00475 }

G4Track * BLRunManager::popOneDeferredTrack (  ) 

popOneDeferredTrack() will return a pointer to one deferred track, removing it from the list of deferred tracks; returns NULL if no more. The returned pointer should eventually be deleted. Note that a loop calling this function can defer other tracks or secondaries, which simply extends the loop. Order is LIFO.

References stackManager.

00478 {
00479         G4VTrajectory *trajectory=0;
00480         return stackManager->PopNextTrack(&trajectory);
00481 }

void BLRunManager::RunInitialization (  )  [inline]

Referenced by beamOnCollective(), and beginRun().

00330 { G4RunManager::RunInitialization(); }

void BLRunManager::DoEventLoop ( G4int  n_event,
const char *  macroFile,
G4int  n_select 
) [inline]

00331                                                                              {
00332                 G4RunManager::DoEventLoop(n_event,macroFile,n_select);
00333         }

void BLRunManager::RunTermination (  )  [inline]

Referenced by beamOnCollective(), and endRun().

00334 { G4RunManager::RunTermination(); }

G4EventManager* BLRunManager::getEventManager (  )  [inline]

Referenced by BLCMDtrackermode::callback().

00335 { return eventManager; }

G4TrackingManager* BLRunManager::getTrackingManager (  )  [inline]

References trackManager.

00336 { return trackManager; }

G4StackManager* BLRunManager::getStackManager (  )  [inline]

References stackManager.

00337 { return stackManager; }

void BLRunManager::setCurrentEvent ( G4Event *  ev  )  [inline]

Referenced by BLCMDtrace::collectiveStep().

00338 { currentEvent = ev; }

G4bool BLRunManager::getRunAborted (  )  [inline]

Referenced by BLCMDtrackermode::callback().

00339 { return runAborted; }


Member Data Documentation

BLRunManager * BLRunManager::singleton = 0 [static, private]

Referenced by BLRunManager(), and getObject().

G4TrackingManager* BLRunManager::trackManager [private]

G4StackManager* BLRunManager::stackManager [private]

std::vector<BLTrackData> BLRunManager::trackVector [private]

G4double BLRunManager::stepTime [private]

G4double BLRunManager::deltaT [private]

bool BLRunManager::rejected [private]


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