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.
void DoEventLoop (G4int n_event, const char *macroFile, G4int n_select)
 DoEventLoop() is REPLACED, so abondonCurrentEvent() works.
void abandonCurrentEvent ()
 abandonCurrentEvent() will abandon the current event by performing a longjmp into DoEventLoop(). Not guaranteed to work, but probably does, even if called from BLAlarm::signal(). NOTE: If the current event can be abandoned, this function never returns; if the event cannot be abandoned, it returns immediately.
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 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
jmp_buf jmpBuf
bool validJmpBuf
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, stepTime, trackManager, and validJmpBuf.

Referenced by getObject().

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

BLRunManager::~BLRunManager (  )  [virtual]

Destructor.

References computeVector, and trackVector.

00131 {
00132         trackVector.clear();
00133         computeVector.clear();
00134 }


Member Function Documentation

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

References nextSecondaryTrackID, and trackVector.

Referenced by beamOnCollective().

00384 {
00385         G4int evId = currentEvent->GetEventID();
00386         if(track->GetTrackID() <= 0) {
00387                 track->SetTrackID(nextSecondaryTrackID++);
00388         }
00389 
00390         G4Event *tmp = currentEvent;
00391         currentEvent = new G4Event(*currentEvent);
00392 
00393         trackVector.push_back(BLTrackData(currentEvent,track));
00394 
00395         currentEvent = tmp;
00396 }

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 Bunch::Bunch(), BLCMDtrackermode::callback(), BLCMDtrace::command(), BLCMDspacechargelw::command(), BLCMDspacecharge::command(), BLCMDcollective::command(), BLCMDtracker::fitTrack(), BLManager::handleCallbacks(), BLCMDtracker::handlePreviousTracks(), BLCMDtracker::operator()(), BLCMDtracker::PostUserTrackingAction(), BLAlarm::sighandler(), TraceNTuple::TraceNTuple(), and TStepLimiter::TStepLimiter().

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

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().

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

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

DoEventLoop() is REPLACED, so abondonCurrentEvent() works.

References jmpBuf, and validJmpBuf.

00144 {
00145         if(verboseLevel>0) { timer->Start(); }
00146 
00147         G4String msg;
00148         if(macroFile!=0) { 
00149                 if(n_select<0) n_select = n_event;
00150                 msg = "/control/execute ";
00151                 msg += macroFile;
00152         } else {
00153                 n_select = -1;
00154         }
00155 
00156         int local_n_select = n_select; // argument can be clobbered by longjmp
00157 
00158         // Event loop
00159         volatile G4int i_event = 0;
00160         if(setjmp(jmpBuf) != 0) ++i_event;
00161         validJmpBuf = true;
00162         while(i_event<n_event) {
00163                 currentEvent = GenerateEvent(i_event);
00164                 eventManager->ProcessOneEvent(currentEvent);
00165                 AnalyzeEvent(currentEvent);
00166                 UpdateScoring();
00167                 if(i_event<local_n_select)
00168                         G4UImanager::GetUIpointer()->ApplyCommand(msg);
00169                 StackPreviousEvent(currentEvent);
00170                 currentEvent = 0;
00171                 if(runAborted) break;
00172                 ++i_event;
00173         }
00174         validJmpBuf = false;
00175 
00176         if(verboseLevel>0) {
00177                 timer->Stop();
00178                 G4cout << "Run terminated." << G4endl;
00179                 G4cout << "Run Summary" << G4endl;
00180                 if(runAborted) {
00181                         G4cout << "  Run Aborted after " << i_event + 1 <<
00182                                                 " events processed." << G4endl;
00183                 } else {
00184                         G4cout << "  Number of events processed : " << n_event
00185                                                                 << G4endl; 
00186                 }
00187                 G4cout << "  "  << *timer << G4endl;
00188         }
00189 }

void BLRunManager::abandonCurrentEvent (  ) 

abandonCurrentEvent() will abandon the current event by performing a longjmp into DoEventLoop(). Not guaranteed to work, but probably does, even if called from BLAlarm::signal(). NOTE: If the current event can be abandoned, this function never returns; if the event cannot be abandoned, it returns immediately.

