BLRunManager.hh

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