#include <BLRunManager.hh>
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 | getNextStepTime () |
getNextStepTime() will get the time for the next step in time NOTE: inside collectiveStep() this function returns the time of the CURRENT step, not the next one. | |
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 BLRunManager * | getObject () |
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< BLTrackData > | trackVector |
std::vector < BLCollectiveComputation * > | computeVector |
bool | collectiveMode |
bool | keepSecondaries |
G4double | nextStepTime |
G4double | deltaT |
G4int | currentTrackIndex |
bool | rejected |
int | nextSecondaryTrackID |
Static Private Attributes | |
static BLRunManager * | singleton = 0 |
BLRunManager::BLRunManager | ( | ) |
Constructor.
References collectiveMode, currentTrackIndex, deltaT, keepSecondaries, nextSecondaryTrackID, nextStepTime, rejected, singleton, and stackManager.
Referenced by getObject().
00111 : G4RunManager(), trackVector(), computeVector() 00112 { 00113 assert(singleton==0); 00114 singleton = this; 00115 00116 stackManager = 0; 00117 collectiveMode = false; 00118 keepSecondaries = true; 00119 nextStepTime = 0.0; 00120 deltaT = 1.0*ns; 00121 currentTrackIndex = -1; 00122 rejected = false; 00123 nextSecondaryTrackID = -9999; 00124 }
BLRunManager::~BLRunManager | ( | ) | [virtual] |
Destructor.
References computeVector, and trackVector.
00127 { 00128 trackVector.clear(); 00129 computeVector.clear(); 00130 }
void BLRunManager::appendTrack | ( | G4Track * | track | ) | [private] |
References nextSecondaryTrackID, and trackVector.
Referenced by beamOnCollective().
00321 { 00322 G4int evId = currentEvent->GetEventID(); 00323 if(track->GetTrackID() <= 0) { 00324 track->SetTrackID(nextSecondaryTrackID++); 00325 } 00326 00327 G4Event *tmp = currentEvent; 00328 currentEvent = new G4Event(*currentEvent); 00329 00330 trackVector.push_back(BLTrackData(currentEvent,track)); 00331 00332 currentEvent = tmp; 00333 }
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().
00100 { 00101 if(!singleton) { 00102 new BLRunManager(); 00103 // startup sequencing prevents registering callback, so it 00104 // is hard-coded in BLManager::handleCallbacks(). This is due 00105 // to the BLManager constructor calling this function. 00106 } 00107 00108 return singleton; 00109 }
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().
00132 { 00133 if(collectiveMode && BLManager::getObject()->getState() == BEAM) 00134 beamOnCollective(nEvents); 00135 else 00136 G4RunManager::BeamOn(nEvents,macroFile,n_select); 00137 }
void BLRunManager::callback | ( | int | type | ) | [virtual] |
callback() from BLCallback.
Reimplemented from BLCallback.
References collectiveMode.
Referenced by BLManager::handleCallbacks().
00140 { 00141 // Only perform actions between reference and Beam 00142 if(type != 1) return; 00143 00144 // in collectiveMode, add a TStepLimiter to every particle 00145 if(collectiveMode) { 00146 printf("BLRunManager: adding TStepLimiter processes to all particles\n"); 00147 G4ParticleTable::G4PTblDicIterator *theParticleIterator 00148 = G4ParticleTable::GetParticleTable()->GetIterator(); 00149 theParticleIterator->reset(); 00150 while((*theParticleIterator)()) { 00151 G4ParticleDefinition *pd = 00152 theParticleIterator->value(); 00153 if(pd->IsShortLived()) continue; 00154 G4ProcessManager *pmgr = pd->GetProcessManager(); 00155 if(!pmgr) continue; 00156 TStepLimiter *tsl = new TStepLimiter(); 00157 pmgr->AddProcess(tsl,-1,-1,3); 00158 pmgr->SetProcessOrdering(tsl,idxPostStep); 00159 pmgr->SetProcessOrdering(tsl,idxAtRest); 00160 } 00161 } 00162 }
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(), computeVector, currentTrackIndex, deltaT, BLTrackData::event, BLManager::getNextSecondaryTrackID(), BLManager::getObject(), TStepLimiter::getPrevStatus(), BLManager::getPrimaryTrackID(), IDLE, keepSecondaries, nextSecondaryTrackID, nextStepTime, processOneTrack(), rejected, RunInitialization(), RunTermination(), BLManager::setEventID(), BLManager::setEventTimeLimit(), BLManager::setState(), stackManager, BLTrackData::track, trackManager, and trackVector.
Referenced by BeamOn().
00165 { 00166 BLManager *manager = BLManager::getObject(); 00167 G4TransportationManager *transportManager = 00168 G4TransportationManager::GetTransportationManager(); 00169 G4StateManager* stateManager = G4StateManager::GetStateManager(); 00170 numberOfEventToBeProcessed = nEvents; 00171 00172 if(!ConfirmBeamOnCondition()) { 00173 G4Exception("BLRunManager","Cannot run beam",FatalException,""); 00174 } 00175 printf("================== In Collective Mode ==================\n"); 00176 ConstructScoringWorlds(); 00177 RunInitialization(); 00178 G4ApplicationState currentState = stateManager->GetCurrentState(); 00179 if(currentState != G4State_GeomClosed) { 00180 G4Exception("BLRunManager","Geometry not closed",FatalException, 00181 ""); 00182 } 00183 00184 manager->setState(BEAM); 00185 stateManager->SetNewState(G4State_EventProc); 00186 trackManager = 0; // unused in collective mode 00187 stackManager = eventManager->GetStackManager(); 00188 runAborted = false; 00189 00190 // the per-event time limit makes no sense in collective mode 00191 manager->setEventTimeLimit(-1); // infinite 00192 00193 // event loop to create trackVector 00194 G4PrimaryTransformer *transformer = new G4PrimaryTransformer(); 00195 G4Navigator* navigator = transportManager->GetNavigatorForTracking(); 00196 nextStepTime = -DBL_MAX; 00197 int nev; 00198 for(nev=0; nev<nEvents; ++nev) { 00199 currentEvent = GenerateEvent(nev); 00200 if(!currentEvent || runAborted) break; 00201 stackManager->PrepareNewEvent(); 00202 manager->BeginOfEventAction(currentEvent); 00203 G4TrackVector *tv=transformer->GimmePrimaries(currentEvent,0); 00204 int evId=currentEvent->GetEventID(); 00205 if(evId <= 0) { 00206 evId = nev; 00207 currentEvent->SetEventID(evId); 00208 } 00209 nextSecondaryTrackID = manager->getNextSecondaryTrackID(); 00210 for(unsigned j=0; j<tv->size(); ++j) { 00211 G4Track *track = (*tv)[j]; 00212 assert(j==0); // cannot handle multiple primaries 00213 track->SetTrackID(manager->getPrimaryTrackID()); 00214 appendTrack(track); 00215 G4double t = track->GetGlobalTime(); 00216 if(t > nextStepTime) nextStepTime = t; 00217 } 00218 if(tv->size() == 0) delete currentEvent; 00219 tv->clear(); 00220 if(runAborted) break; 00221 transportManager->SetNavigatorForTracking(navigator); 00222 manager->setEventID(evId+1); 00223 } 00224 00225 printf("=========== Collective: %d Events, %ld Tracks ==============\n", 00226 nev,(long)trackVector.size()); 00227 00228 for(unsigned i=0; i<computeVector.size(); ++i) 00229 computeVector[i]->beginCollectiveTracking(trackVector); 00230 00231 // Step loop 00232 runAborted = false; 00233 int nActive; 00234 bool first=true; 00235 G4double saveNextStepTime = nextStepTime; 00236 std::vector<G4Track*> saveVector; 00237 do { 00238 if(runAborted) goto RunAborted; 00239 nActive = 0; 00240 rejected = false; 00241 // note: secondaries can be appended to trackVector during loop 00242 for(unsigned i=0; i<trackVector.size(); ++i) { 00243 if(rejected && !first) break; 00244 BLTrackData *td = &trackVector[i]; 00245 G4Track *track = td->track; 00246 G4TrackStatus trackStatus = track->GetTrackStatus(); 00247 //printf("Collective Track loop: ev=%d trk=%d status=%d fAlive=%d fSuspend=%d\n",td->event->GetEventID(),track->GetTrackID(),(int)trackStatus, fAlive,fSuspend); 00248 if(trackStatus != fAlive && 00249 trackStatus != fStopButAlive) 00250 continue; 00251 ++nActive; 00252 currentEvent = td->event; 00253 manager->setEventID(td->event->GetEventID()); 00254 currentTrackIndex = i; 00255 for(unsigned j=0; j<computeVector.size(); ++j) 00256 computeVector[j]->beginTrack(trackVector,i); 00257 if(runAborted) goto RunAborted; 00258 processOneTrack(track); 00259 currentTrackIndex = -1; 00260 if(td->event->IsAborted()) 00261 track->SetTrackStatus(fKillTrackAndSecondaries); 00262 if(runAborted) goto RunAborted; 00263 trackStatus = track->GetTrackStatus(); 00264 if(trackStatus == fSuspend) 00265 track->SetTrackStatus(TStepLimiter::getPrevStatus()); 00266 G4TrackVector *tv = trackManager->GimmeSecondaries(); 00267 if(keepSecondaries && trackStatus != fKillTrackAndSecondaries) { 00268 for(unsigned j=0; j<tv->size(); ++j) 00269 appendTrack((*tv)[j]); 00270 tv->clear(); 00271 } else { 00272 for(unsigned j=0; j<tv->size(); ++j) 00273 delete (*tv)[j]; 00274 tv->clear(); 00275 } 00276 } 00277 if(runAborted) goto RunAborted; 00278 // don't call collectiveStep() if rejected during tracking 00279 if(!rejected || first) { 00280 for(unsigned i=0; i<computeVector.size(); ++i) 00281 computeVector[i]->collectiveStep(trackVector); 00282 } 00283 if(rejected && !first) { 00284 // discard current tracks, and restore saved ones 00285 unsigned n=saveVector.size(); 00286 // first, discard any added secondaries 00287 while(trackVector.size() > n) { 00288 delete trackVector.back().track; 00289 delete trackVector.back().event; 00290 trackVector.pop_back(); 00291 } 00292 assert(trackVector.size() == saveVector.size()); 00293 for(unsigned i=0; i<n; ++i) { 00294 delete trackVector[i].track; 00295 trackVector[i].track = saveVector[i]; 00296 } 00297 nextStepTime = saveNextStepTime; 00298 } else { 00299 // delete saved tracks, and save current ones 00300 for(unsigned i=0; i<saveVector.size(); ++i) 00301 delete saveVector[i]; 00302 saveVector.clear(); 00303 for(unsigned i=0; i<trackVector.size(); ++i) 00304 saveVector.push_back( 00305 new G4Track(*trackVector[i].track)); 00306 saveNextStepTime = nextStepTime; 00307 nextStepTime += deltaT; 00308 } 00309 first = false; 00310 } while(nActive > 0); 00311 00312 RunAborted: 00313 for(unsigned i=0; i<computeVector.size(); ++i) 00314 computeVector[i]->endCollectiveTracking(trackVector); 00315 00316 RunTermination(); 00317 manager->setState(IDLE); 00318 }
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().
00178 { 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.
00186 { 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.
00191 { return currentTrackIndex; }
bool BLRunManager::getCollectiveMode | ( | ) | [inline] |
getCollectiveMode() returns true if in collective mode.
References collectiveMode.
Referenced by BLManager::BeginOfEventAction(), BLManager::getExternalParentID(), BLManager::getExternalTrackID(), BLManager::incrEventsProcessed(), BLManager::PostUserTrackingAction(), BLCMDtrace::PostUserTrackingAction(), BLManager::PreUserTrackingAction(), BLCMDtrace::PreUserTrackingAction(), BLManager::trackTuneAndReferenceParticles(), and BLCMDtrace::UserSteppingAction().
00194 { return collectiveMode; }
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().
00205 { collectiveMode=flag; 00206 keepSecondaries=_keepSecondaries; 00207 setDeltaT(_deltaT); 00208 }
G4double BLRunManager::getNextStepTime | ( | ) | [inline] |
getNextStepTime() will get the time for the next step in time NOTE: inside collectiveStep() this function returns the time of the CURRENT step, not the next one.
References nextStepTime.
Referenced by TStepLimiter::AtRestGetPhysicalInteractionLength(), BLCMDcollective::collectiveStep(), and TStepLimiter::PostStepGetPhysicalInteractionLength().
00213 { return nextStepTime; }
G4double BLRunManager::getDeltaT | ( | ) | [inline] |
getDeltaT() will get the time interval between time steps
References deltaT.
00220 { 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 deltaT.
Referenced by BLCMDcollective::beginCollectiveTracking(), BLCMDcollective::command(), and setCollectiveMode().
00227 { assert(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.
00241 { 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().
00336 { 00337 BLManager::getObject()->setState(BEAM); 00338 ConstructScoringWorlds(); 00339 RunInitialization(); 00340 SetRunIDCounter(runid); 00341 runAborted = false; 00342 }
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().
00345 { 00346 RunTermination(); 00347 BLManager::getObject()->setState(IDLE); 00348 }
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().
00351 { 00352 if(currentEvent) delete currentEvent; 00353 currentEvent = new G4Event(evid); 00354 trackManager = eventManager->GetTrackingManager(); 00355 stackManager = eventManager->GetStackManager(); 00356 BLManager::getObject()->BeginOfEventAction(currentEvent); 00357 }
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().
00360 { 00361 BLManager::getObject()->EndOfEventAction(currentEvent); 00362 if(currentEvent) delete currentEvent; 00363 currentEvent = 0; 00364 trackManager = 0; 00365 stackManager = 0; 00366 }
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, BLManager::getObject(), BLManager::getPrimaryParentID(), BLManager::getPrimaryTrackID(), BLManager::setState(), and SPECIAL.
00369 { 00370 static bool more=true; 00371 static int nev=0; 00372 static BLManager *manager=0; 00373 static G4PrimaryTransformer *transformer=0; 00374 if(!manager) { 00375 manager = BLManager::getObject(); 00376 transformer = new G4PrimaryTransformer(); 00377 } 00378 00379 if(!more) return false; 00380 00381 manager->setState(BEAM); 00382 currentEvent = GenerateEvent(nev); 00383 if(!currentEvent || runAborted) { 00384 more = false; 00385 return false; 00386 } 00387 *pevent = currentEvent; 00388 G4TrackVector *tv=transformer->GimmePrimaries(currentEvent,0); 00389 int evId=currentEvent->GetEventID(); 00390 if(evId <= 0) { 00391 evId = nev; 00392 currentEvent->SetEventID(evId); 00393 } 00394 for(unsigned j=0; j<tv->size(); ++j) { 00395 G4Track *track = (*tv)[j]; 00396 assert(j==0); // cannot handle multiple primaries 00397 track->SetTrackID(manager->getPrimaryTrackID()); 00398 track->SetParentID(manager->getPrimaryParentID()); 00399 *ptrack = track; 00400 } 00401 tv->clear(); 00402 if(runAborted) { 00403 more = false; 00404 return false; 00405 } 00406 ++nev; 00407 manager->setState(SPECIAL); 00408 return true; 00409 }
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 stackManager, and trackManager.
Referenced by beamOnCollective(), BLCMDtracker::fitTrack(), BLCMDtracker::operator()(), and processAllDeferredTracksAndTheirSecondaries().
00412 { 00413 if(!trackManager || !stackManager) { 00414 trackManager = eventManager->GetTrackingManager(); 00415 stackManager = eventManager->GetStackManager(); 00416 } 00417 assert(track != 0); 00418 trackManager->ProcessOneTrack(track); 00419 }
void BLRunManager::discardAllSecondaries | ( | ) |
discardAllSecondaries() will discard all secondaries.
References trackManager.
00422 { 00423 G4TrackVector *secondaries = trackManager->GimmeSecondaries(); 00424 for(unsigned j=0; j<secondaries->size(); ++j) 00425 delete (*secondaries)[j]; 00426 secondaries->clear(); 00427 }
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().
00430 { 00431 int n=0; 00432 G4TrackVector *secondaries = trackManager->GimmeSecondaries(); 00433 for(unsigned j=0; j<secondaries->size(); ++j) { 00434 G4Track *track = (*secondaries)[j]; 00435 track->SetTrackID(secondaryid++); 00436 track->SetParentID(parentid); 00437 deferOneTrack(track); 00438 ++n; 00439 } 00440 secondaries->clear(); 00441 00442 return n; 00443 }
void BLRunManager::deferOneTrack | ( | G4Track * | track | ) |
deferOneTrack() will defer a single track until deferred tracks are processed.
References stackManager.
Referenced by deferAllSecondaries().
00446 { 00447 stackManager->PushOneTrack(track); 00448 }
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().
00297 { 00298 int first = secondaryid; 00299 secondaryid += deferAllSecondaries(secondaryid,parentid); 00300 secondaryid += 00301 processAllDeferredTracksAndTheirSecondaries(secondaryid); 00302 return secondaryid-first; 00303 }
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().
00451 { 00452 int first=secondaryid; 00453 for(;;) { 00454 G4VTrajectory *trajectory=0; 00455 G4Track *track = stackManager->PopNextTrack(&trajectory); 00456 if(!track) break; 00457 processOneTrack(track); 00458 int parentid = track->GetTrackID(); 00459 delete track; 00460 secondaryid += deferAllSecondaries(secondaryid,parentid); 00461 } 00462 00463 return secondaryid-first; 00464 }
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.
00467 { 00468 G4TrackVector *secondaries = trackManager->GimmeSecondaries(); 00469 int i = secondaries->size(); 00470 if(i == 0) return 0; 00471 G4Track *track = secondaries->back(); 00472 secondaries->pop_back(); 00473 return track; 00474 }
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.
00477 { 00478 G4VTrajectory *trajectory=0; 00479 return stackManager->PopNextTrack(&trajectory); 00480 }
void BLRunManager::RunInitialization | ( | ) | [inline] |
void BLRunManager::DoEventLoop | ( | G4int | n_event, | |
const char * | macroFile, | |||
G4int | n_select | |||
) | [inline] |
void BLRunManager::RunTermination | ( | ) | [inline] |
G4EventManager* BLRunManager::getEventManager | ( | ) | [inline] |
G4TrackingManager* BLRunManager::getTrackingManager | ( | ) | [inline] |
G4StackManager* BLRunManager::getStackManager | ( | ) | [inline] |
void BLRunManager::setCurrentEvent | ( | G4Event * | ev | ) | [inline] |
G4bool BLRunManager::getRunAborted | ( | ) | [inline] |
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] |
Referenced by appendTrack(), beamOnCollective(), getTrackVector(), and ~BLRunManager().
std::vector<BLCollectiveComputation *> BLRunManager::computeVector [private] |
Referenced by beamOnCollective(), registerCollectiveComputation(), and ~BLRunManager().
bool BLRunManager::collectiveMode [private] |
Referenced by BeamOn(), BLRunManager(), callback(), getCollectiveMode(), and setCollectiveMode().
bool BLRunManager::keepSecondaries [private] |
Referenced by beamOnCollective(), BLRunManager(), and setCollectiveMode().
G4double BLRunManager::nextStepTime [private] |
Referenced by beamOnCollective(), BLRunManager(), and getNextStepTime().
G4double BLRunManager::deltaT [private] |
Referenced by beamOnCollective(), BLRunManager(), getDeltaT(), and setDeltaT().
G4int BLRunManager::currentTrackIndex [private] |
Referenced by beamOnCollective(), BLRunManager(), and getCurrentTrackIndex().
bool BLRunManager::rejected [private] |
Referenced by beamOnCollective(), BLRunManager(), and rejectCollectiveStep().
int BLRunManager::nextSecondaryTrackID [private] |
Referenced by appendTrack(), beamOnCollective(), and BLRunManager().