References jmpBuf, and validJmpBuf.

Referenced by BLAlarm::sighandler().

00192 {
00193         if(validJmpBuf) {
00194                 G4Exception("BLRunManager","Event Time Limit Alarm",
00195                                                         EventMustBeAborted,"");
00196                 G4StateManager::GetStateManager()->
00197                                                 SetNewState(G4State_GeomClosed);
00198                 longjmp(jmpBuf,1);
00199         }
00200 }

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

callback() from BLCallback.

Reimplemented from BLCallback.

References collectiveMode.

Referenced by BLManager::handleCallbacks().

00203 {
00204         // Only perform actions between reference and Beam
00205         if(type != 1) return;
00206 
00207         // in collectiveMode, add a TStepLimiter to every particle
00208         if(collectiveMode) {
00209                 printf("BLRunManager: adding TStepLimiter processes to all particles\n");
00210                 G4ParticleTable::G4PTblDicIterator *theParticleIterator
00211                         = G4ParticleTable::GetParticleTable()->GetIterator();
00212                 theParticleIterator->reset();
00213                 while((*theParticleIterator)()) {
00214                         G4ParticleDefinition *pd = 
00215                                         theParticleIterator->value();
00216                         if(pd->IsShortLived()) continue;
00217                         G4ProcessManager *pmgr = pd->GetProcessManager();
00218                         if(!pmgr) continue;
00219                         TStepLimiter *tsl = new TStepLimiter();
00220                         pmgr->AddProcess(tsl,-1,-1,3);
00221                         pmgr->SetProcessOrdering(tsl,idxPostStep);
00222                         pmgr->SetProcessOrdering(tsl,idxAtRest);
00223                 }
00224         }
00225 }

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().

