00001 // BLRunManager.hh 00002 /* 00003 This source file is part of G4beamline, http://g4beamline.muonsinc.com 00004 Copyright (C) 2003,2004,2005,2006 by Tom Roberts, all rights reserved. 00005 00006 This program is free software; you can redistribute it and/or 00007 modify it under the terms of the GNU General Public License 00008 as published by the Free Software Foundation; either version 2 00009 of the License, or (at your option) any later version. 00010 00011 This program is distributed in the hope that it will be useful, 00012 but WITHOUT ANY WARRANTY; without even the implied warranty of 00013 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00014 GNU General Public License for more details. 00015 00016 http://www.gnu.org/copyleft/gpl.html 00017 */ 00018 00019 #ifndef BLRUNMANAGER_H 00020 #define BLRUNMANAGER_H 00021 00022 #include <vector> 00023 #include <map> 00024 00025 #include "G4RunManager.hh" 00026 #include "G4EventManager.hh" 00027 #include "G4TrackingManager.hh" 00028 #include "G4StackManager.hh" 00029 #include "G4SteppingManager.hh" 00030 #include "G4Navigator.hh" 00031 #include "G4PropagatorInField.hh" 00032 00033 #include "BLAssert.hh" 00034 #include "BLCollectiveComputation.hh" 00035 #include "BLCallback.hh" 00036 00037 00038 /** class BLRunManager -- enhanced version of G4RunManager. 00039 * 00040 * There are three basic modes of BLRunManager: normal, collective, 00041 * and enhanced. 00042 * 00043 * Normal mode: 00044 * This class inherents from G4RunManager, and can thus be used 00045 * just like it, for standard control of the Geant4 simulation 00046 * (via function BeamOn()). 00047 * 00048 * Collective mode: 00049 * 00050 * In collective mode a vector<G4TrackData*> is generated from the 00051 * input file's beam commands, and then all tracks are stepped for 00052 * one time step after which the BLCollectiveComputation::collectiveStep() 00053 * is called; this is continued until all tracks are killed. Any 00054 * secondaries generated can either be kept or discarded (see 00055 * setCollectiveMode() below). The TrackID-s of primary tracks are 00056 * preserved, and TrackIDs for secondary tracks are assigned sequentially 00057 * from 1001, ignoring EventID. 00058 * 00059 * Collective mode is selected by calling setCollectiveMode(), and 00060 * then behaving externally just like normal mode (i.e. just call 00061 * BeamOn()). 00062 * 00063 * Enhanced mode: 00064 * This class also implements functions that permit enhanced control 00065 * of the simulation -- in particular there is a processOneTrack() 00066 * function that simply processes a single track. 00067 * 00068 * NOTE: the word "track" is ambiguous, and is commonly used as both 00069 * a noun and a verb. Here it is a noun only (referring to a G4Track), 00070 * and the verb is "process". Processing a track consists of stepping 00071 * it through the simulation until it is killed or suspended, while 00072 * holding all secondaries it generates. The secondaries are themselves 00073 * tracks, and must be dealt with before another call to processOneTrack() 00074 * (see below). 00075 * 00076 * Enhanced run control: All of the usual G4beamline initialization 00077 * must occur, so the following code is best put into a BLCallback 00078 * function registered to replace the main loop -- that means it gets 00079 * called after the post-Reference callbacks are called, and that after 00080 * all replace-main-loop callbacks are called the program exits (thus 00081 * avoiding the main loop). 00082 * 00083 * BLRunManager *runmgr = BLRunManager::getObject(); 00084 * runmgr->beginRun(int runid); 00085 * --- loop over events --- 00086 * runmgr->beginEvent(int evid); 00087 * --- loop over tracks --- 00088 * --- get a track 00089 * runmgr->processOneTrack(track); 00090 * --- handle secondaries (see below) 00091 * runmgr->endEvent(); 00092 * runmgr->endRun(); 00093 * // write out NTuples, etc. (if appropriate) 00094 * 00095 * There can be a loop over runs, if necessary. 00096 * 00097 * Note that these enhanced run control routines do not create or handle 00098 * trajectories, and ignore them if present. Hence no visualization 00099 * of tracks is possible. 00100 * 00101 * Note that secondaries are tracks, but they have not yet been assigned 00102 * either a trackid or a parentid. Those should be assigned before the 00103 * track is deferred, as all deferred tracks should have those id-s. 00104 * 00105 * Handling secondaries -- do one of: 00106 * runmgr->discardAllSecondaries(); 00107 * or 00108 * runmgr->processAllSecondariesAndDeferredTracks(); 00109 * or 00110 * runmgr->deferAllSecondaries(); 00111 * // followed sometime later by: 00112 * runmgr->processAllDeferredTracksAndTheirSecondaries(); 00113 * or 00114 * while((track=runmgr->popOneSecondary()) != 0) { 00115 * // NOTE: you cannot call processOneTrack() here! 00116 * // -- that would discard any remaining secondaries. 00117 * if(--- condition ---) 00118 * runmgr->deferOneTrack(track); 00119 * else 00120 * delete track; 00121 * } 00122 * // followed sometime later by: 00123 * runmgr->processAllDeferredTracksAndTheirSecondaries(); 00124 **/ 00125 class BLRunManager : public G4RunManager, public BLCallback { 00126 static BLRunManager *singleton; 00127 G4TrackingManager *trackManager; 00128 G4StackManager *stackManager; 00129 std::vector<BLTrackData> trackVector; 00130 std::vector<BLCollectiveComputation *> computeVector; 00131 bool collectiveMode; 00132 bool keepSecondaries; 00133 G4double stepTime; // collective mode only 00134 G4double deltaT; // collective mode only 00135 G4int currentTrackIndex; // collective mode only 00136 bool rejected; // collective mode only 00137 int nextSecondaryTrackID; // collective mode only 00138 void appendTrack(G4Track *track); // collective mode only 00139 public: 00140 /// returns the (singleton) BLRunManager. 00141 /// Note that G4RunManager::GetRunManager() returns a pointr to 00142 /// the base class, not this class. 00143 static BLRunManager *getObject(); 00144 00145 /// Constructor 00146 BLRunManager(); 00147 00148 /// Destructor 00149 virtual ~BLRunManager(); 00150 00151 /// BeamOn() will track beam events, in normal or collective mode. 00152 /// See setCollectiveMode() below. 00153 virtual void BeamOn(int nEvents, const char *macroFile=0, 00154 G4int n_select=-1); 00155 00156 /// callback() from BLCallback. 00157 virtual void callback(int type); 00158 00159 // Collective run control: 00160 // NONE of these routines are called in normal or enhanced mode. 00161 00162 /// beamOnCollective() processes an entire run in collective mode. 00163 /// In collective mode, a vector of BLTrackData is constructed, and 00164 /// all tracks in it are stepped one time step at a time. This is 00165 /// essentially an inversion of the usual loop over steps inside 00166 /// loops over tracks and events. 00167 /// All secondaries are immediately added to the end of trackVector, 00168 /// if keepSecondaries is true. 00169 /// The Tune and Reference particles are always tracked in normal mode. 00170 void beamOnCollective(int nEvents); 00171 00172 /// registerCollectiveComputation() will register a computation to 00173 /// be called after every time step in collective mode. Multiple 00174 /// computations are called in the order they are registered. 00175 /// Must be called before beamOnCollective() or steps may be missed 00176 /// (usually called in the command() function of a BLCommand, so the 00177 /// order of commands in the input file is preserved). 00178 void registerCollectiveComputation(BLCollectiveComputation *compute) 00179 { computeVector.push_back(compute); } 00180 00181 /// getTrackVector() returns the trackVector. 00182 /// Note: the track vector and its tracks can be modified inside the 00183 /// collective computation (takes effect immediately). 00184 /// As BLRunManager is a singleton, the trackVector never changes. 00185 /// No track is ever removed from the track vector, but they can 00186 /// have their status set to fStopAndKill which effectively does that. 00187 std::vector<BLTrackData> &getTrackVector() { return trackVector; } 00188 00189 /// getCurrentTrackIndex() returns the index in the track vector 00190 /// of the track currently being tracked. Returns -1 if not in 00191 /// collective mode, or if tracking is not in progress. 00192 int getCurrentTrackIndex() { return currentTrackIndex; } 00193 00194 /// getCollectiveMode() returns true if in collective mode. 00195 bool getCollectiveMode() { return collectiveMode; } 00196 00197 /// setCollectiveMode() arranges for calls to BeamOn() to call 00198 /// beamOnCollective() rather than the usual G4RunManager function. 00199 /// This is for BEAM only; TUNE and REFERENCE particles are always 00200 /// tracked normally via G4RunManager::BeamOn(). 00201 /// _keepSecondaries applies ONLY in collective mode (in normal mode 00202 /// secondaries are always kept, in enhanced mode that is up to 00203 /// the user). 00204 void setCollectiveMode(bool flag=true, bool _keepSecondaries=true, 00205 G4double _deltaT=1.0*ns) 00206 { collectiveMode=flag; 00207 keepSecondaries=_keepSecondaries; 00208 setDeltaT(_deltaT); 00209 } 00210 00211 /// getStepTime() will get the ste time. 00212 /// NOTE: inside collectiveStep() this function returns the time 00213 /// of the current step; outside collectiveStep() it returns the 00214 /// next step time (i.e. the time at which tracks should be suspended). 00215 G4double getStepTime() { return stepTime; } 00216 00217 /// setStepTime() will set the time for the next step in time. 00218 /// If called from within collectiveStep(), it may confuse other 00219 /// classes collectiveStep()-s. If called from beginCollectiveTracking() 00220 /// it sets the time of the first step. 00221 void setStepTime(G4double t) { stepTime = t; } 00222 00223 /// getDeltaT() will get the time interval between time steps 00224 G4double getDeltaT() { return deltaT; } 00225 00226 /// setDeltaT() will set the time interval between time steps. 00227 /// All tracks will be stepped to the same value of global time. 00228 /// Initiall value is 1.0*ns. deltaT <= 0.0 is an error. 00229 /// Can be called at any time; if called within collectiveStep() 00230 /// it will affect the immediately following time step. 00231 void setDeltaT(G4double dt) { BLAssert(dt>0.0); deltaT = dt; } 00232 00233 /// rejectCollectiveStep() will cause the run manager to abandon 00234 /// tracking (if in progress) and to discard all of the current 00235 /// tracks, replacing the vector with a copy saved immediately 00236 /// after the previous call to collectiveStep() returned. 00237 /// Thus the current collective step will be repeated as if 00238 /// no tracking had been performed since the previous return from 00239 /// collectiveStep(). For this to make sense, setDeltaT() 00240 /// should have been called with a smaller value 00241 /// -- this is intended for a collective algorithm to dynamically 00242 /// adapt the value of deltaT. 00243 /// Note: the first step is to bring all beam tracks to a common 00244 /// value of global time, and rejecting this step is ignored. 00245 void rejectCollectiveStep() { rejected = true; } 00246 00247 00248 00249 // Enhanced run control: 00250 // NONE of these routines are called in normal run-control mode. 00251 // NONE of these routines are called in collective run-control mode. 00252 // In most cases these routines should be called from a callback() 00253 // function registered to replace the main loop (exception: calling 00254 // processOneTrack() as "pre-tracking" for a collective computation). 00255 00256 /// beginRun() begins a run in enhanced run-control mode. 00257 void beginRun(int runid=0); 00258 00259 /// endRun() ends a run in enhanced run-control mode. 00260 void endRun(); 00261 00262 /// beginEvent() begins an event in enhanced run-control mode. 00263 void beginEvent(int evid=0); 00264 00265 /// endEvent() ends an event in enhanced run-control mode. 00266 void endEvent(); 00267 00268 /// getNextBeamEventAndTrack() will get the next beam event and track. 00269 /// returns true if one is returned, false if none are left. 00270 /// Neithe event actions nor track actions are performed. 00271 /// Both *pevent and *ptrack must be deleted. 00272 bool getNextBeamEventAndTrack(G4Event **pevent, G4Track **ptrack); 00273 00274 /// processOneTrack() tracks a single track until it is suspended or 00275 /// killed. 00276 /// Note that before tracking begins it deletes any held secondaries, 00277 /// so if secondaries are not to be discarded they must be deferred 00278 /// or processed before the next call. 00279 /// You may want to call BLManager::setState(SPECIAL) to prevent 00280 /// this track from being entered into NTuples for virtualdetectors 00281 /// (etc.) encountered during processing. 00282 /// This function MUST NOT BE CALLED if any track is being processed 00283 /// via any means (e.g. the normal BeamOn() is executing). 00284 void processOneTrack(G4Track *track); 00285 00286 /// discardAllSecondaries() will discard all secondaries. 00287 void discardAllSecondaries(); 00288 00289 /// deferAllSecondaries() will defer all secondaries until 00290 /// deferred tracks are processed. Returns # tracks deferred. 00291 int deferAllSecondaries(int secondaryid=10000, int parentid=-1); 00292 00293 /// deferOneTrack() will defer a single track until deferred tracks 00294 /// are processed. 00295 void deferOneTrack(G4Track *track); 00296 00297 /// processAllSecondariesAndDeferredTracks() will defer all secondaries 00298 /// and then process all deferred tracks and their secondaries. 00299 /// returns # tracks processed. 00300 int processAllSecondariesAndDeferredTracks(int secondaryid=10000, 00301 int parentid=-1) { 00302 int first = secondaryid; 00303 secondaryid += deferAllSecondaries(secondaryid,parentid); 00304 secondaryid += 00305 processAllDeferredTracksAndTheirSecondaries(secondaryid); 00306 return secondaryid-first; 00307 } 00308 00309 /// processAllDeferredTracksAndTheirSecondaries() will process all 00310 /// deferred tracks, including processing their secondaries. 00311 /// Returns # tracks processed. 00312 int processAllDeferredTracksAndTheirSecondaries(int trackid=10000); 00313 00314 /// popOneSecondary() will return a pointer to one secondary, removing 00315 /// it from the list of secondaries; returns NULL if no more. 00316 /// Intended to be called in a loop immediately after a call to 00317 /// processOneTrack(). The returned pointer should eventually be 00318 /// deleted. Order is LIFO. 00319 G4Track *popOneSecondary(); 00320 00321 /// popOneDeferredTrack() will return a pointer to one deferred 00322 /// track, removing it from the list of deferred tracks; returns 00323 /// NULL if no more. The returned pointer should eventually be 00324 /// deleted. Note that a loop calling this function can defer other 00325 /// tracks or secondaries, which simply extends the loop. 00326 /// Order is LIFO. 00327 G4Track *popOneDeferredTrack(); 00328 00329 // protected functions made public (for testing and special purposes): 00330 void RunInitialization() { G4RunManager::RunInitialization(); } 00331 void DoEventLoop(G4int n_event,const char* macroFile,G4int n_select) { 00332 G4RunManager::DoEventLoop(n_event,macroFile,n_select); 00333 } 00334 void RunTermination() { G4RunManager::RunTermination(); } 00335 G4EventManager *getEventManager() { return eventManager; } 00336 G4TrackingManager *getTrackingManager() { return trackManager; } 00337 G4StackManager *getStackManager() { return stackManager; } 00338 void setCurrentEvent(G4Event *ev) { currentEvent = ev; } 00339 G4bool getRunAborted() { return runAborted; } 00340 }; 00341 00342 #endif // BLRUNMANAGER_H