00228 {
00229         BLManager *manager = BLManager::getObject();
00230         G4TransportationManager *transportManager =
00231                         G4TransportationManager::GetTransportationManager();
00232         G4StateManager* stateManager = G4StateManager::GetStateManager();
00233         numberOfEventToBeProcessed = nEvents;
00234 
00235         if(!ConfirmBeamOnCondition()) {
00236                 G4Exception("BLRunManager","Cannot run beam",FatalException,"");
00237         }
00238         printf("================== In Collective Mode ==================\n");
00239         ConstructScoringWorlds();
00240         RunInitialization();
00241         G4ApplicationState currentState = stateManager->GetCurrentState();
00242         if(currentState != G4State_GeomClosed) {
00243                 G4Exception("BLRunManager","Geometry not closed",FatalException,
00244                                                                         "");
00245         }
00246 
00247         manager->setState(BEAM);
00248         stateManager->SetNewState(G4State_EventProc);
00249         trackManager = 0; // unused in collective mode
00250         stackManager = eventManager->GetStackManager();
00251         runAborted = false;
00252 
00253         // the per-event time limit makes no sense in collective mode
00254         manager->setEventTimeLimit(-1); // infinite
00255 
00256         // event loop to create trackVector
00257         G4PrimaryTransformer *transformer = new G4PrimaryTransformer();
00258         G4Navigator* navigator = transportManager->GetNavigatorForTracking();
00259         stepTime = -DBL_MAX;
00260         int nev;
00261         for(nev=0; nev<nEvents; ++nev) {
00262                 currentEvent = GenerateEvent(nev);
00263                 if(!currentEvent || runAborted) break;
00264                 stackManager->PrepareNewEvent();
00265                 manager->BeginOfEventAction(currentEvent);
00266                 G4TrackVector *tv=transformer->GimmePrimaries(currentEvent,0);
00267                 int evId=currentEvent->GetEventID();
00268                 if(evId <= 0) {
00269                         evId = nev;
00270                         currentEvent->SetEventID(evId);
00271                 }
00272                 nextSecondaryTrackID = manager->getNextSecondaryTrackID();
00273                 for(unsigned j=0; j<tv->size(); ++j) {
00274                         G4Track *track = (*tv)[j];
00275                         BLAssert(j==0); // cannot handle multiple primaries 
00276                         track->SetTrackID(manager->getPrimaryTrackID());
00277                         appendTrack(track);
00278                         G4double t = track->GetGlobalTime();
00279                         if(t > stepTime) stepTime = t;
00280                 }
00281                 if(tv->size() == 0) delete currentEvent;
00282                 tv->clear();
00283                 if(runAborted) break;
00284                 transportManager->SetNavigatorForTracking(navigator);
00285                 manager->setEventID(evId+1);
00286         }
00287 
00288         printf("=========== Collective: %d Events, %ld Tracks ==============\n",
00289                                 nev,(long)trackVector.size());
00290 
00291         for(unsigned i=0; i<computeVector.size(); ++i)
00292                 computeVector[i]->beginCollectiveTracking(trackVector);
00293 
00294         // Step loop
00295         runAborted = false;
00296         int nActive;
00297         bool first=true;
00298         G4double saveStepTime = stepTime;
00299         std::vector<G4Track*> saveVector;
00300         do {
00301                 if(runAborted) goto RunAborted;
00302                 nActive = 0;
00303                 rejected = false;
00304                 // note: secondaries can be appended to trackVector during loop
00305                 for(unsigned i=0; i<trackVector.size(); ++i) {
00306                         if(rejected && !first) break;
00307                         BLTrackData *td = &trackVector[i];
00308                         G4Track *track = td->track;
00309                         G4TrackStatus trackStatus = track->GetTrackStatus();
00310 //printf("Collective Track loop: ev=%d trk=%d status=%d fAlive=%d fSuspend=%d\n",td->event->GetEventID(),track->GetTrackID(),(int)trackStatus, fAlive,fSuspend);
00311                         if(trackStatus != fAlive && 
00312                                                 trackStatus != fStopButAlive)
00313                                 continue;
00314                         ++nActive;
00315                         currentEvent = td->event;
00316                         manager->setEventID(td->event->GetEventID());
00317                         currentTrackIndex = i;
00318                         for(unsigned j=0; j<computeVector.size(); ++j)
00319                                 computeVector[j]->beginTrack(trackVector,i);
00320                         if(runAborted) goto RunAborted;
00321                         processOneTrack(track);
00322                         currentTrackIndex = -1;
00323                         if(td->event->IsAborted())
00324                                 track->SetTrackStatus(fKillTrackAndSecondaries);
00325                         if(runAborted) goto RunAborted;
00326                         trackStatus = track->GetTrackStatus();
00327                         if(trackStatus == fSuspend)
00328                                 track->SetTrackStatus(TStepLimiter::getPrevStatus());
00329                         G4TrackVector *tv = trackManager->GimmeSecondaries();
00330                         if(keepSecondaries && trackStatus != fKillTrackAndSecondaries) {
00331                                 for(unsigned j=0; j<tv->size(); ++j)
00332                                         appendTrack((*tv)[j]);
00333                                 tv->clear();
00334                         } else {
00335                                 for(unsigned j=0; j<tv->size(); ++j)
00336                                         delete (*tv)[j];
00337                                 tv->clear();
00338                         }
00339                 }
00340                 if(runAborted) goto RunAborted;
00341                 // don't call collectiveStep() if rejected during tracking
00342                 if(!rejected || first) {
00343                         for(unsigned i=0; i<computeVector.size(); ++i)
00344                                 computeVector[i]->collectiveStep(trackVector);
00345                 }
00346                 if(rejected && !first) {
00347                         // discard current tracks, and restore saved ones
00348                         unsigned n=saveVector.size();
00349                         // first, discard any added secondaries
00350                         while(trackVector.size() > n) {
00351                                 delete trackVector.back().track;
00352                                 delete trackVector.back().event;
00353                                 trackVector.pop_back();
00354                         }
00355                         BLAssert(trackVector.size() == saveVector.size());
00356                         for(unsigned i=0; i<n; ++i) {
00357                                 delete trackVector[i].track;
00358                                 trackVector[i].track = saveVector[i];
00359                         }
00360                         stepTime = saveStepTime;
00361                 } else {
00362                         // delete saved tracks, and save current ones
00363                         for(unsigned i=0; i<saveVector.size(); ++i)
00364                                 delete saveVector[i];
00365                         saveVector.clear();
00366                         for(unsigned i=0; i<trackVector.size(); ++i)
00367                                 saveVector.push_back(
00368                                         new G4Track(*trackVector[i].track));
00369                         saveStepTime = stepTime;
00370                         stepTime += deltaT;
00371                 }
00372                 first = false;
00373         } while(nActive > 0);
00374 
00375 RunAborted:
00376         for(unsigned i=0; i<computeVector.size(); ++i)
00377                 computeVector[i]->endCollectiveTracking(trackVector);
00378 
00379         RunTermination();
00380         manager->setState(IDLE);
00381 }

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(), BLCMDspacechargelw::command(), BLCMDspacecharge::command(), and BLCMDcollective::command().

00192                 { 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.

00200 { 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.

00205 { 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 BLCMDspacechargelw::command(), BLCMDspacecharge::command(), BLCMDcollective::command(), and BLManager::trackTuneAndReferenceParticles().

00219                 { collectiveMode=flag;
00220                   keepSecondaries=_keepSecondaries;
00221                   setDeltaT(_deltaT);
00222                 }

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 Bunch::addFieldValue(), TStepLimiter::AtRestGetPhysicalInteractionLength(), BLCMDcollective::collectiveStep(), Bunch::computeField(), Bunch::init(), TStepLimiter::PostStepGetPhysicalInteractionLength(), and Bunch::resizeGrid().

00228 { 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.

00234 { stepTime = t; }

G4double BLRunManager::getDeltaT (  )  [inline]

getDeltaT() will get the time interval between time steps

References deltaT.

Referenced by Bunch::addFieldValue(), Bunch::computeField(), and Bunch::resizeGrid().

00237 { 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 BLCMDspacechargelw::beginCollectiveTracking(), BLCMDspacecharge::beginCollectiveTracking(), BLCMDspacechargelw::command(), BLCMDspacecharge::command(), BLCMDcollective::command(), and setCollectiveMode().

00244 { 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.

00258 { 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().

00399 {
00400         BLManager::getObject()->setState(BEAM);
00401         ConstructScoringWorlds();
00402         RunInitialization();
00403         SetRunIDCounter(runid);
00404         runAborted = false;
00405 }

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().

00408 {
00409         RunTermination();
00410         BLManager::getObject()->setState(IDLE);
00411 }

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().

00414 {
00415         if(currentEvent) delete currentEvent;
00416         currentEvent = new G4Event(evid);
00417         trackManager = eventManager->GetTrackingManager();
00418         stackManager = eventManager->GetStackManager();
00419         BLManager::getObject()->BeginOfEventAction(currentEvent);
00420 }

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().

00423 {
00424         BLManager::getObject()->EndOfEventAction(currentEvent);
00425         if(currentEvent) delete currentEvent;
00426         currentEvent = 0;
00427         trackManager = 0;
00428         stackManager = 0;
00429 }

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.

00432 {
00433         static bool more=true;
00434         static int nev=0;
00435         static BLManager *manager=0;
00436         static G4PrimaryTransformer *transformer=0;
00437         if(!manager) {
00438                 manager = BLManager::getObject();
00439                 transformer = new G4PrimaryTransformer();
00440         }
00441 
00442         if(!more) return false;
00443 
00444         manager->setState(BEAM);
00445         currentEvent = GenerateEvent(nev);
00446         if(!currentEvent || runAborted) {
00447                 more = false;
00448                 return false;
00449         }
00450         *pevent = currentEvent;
00451         G4TrackVector *tv=transformer->GimmePrimaries(currentEvent,0);
00452         int evId=currentEvent->GetEventID();
00453         if(evId <= 0) {
00454                 evId = nev;
00455                 currentEvent->SetEventID(evId);
00456         }
00457         for(unsigned j=0; j<tv->size(); ++j) {
00458                 G4Track *track = (*tv)[j];
00459                 BLAssert(j==0); // cannot handle multiple primaries 
00460                 track->SetTrackID(manager->getPrimaryTrackID());
00461                 track->SetParentID(manager->getPrimaryParentID());
00462                 *ptrack = track;
00463         }
00464         tv->clear();
00465         if(runAborted) {
00466                 more = false;
00467                 return false;
00468         }
00469         ++nev;
00470         manager->setState(SPECIAL);
00471         return true;
00472 }

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(), BLCMDspacechargelw::beginTrack(), BLCMDtracker::fitTrack(), BLCMDtracker::operator()(), and processAllDeferredTracksAndTheirSecondaries().

00475 {
00476         if(!trackManager || !stackManager) {
00477                 trackManager = eventManager->GetTrackingManager();
00478                 stackManager = eventManager->GetStackManager();
00479         }
00480         BLAssert(track != 0);
00481         trackManager->ProcessOneTrack(track);
00482 }

void BLRunManager::discardAllSecondaries (  ) 

discardAllSecondaries() will discard all secondaries.

References trackManager.

00485 {
00486         G4TrackVector *secondaries = trackManager->GimmeSecondaries();
00487         for(unsigned j=0; j<secondaries->size(); ++j)
00488                 delete (*secondaries)[j];
00489         secondaries->clear();
00490 }

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().

00493 {
00494         int n=0;
00495         G4TrackVector *secondaries = trackManager->GimmeSecondaries();
00496         for(unsigned j=0; j<secondaries->size(); ++j) {
00497                 G4Track *track = (*secondaries)[j];
00498                 track->SetTrackID(secondaryid++);
00499                 track->SetParentID(parentid); 
00500                 deferOneTrack(track);
00501                 ++n;
00502         }
00503         secondaries->clear();
00504 
00505         return n;
00506 }

void BLRunManager::deferOneTrack ( G4Track *  track  ) 

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

References stackManager.

Referenced by deferAllSecondaries().

00509 {
00510         stackManager->PushOneTrack(track);
00511 }

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().

00314                                                                          {
00315                 int first = secondaryid;
00316                 secondaryid += deferAllSecondaries(secondaryid,parentid);
00317                 secondaryid +=
00318                        processAllDeferredTracksAndTheirSecondaries(secondaryid);
00319                 return secondaryid-first;
00320         }

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().

00514 {
00515         int first=secondaryid;
00516         for(;;) {
00517                 G4VTrajectory *trajectory=0;
00518                 G4Track *track = stackManager->PopNextTrack(&trajectory);
00519                 if(!track) break;
00520                 processOneTrack(track);
00521                 int parentid = track->GetTrackID();
00522                 delete track;
00523                 secondaryid += deferAllSecondaries(secondaryid,parentid);
00524         }
00525 
00526         return secondaryid-first;
00527 }

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.

00530 {
00531         G4TrackVector *secondaries = trackManager->GimmeSecondaries();
00532         int i = secondaries->size();
00533         if(i == 0) return 0;
00534         G4Track *track = secondaries->back();
00535         secondaries->pop_back();
00536         return track;
00537 }

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.

00540 {
00541         G4VTrajectory *trajectory=0;
00542         return stackManager->PopNextTrack(&trajectory);
00543 }

void BLRunManager::RunInitialization (  )  [inline]

Referenced by beamOnCollective(), and beginRun().

00343 { G4RunManager::RunInitialization(); }

void BLRunManager::RunTermination (  )  [inline]

Referenced by beamOnCollective(), and endRun().

00344 { G4RunManager::RunTermination(); }

G4EventManager* BLRunManager::getEventManager (  )  [inline]

Referenced by BLCMDtrackermode::callback().

00345 { return eventManager; }

G4TrackingManager* BLRunManager::getTrackingManager (  )  [inline]

References trackManager.

00346 { return trackManager; }

G4StackManager* BLRunManager::getStackManager (  )  [inline]

References stackManager.

00347 { return stackManager; }

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

Referenced by BLCMDtrace::collectiveStep().

00348 { currentEvent = ev; }

G4bool BLRunManager::getRunAborted (  )  [inline]

Referenced by BLCMDtrackermode::callback().

00349 { return runAborted; }


Member Data Documentation

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

Referenced by BLRunManager(), and getObject().

G4TrackingManager* BLRunManager::trackManager [private]

G4StackManager* BLRunManager::stackManager [private]

jmp_buf BLRunManager::jmpBuf [private]

Referenced by abandonCurrentEvent(), and DoEventLoop().

bool BLRunManager::validJmpBuf [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