BLManager Class Reference

#include <BLManager.hh>

List of all members.


Detailed Description

BLManager is the overall manager for g4beamline, managing all aspects of the run.

Note it is uses all the G4 user action classes, implemented as classes under BLManager.

TrackID-s: It is non-trivial to preserve TrackID-s from input files, because Geant4 considers TrackID to be an internal variable, and it assigns them in a non-customizable manner. So BLManager keeps a trackIDMap[] that converts from internal (Geant4) to external (User) TrackID-s. It also uses BLTrackInfo to associate the external IDs with the track -- the external trackID is determined either in BLBeam (when the track is created), or in BLManager::PreUserTrackingAction for secondaries.

Public Member Functions

void delayedConstruction ()
 delayedConstruction() performs things which must wait until all static initializers have executed (e.g. in Geant4 routines).
 ~BLManager ()
 Destructor.
BLManagerState getState ()
 getState() returns the current state.
void setState (BLManagerState _state)
 setState sets the state
int getSteppingVerbose ()
 getSteppingVerbose() returns steppingVerbose. NOTE: use this during tracking instead of Param.getInt("steppingVerbose");
void setSteppingVerbose (int v)
 setSteppingVerbose() updates steppingVerbose. -- NOTE: many other classes relay on the Parameter, not the valud in this class.
int getEventTimeLimit ()
 getEventTimeLimit() returns the CPU time limit for events (seconds). -1 mean infinite.
void setEventTimeLimit (int sec)
 setEventTimeLimit() sets the CPU time limit for events (seconds). -1 mean infinite.
G4int getEventID () const
 getEventID() gets the current eventID;
void setEventID (int evId)
 setEventID() sets the current eventID;
bool skipEvent (int EventID)
 skipEvent() will check the eventList vector (which comes from the eventCut file).
void incrEventsProcessed (int eventID)
 incrEventsProcessed() will increment eventsProcessed. For special uses only (e.g.MPI).
void registerSteppingAction (SteppingAction *sa)
 registerSteppingAction() registers a SteppingAction to be called for each step (regardless of state).
void registerSteppingAction (G4VPhysicalVolume *physicalVol, SteppingAction *sa)
 registerSteppingAction() registers a SteppingAction to be called for each step (regardless of state), for every step involving the physicalVol. The callback will be called if physicalVol is either the pre- or post-step physical volume (once if both). If physicalVol==0 it is called every reference particle step. LIMITATION: only one callback can be registered for a given physicalVol (except 0).
void registerTuneParticleStep (G4VPhysicalVolume *physicalVol, SteppingAction *sa)
 registerTuneParticleStep() registers a SteppingAction to be called for every tune particle step involving the physicalVol. The callback will be called if physicalVol is either the pre- or post-step physical volume (once if both). If physicalVol==0 it is called every tune particle step. LIMITATION: only one callback can be registered for a given physicalVol (except 0).
void registerReferenceParticleStep (G4VPhysicalVolume *physicalVol, SteppingAction *sa)
 registerReferenceParticleStep() registers a SteppingAction to be called for every reference particle step involving the physicalVol. The callback will be called if physicalVol is either the pre- or post-step physical volume (once if both). If physicalVol==0 it is called every reference particle step. LIMITATION: only one callback can be registered for a given physicalVol (except 0).
void registerBeamStep (G4VPhysicalVolume *physicalVol, SteppingAction *sa)
 registerBeamStep() registers a SteppingAction to be called for every beam step involving the physicalVol. (beam particles are everything except reference and tune.) The callback will be called if physicalVol is either the pre- or post-step physical volume (once if both). if physicalVol==0 it is called every beam step. LIMITATION: only one callback can be registered for a given physicalVol (except 0).
void registerZStep (G4double z, ZSteppingAction *sa, G4int when=7)
 registerZStep() will force a step to occur near the given z position, and will call the ZSteppingAction for it, interpolating to the desired z value (Centerline coords). when is a bitwise OR of 1=tune, 2=reference, 4=beam.
void registerStackingAction (StackingAction *sa)
 registerStackingAction() registers a StackingAction to be called by the Geant4 stacking manager.
void setSteppingFormat ()
 setSteppingFormat() sets the verbose printing format according to parameter "steppingFormat".
G4String getFormatHelp ()
 getFormatHelp() returns a string with help text about valid format items.
void appendVerboseFormat (G4String fmt)
 appendVerboseFormat() appends fmt to the format for printing when param steppingVerbose is nonzero
void steppingVerbosePrint (const G4Step *step, const G4Track *track, int header=0)
 steppingVerbosePrint() will print this step according to the current verboseFormat. Prints header if header != 0).
void setSteppingManager (G4SteppingManager *p)
 setSteppingManager() sets the pointer to the current G4SteppingManager.
G4SteppingManager * getSteppingManager ()
 getSteppingManager() returns a pointer to the current G4SteppingManager.
void setTrackingManager (G4TrackingManager *p)
 setTrackingManager() sets the pointer to the current G4TrackingManager.
void initialize ()
 initialize() will initialize the BLManager object, and the geant4 kernel, thus constructing the geometry in the world group. Note that registerPhysics() must be called before initialize() (normally done by a "physics" command in the input file).
void trackTuneAndReferenceParticles ()
 trackTuneAndReferenceParticles() will generate and track the tune particle and then the reference particle.
void trackBeam ()
 trackBeam() will generate the defined beam and track each event.
void displayVisual ()
 displayVisual() will display the detector visually.
void displayGeometry (G4VPhysicalVolume *phys=0, int level=0)
 displayGeometry() will display the geant4 geometry. This is a hierarchical ASCII listing of all volumes. if phys==0 then use the worldPhysicalVolume.
void registerPhysics (BLPhysics *_physics)
 registerPhysics() registers the BLPhysics object. It also sets the physics list to the BLRunManager, so following commands can find particle by name.
BLPhysicsgetPhysics ()
 getPhysics returns the registered BLPhysics object.
void registerBeam (BLBeam *_beam)
 registerBeam() registers a BLBeam object for beam generation. Multiple BLBeam objects can be registered, used in order.
void registerReference (BLBeam *_beam)
 registerReference() registers a BLBeam object for reference particle generation. Multiple BLBeam objects can be registered, used in order.
int nReference ()
 nReference() returns the number of reference particles registered.
void registerRunAction (RunAction *a, G4bool beamOnly=true)
 registerRunAction() registers a UserRunAction. If beamOnly is true (the default), the callback is made only if state==BEAM.
void registerEventAction (EventAction *a, G4bool beamOnly=true)
 registerEventAction() registers a UserEventAction. If beamOnly is true (the default), the callback is made only if state==BEAM.
void registerTrackingAction (TrackingAction *a)
 registerTrackingAction() registers a UserTackingAction.
void registerUserCode (BLUserCode *instance)
 registerUserCode() registers a BLUserCode instance.
std::vector< BLUserCode * > getUserCodeInstances (G4String type)
 getUserCodeInstances() returns a vector of all registered instances of BLUserCode with the specified type.
void registerCallback (BLCallback *cb, int type)
 registerCallback() registers a BLCallback. type=0 for pre-Tune particle, type=1 for post-Reference (pre-beam tracking), type=2 for post-beam tracking. type=3 for replacing the main program loop. type=4 for visualization. NOTE: if there are type=3 callbacks, when the last one returns the closes up by summarizing NTuples and callng handleCallbacks(2), and then the program exits. This prevents the main program loop from executing. type=3 callbacks are called just after the type=1 callbacks (i.e. after the post-Reference callbacks).
void handleCallbacks (int type)
 handleCallbacks() calls all applicable registered callbacks. type=0 for pre-reference particle, type=1 for post-center (pre-beam tracking), type=2 for post-beam tracking (just before program exit). type=3 for replacing the main program loop. type=4 for visualization.
G4VPhysicalVolume * getWorldPhysicalVolume ()
 getWorldPhysicalVolume() returns a pointer to the world PV. Note that it must already have been consructed, so this function returns NULL before construct() is called (by main).
void setPRNGSeedMethod (PRNGSeedMethod method)
 setPRNGSeedMethod() will set the method used to seed the pseudo random number generator at the start of each event.
PRNGSeedMethod getPRNGSeedMethod ()
 getPRNGSeedMethod() will return the method used to seed the pseudo random number generator at the start of each event.
void setEventCutFile (G4String &f)
 setEventCutFile() will set the file containing event s for cutting. The file is ASCII, with one event number per line; lines beginning with # are ignored.
void UserSteppingAction (const G4Step *step)
 UserSteppingAction() from G4UserSteppingAction.
G4VPhysicalVolume * Construct ()
 Construct() from G4VUserDetectorConstruction.
void PreUserTrackingAction (const G4Track *track)
 PreUserTrackingAction() from G4UserTrackingAction.
void PostUserTrackingAction (const G4Track *track)
 PostUserTrackingAction() from G4UserTrackingAction.
void BeginOfRunAction (const G4Run *run)
 BeginOfRunAction() from G4UserRunAction.
void EndOfRunAction (const G4Run *run)
 EndOfRunAction() from G4UserRunAction.
void BeginOfEventAction (const G4Event *event)
 BeginOfEventAction() from G4UserEventAction.
void EndOfEventAction (const G4Event *event)
 EndOfEventAction() from G4UserEventAction.
void GeneratePrimaries (G4Event *event)
 GeneratePrimaries() from G4VUserPrimaryGeneratorAction.
G4ClassificationOfNewTrack ClassifyNewTrack (const G4Track *)
 ClassifyNewTrack() from G4UserStackingAction.
void NewStage ()
 NewStage() from G4userStackingAction.
void PrepareNewEvent ()
 PrepareNewEvent() from G4StackingAction.
void clearTrackIDMap ()
 clearTrackIDMap() clears the TrackID map.
void setNextSecondaryTrackID (int next)
 setNextSecondaryTrackID() sets the external TrackID for the next secondary track. Automatically incremented for subsequent secondaries.
int getNextSecondaryTrackID ()
int getExternalTrackID (const G4Track *track)
 getExternalTrackID() returns the external TrackID for the given track. In collective mode, internal and external trackID-s are the same.
int getExternalParentID (const G4Track *track)
 getExternalParentID() returns the external ParentID for the given track. In collective mode, internal and external trackID-s are the same.
void setExternalTrackID (G4Track *track, int trackID, int parentID)
 setExternalTrackID() will set the external trackID and parentID. if trackID<0 uses nextSecondarTrackID++.
int getPrimaryTrackID ()
 getPrimaryTrackID() returns the primaryTrackID set by the beam command
int getPrimaryParentID ()
 getPrimaryParentID() returns the primaryParentID set by the beam command
void setPrimaryTrackID (int t, int p)
 setPrimaryTrackID() sets track and parent IDs for a primary track.
G4bool Notify (const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
 Notify() from G4VExceptionHandler.
void exceptionSummary ()
 exceptionSummary() prints a summary of all exceptions

Static Public Member Functions

static BLManagergetObject ()
 getObject() will return a pointer to the single BLManager object, creating it if necessary. Does only the immediate constructor, not delayedConstruction(). that means it is Ok to register capabilities, but not much else.
static bool isInitialized ()
 isInitialized() returns true if the BLManager has been initialized.

Private Member Functions

 BLManager ()
 private constructor -- immediate construction only.
void insertZStep (std::vector< ZStep > &vector, G4double z, ZSteppingAction *action)
void readEventCutFile ()

Private Attributes

BLRunManagerrunManager
BLManagerState state
G4int steppingVerbose
unsigned int beamIndex
G4int histoUpdate
time_t startRun
time_t startEvent
G4int eventTimeLimit
BLPhysicsphysics
G4VPhysicalVolume * worldPhysicalVolume
G4int eventID
G4int trackID
G4ParticleDefinition * trackParticleDef
G4int fatalExceptions
G4int eventsAborted
G4int tracksAborted
G4int warnings
G4int prevEventID
G4int eventsProcessed
bool endRun
PRNGSeedMethod seedMethod
G4String eventCutFile
std::set< int > eventList
G4SteppingManager * fpSteppingManager
G4TrackingManager * fpTrackingManager
std::vector< BLBeam * > beamVector
std::vector< BLBeam * > referenceVector
std::vector< RunAction * > runActionVector
std::vector< RunAction * > beamRunActionVector
std::vector< EventAction * > eventActionVector
std::vector< EventAction * > beamEventActionVector
std::vector< TrackingAction * > trackingActionVector
std::vector< BLCallback * > preReferenceCallbackVector
std::vector< BLCallback * > postReferenceCallbackVector
std::vector< BLCallback * > postTrackingCallbackVector
std::vector< BLCallback * > replaceMainLoopCallbackVector
std::vector< BLCallback * > visualizationCallbackVector
std::vector< SteppingAction * > allStepVector
std::map< G4VPhysicalVolume
*, SteppingAction * > 
allStepMap
std::map< G4VPhysicalVolume
*, SteppingAction * > 
tpStepMap
std::map< G4VPhysicalVolume
*, SteppingAction * > 
rpStepMap
std::vector< SteppingAction * > tpStepVector
std::vector< SteppingAction * > rpStepVector
std::map< G4VPhysicalVolume
*, SteppingAction * > 
beamStepMap
std::vector< SteppingAction * > beamStepVector
std::vector< int > verboseFormat
G4double zTolerance
std::vector< ZSteptuneZStep
std::vector< ZStepreferenceZStep
std::vector< ZStepbeamZStep
std::vector< ZStep > * currentZStep
unsigned indexZStep
G4double prevZ
G4int nStuck
std::vector< StackingAction * > stackingActionVector
int nextSecondaryTrackID
std::map< G4int, G4int > trackIDMap
int primaryTrackID
int primaryParentID
std::vector< BLUserCode * > userCodeVector
std::map< G4String, int > exceptionCount

Static Private Attributes

static BLManagerblManager = 0
static bool initialized = false

Classes

class  EventAction
class  PrimaryGeneratorAction
class  RunAction
class  StackingAction
class  SteppingAction
class  TrackingAction
struct  ZStep
class  ZSteppingAction


Constructor & Destructor Documentation

BLManager::BLManager (  )  [private]

private constructor -- immediate construction only.

References beamIndex, beamZStep, blManager, CL, currentZStep, endRun, EVENT_NUMBER, eventCutFile, eventID, eventsAborted, eventsProcessed, eventTimeLimit, fatalExceptions, BLParam::getDouble(), GLOBAL, histoUpdate, IDLE, indexZStep, BLAlarm::init(), KE, NEWLINE, nextSecondaryTrackID, NSTEP, nStuck, Param, physics, prevEventID, prevZ, primaryParentID, primaryTrackID, PROCESS, referenceZStep, runManager, seedMethod, startEvent, startRun, state, STEP, steppingVerbose, trackID, trackParticleDef, tracksAborted, tuneZStep, verboseFormat, VOL, warnings, worldPhysicalVolume, and zTolerance.

Referenced by getObject().

00199                      : G4VUserDetectorConstruction(), 
00200                 G4VExceptionHandler(), eventList(),
00201                 beamVector(), referenceVector(),
00202                 runActionVector(), beamRunActionVector(),
00203                 eventActionVector(), beamEventActionVector(),
00204                 trackingActionVector(), preReferenceCallbackVector(),
00205                 postReferenceCallbackVector(), postTrackingCallbackVector(),
00206                 replaceMainLoopCallbackVector(), visualizationCallbackVector(),
00207                 allStepVector(), allStepMap(), tpStepMap(), rpStepMap(),
00208                 tpStepVector(), rpStepVector(),
00209                 beamStepMap(), beamStepVector(), verboseFormat(),
00210                 tuneZStep(), referenceZStep(), beamZStep(), 
00211                 stackingActionVector(), trackIDMap(), userCodeVector(),
00212                 exceptionCount()
00213 {
00214         // Parameter definitions (variables are all unused):
00215         BLSetParam unused_1("histoFile","g4beamline",
00216                                 "Default (Root) NTuple output filename");
00217         BLSetParam unused_2("histoUpdate","0",
00218                                 "Output update interval (events)");
00219         BLSetParam unused_3("viewer","none",
00220                                 "Visualization driver selected (default=none)");
00221         BLSetParam unused_4("eventTimeLimit","30",
00222                                 "CPU Time Limit (sec)");
00223         BLSetParam unused_5("steppingVerbose","0",
00224                                 "Set nonzero to print each step");
00225         BLSetParam unused_6("steppingFormat", "N GLOBAL CL KE STEP VOL PROCESS",
00226                                 "Format for printing steps");
00227         BLSetParam unused_7("zTolerance","2.0",
00228                                 "Tolerance for Z steps (mm)");
00229 
00230         if(blManager)
00231                 G4Exception("BLManager","Object Already Exists",FatalException,
00232                                                                         "");
00233         blManager = this;
00234 
00235         state = IDLE;
00236         steppingVerbose = 0;
00237 
00238         // ensure ZStep vectors have bookends
00239         tuneZStep.push_back(ZStep(-DBL_MAX,0));
00240         tuneZStep.push_back(ZStep(DBL_MAX,0));
00241         referenceZStep.push_back(ZStep(-DBL_MAX,0));
00242         referenceZStep.push_back(ZStep(DBL_MAX,0));
00243         beamZStep.push_back(ZStep(-DBL_MAX,0));
00244         beamZStep.push_back(ZStep(DBL_MAX,0));
00245         currentZStep = &beamZStep;
00246         indexZStep = 1;
00247         zTolerance = Param.getDouble("zTolerance");
00248         prevZ = -DBL_MAX;
00249 
00250         // default verbose format
00251         verboseFormat.push_back(NSTEP);
00252         verboseFormat.push_back(GLOBAL);
00253         verboseFormat.push_back(CL);
00254         verboseFormat.push_back(KE);
00255         verboseFormat.push_back(STEP);
00256         verboseFormat.push_back(VOL);
00257         verboseFormat.push_back(PROCESS);
00258         verboseFormat.push_back(NEWLINE);
00259 
00260         runManager = 0;
00261         state = IDLE;
00262         beamIndex = 0;
00263         histoUpdate = 0;
00264         startRun = 0;
00265         startEvent = 0;
00266         eventTimeLimit = 30;
00267         physics = 0;
00268         worldPhysicalVolume = 0;
00269         eventID = 1;
00270         trackID=-9999;
00271         trackParticleDef = 0;
00272         fatalExceptions = 0;
00273         eventsAborted=0;
00274         tracksAborted = 0;
00275         warnings=0;
00276         prevEventID = 0;
00277         eventsProcessed = 0;
00278         endRun = false;
00279         seedMethod = EVENT_NUMBER;
00280         eventCutFile = "";
00281         nStuck = 0;
00282         primaryTrackID = -1;
00283         primaryParentID = -1;
00284         nextSecondaryTrackID = 1001;
00285 
00286         BLAlarm::init();
00287 
00288         G4StateManager * stateManager = G4StateManager::GetStateManager();
00289         stateManager->SetExceptionHandler(this);
00290 }

BLManager::~BLManager (  ) 

Destructor.

References blManager.

00300 {
00301         // cannot delete runManager, as it deletes the detector, 
00302         // steppingaction, eventaction, ..., and that becomes multiple
00303         // deletes of this object. But we are exiting, so we don't really
00304         // need to delete anything.
00305         //delete runManager;
00306 
00307         // Instead of deleting things, just open the geometry to avoid 
00308         // warnings when exit() is called.
00309         G4GeometryManager::GetInstance()->OpenGeometry();
00310 
00311         blManager = 0;
00312 }


Member Function Documentation

void BLManager::insertZStep ( std::vector< ZStep > &  vector,
G4double  z,
ZSteppingAction action 
) [private]

Referenced by registerZStep().

00386 {
00387         // the bookends in vector (see BLManager()) ensure this works
00388         std::vector<ZStep>::iterator i;
00389         for(i=vector.begin(); i<vector.end(); ++i) {
00390                 if(z < i->z) {  // preserve order for ==
00391                         ZStep zs(z,action);
00392                         vector.insert(i,zs);
00393                         break;
00394                 }
00395         }
00396 }

void BLManager::readEventCutFile (  )  [private]

References eventCutFile, and eventList.

Referenced by BeginOfRunAction().

01506 {
01507         eventList.clear();
01508         if(eventCutFile.size() == 0) return;
01509 
01510         FILE *in = fopen(eventCutFile.c_str(),"r");
01511         if(!in) {
01512                 printf("WARNING: Cannot read eventCutFile '%s' -- IGNORED\n",
01513                                                         eventCutFile.c_str());
01514                 return;
01515         }
01516 
01517         char line[1025];
01518         while(fgets(line,sizeof(line),in)) {
01519                 assert(strlen(line) < sizeof(line)-1);
01520                 if(line[0] == '#') continue;
01521                 long ev = strtol(line,0,0);
01522                 eventList.insert(ev);
01523         }
01524 
01525         fclose(in);
01526 
01527         // eventList is ignored if empty, so insert a crazy entry
01528         if(eventList.size() == 0) eventList.insert(0x80000000);
01529 }

BLManager * BLManager::getObject (  )  [static]

getObject() will return a pointer to the single BLManager object, creating it if necessary. Does only the immediate constructor, not delayedConstruction(). that means it is Ok to register capabilities, but not much else.

References BLManager(), and blManager.

Referenced by FOR009NTuple::appendRow(), TrackFileNTuple::appendRow(), BLRunManager::BeamOn(), BLRunManager::beamOnCollective(), BLCMDmovie::beamStep(), BLRunManager::beginEvent(), BLRunManager::beginRun(), BLVirtualDetectorNTuple::BLVirtualDetectorNTuple(), BLCMDzntuple::callback(), BLCMDtrackermode::callback(), BLCMDtotalenergy::callback(), BLCMDshowmaterial::callback(), BLCMDgeometry::callback(), BLCMDtrackcuts::ClassifyNewTrack(), BLCMDcollective::collectiveStep(), BLCMDzntuple::command(), BLCMDusertrackfilter::command(), BLCMDtune::command(), BLCMDtrackermode::command(), BLCMDtracker::command(), BLCMDtrackcuts::command(), BLCMDtrace::command(), BLCMDtotalenergy::command(), BLCMDtimentuple::command(), BLCMDshowmaterial::command(), BLCMDsetdecay::command(), BLCMDreference::command(), BLCMDrandomseed::command(), BLCMDprofile::command(), BLCMDprobefield::command(), BLCMDprintfield::command(), BLCMDprintf::command(), BLCMDphysics::command(), BLCMDparticlesource::command(), BLCMDparticlecolor::command(), BLCMDntuple::command(), BLCMDnewparticlentuple::command(), BLCMDmovie::command(), BLCMDlist::command(), BLCMDhelp::command(), BLCMDgeometry::command(), BLCMDg4ui::command(), BLCMDfieldntuple::command(), BLCMDeventcuts::command(), BLCMDcosmicraybeam::command(), BLCMDbug1021::command(), BLCMDbeamlossntuple::command(), BLCMDbeam::command(), BLCMDvirtualdetector::construct(), BLCMDusertrackfilter::construct(), BLCMDtubs::construct(), BLCMDtrap::construct(), BLCMDtrackerplane::construct(), BLCMDtorus::construct(), BLCMDsphere::construct(), BLCMDsolenoid::construct(), BLCMDpolycone::construct(), BLCMDpillbox::construct(), BLCMDparticlefilter::construct(), BLCMDmultipole::construct(), BLCMDidealsectorbend::construct(), BLCMDhelicaldipole::construct(), BLCMDgenericquad::construct(), BLCMDgenericbend::construct(), BLCMDextrusion::construct(), BLCMDcorner::construct(), BLCMDbox::construct(), BLRunManager::endEvent(), BLCMDparticlecolor::EndOfEventAction(), BLRunManager::endRun(), SetDecayInstance::fatalError(), BLCMDtracker::fitTrack(), BLCMDreference::generateReferenceParticle(), BLRunManager::getNextBeamEventAndTrack(), BLCMDparticlecolor::getVisAttributes(), BLCMDtracker::handlePreviousTracks(), BLCMDparam::help(), BLCoordinates::init(), BLCMDlist::listParticles(), main(), BLCMDparticlesource::nextBeamEvent(), BLCMDbeam::nextBeamEvent(), BLCMDtracker::operator()(), ParticleFilterPlacement::ParticleFilterPlacement(), BLCoordinates::ReferenceCoordinates::PostUserTrackingAction(), BLCMDtracker::PostUserTrackingAction(), BLCMDreference::PostUserTrackingAction(), BLCMDparticlecolor::PostUserTrackingAction(), BLCMDbeamlossntuple::PostUserTrackingAction(), BLCoordinates::ReferenceCoordinates::PreUserTrackingAction(), BLCMDtracker::PreUserTrackingAction(), BLCMDreference::PreUserTrackingAction(), BLCMDparticlecolor::PreUserTrackingAction(), BLCMDnewparticlentuple::PreUserTrackingAction(), BLCoordinates::ReferenceCoordinates::ReferenceCoordinates(), BLUserCode::registerUserCode(), BLBeam::setRandomSeedToGenerate(), BLEvaluator::setTrackVariables(), BLCoordinates::start(), steppingVerbosePrint(), Surface::Surface(), TraceNTuple::TraceNTuple(), BLCoordinates::ReferenceCoordinates::UserSteppingAction(), BLVirtualDetectorNTuple::UserSteppingAction(), BLCMDusertrackfilter::UserSteppingAction(), TrackerPlaneInstance::UserSteppingAction(), TimeNTuple::UserSteppingAction(), PillboxField::UserSteppingAction(), ParticleFilterPlacement::UserSteppingAction(), BLCMDzntuple::Entry::UserZSteppingAction(), BLCMDtune::UserZSteppingAction(), BLCMDtracker::UserZSteppingAction(), BLCMDreference::UserZSteppingAction(), BLCMDprofile::Entry::UserZSteppingAction(), and BLCMDprintf::UserZSteppingAction().

00379 {
00380         if(!blManager) new BLManager();
00381         return blManager;
00382 }

void BLManager::delayedConstruction (  ) 

delayedConstruction() performs things which must wait until all static initializers have executed (e.g. in Geant4 routines).

References runManager.

Referenced by main().

00293 {
00294         runManager = new BLRunManager();
00295         G4StateManager * stateManager = G4StateManager::GetStateManager();
00296         stateManager->SetExceptionHandler(this);
00297 }

BLManagerState BLManager::getState (  )  [inline]

void BLManager::setState ( BLManagerState  _state  )  [inline]

int BLManager::getSteppingVerbose (  )  [inline]

getSteppingVerbose() returns steppingVerbose. NOTE: use this during tracking instead of Param.getInt("steppingVerbose");

References steppingVerbose.

Referenced by BLVirtualDetectorNTuple::BLVirtualDetectorNTuple(), BLCoordinates::init(), and BLCoordinates::start().

00194 { return steppingVerbose; }

void BLManager::setSteppingVerbose ( int  v  )  [inline]

setSteppingVerbose() updates steppingVerbose. -- NOTE: many other classes relay on the Parameter, not the valud in this class.

References steppingVerbose.

00198 { steppingVerbose = v; }

int BLManager::getEventTimeLimit (  )  [inline]

getEventTimeLimit() returns the CPU time limit for events (seconds). -1 mean infinite.

References eventTimeLimit.

00202 { return eventTimeLimit; }

void BLManager::setEventTimeLimit ( int  sec  )  [inline]

setEventTimeLimit() sets the CPU time limit for events (seconds). -1 mean infinite.

References eventTimeLimit.

Referenced by BLRunManager::beamOnCollective().

00206 { eventTimeLimit = sec; }

G4int BLManager::getEventID (  )  const [inline]

void BLManager::setEventID ( int  evId  )  [inline]

bool BLManager::skipEvent ( int  EventID  )  [inline]

skipEvent() will check the eventList vector (which comes from the eventCut file).

References eventID, and eventList.

Referenced by GeneratePrimaries().

00217             { return eventList.size() > 0 && eventList.count(eventID) == 0; }

void BLManager::incrEventsProcessed ( int  eventID  ) 

incrEventsProcessed() will increment eventsProcessed. For special uses only (e.g.MPI).

References eventsProcessed, BLRunManager::getCollectiveMode(), runManager, startRun, and BLTime::time().

Referenced by EndOfEventAction().

01087 {
01088         // print event number, if appropriate
01089         ++eventsProcessed;
01090         if(!runManager->getCollectiveMode()) {
01091                 if(eventsProcessed <= 10 ||
01092                    (eventsProcessed < 100 && eventsProcessed%10 == 0) ||
01093                    (eventsProcessed < 1000 && eventsProcessed%100 == 0) ||
01094                    eventsProcessed%1000 == 0) {
01095                         printf("Event %d Completed",eventID);
01096                         int t = BLTime::time() - startRun;
01097                         if(t <= 0) t = 1;
01098                         printf("  %d events  realTime=%d sec  %.1f ev/sec",
01099                                 eventsProcessed,t,(double)eventsProcessed/t);
01100                         printf("\n");
01101                         fflush(stdout);
01102                 }
01103         }
01104 }

void BLManager::registerSteppingAction ( SteppingAction sa  )  [inline]

void BLManager::registerSteppingAction ( G4VPhysicalVolume *  physicalVol,
SteppingAction sa 
) [inline]

registerSteppingAction() registers a SteppingAction to be called for each step (regardless of state), for every step involving the physicalVol. The callback will be called if physicalVol is either the pre- or post-step physical volume (once if both). If physicalVol==0 it is called every reference particle step. LIMITATION: only one callback can be registered for a given physicalVol (except 0).

References allStepMap, and allStepVector.

00238                 { if(physicalVol != 0) allStepMap[physicalVol] = sa;
00239                   else allStepVector.push_back(sa); }

void BLManager::registerTuneParticleStep ( G4VPhysicalVolume *  physicalVol,
SteppingAction sa 
) [inline]

registerTuneParticleStep() registers a SteppingAction to be called for every tune particle step involving the physicalVol. The callback will be called if physicalVol is either the pre- or post-step physical volume (once if both). If physicalVol==0 it is called every tune particle step. LIMITATION: only one callback can be registered for a given physicalVol (except 0).

References tpStepMap, and tpStepVector.

Referenced by BLCMDpillbox::construct().

00250                 { if(physicalVol != 0) tpStepMap[physicalVol] = sa; 
00251                   else tpStepVector.push_back(sa); }

void BLManager::registerReferenceParticleStep ( G4VPhysicalVolume *  physicalVol,
SteppingAction sa 
) [inline]

registerReferenceParticleStep() registers a SteppingAction to be called for every reference particle step involving the physicalVol. The callback will be called if physicalVol is either the pre- or post-step physical volume (once if both). If physicalVol==0 it is called every reference particle step. LIMITATION: only one callback can be registered for a given physicalVol (except 0).

References rpStepMap, and rpStepVector.

Referenced by BLCMDmovie::command(), BLCMDvirtualdetector::construct(), BLCMDtrackerplane::construct(), and BLCoordinates::ReferenceCoordinates::ReferenceCoordinates().

00263                 { if(physicalVol != 0) rpStepMap[physicalVol] = sa; 
00264                   else rpStepVector.push_back(sa); }

void BLManager::registerBeamStep ( G4VPhysicalVolume *  physicalVol,
SteppingAction sa 
) [inline]

registerBeamStep() registers a SteppingAction to be called for every beam step involving the physicalVol. (beam particles are everything except reference and tune.) The callback will be called if physicalVol is either the pre- or post-step physical volume (once if both). if physicalVol==0 it is called every beam step. LIMITATION: only one callback can be registered for a given physicalVol (except 0).

References beamStepMap, and beamStepVector.

Referenced by BLCMDtotalenergy::command(), BLCMDtimentuple::command(), BLCMDmovie::command(), BLCMDvirtualdetector::construct(), and BLCMDtrackerplane::construct().

00276                 { if(physicalVol != 0) beamStepMap[physicalVol] = sa; 
00277                   else beamStepVector.push_back(sa); }

void BLManager::registerZStep ( G4double  z,
ZSteppingAction sa,
G4int  when = 7 
)

registerZStep() will force a step to occur near the given z position, and will call the ZSteppingAction for it, interpolating to the desired z value (Centerline coords). when is a bitwise OR of 1=tune, 2=reference, 4=beam.

References beamZStep, insertZStep(), referenceZStep, and tuneZStep.

Referenced by BLCMDzntuple::callback(), BLCMDtune::command(), BLCMDtracker::command(), BLCMDreference::command(), BLCMDprofile::command(), and BLCMDprintf::command().

00399 {
00400         if(when & 1) insertZStep(tuneZStep,z,sa);
00401         if(when & 2) insertZStep(referenceZStep,z,sa);
00402         if(when & 4) insertZStep(beamZStep,z,sa);
00403 }

void BLManager::registerStackingAction ( StackingAction sa  )  [inline]

registerStackingAction() registers a StackingAction to be called by the Geant4 stacking manager.

References stackingActionVector.

Referenced by BLCMDtrackcuts::command().

00288                 { stackingActionVector.push_back(sa); }

void BLManager::setSteppingFormat (  ) 

setSteppingFormat() sets the verbose printing format according to parameter "steppingFormat".

References appendVerboseFormat(), BLParam::getString(), Param, and BLCommand::parseArgs().

Referenced by initialize().

01212 {
01213         G4String fmt = Param.getString("steppingFormat");
01214         BLArgumentVector args;
01215         BLArgumentMap namedArgs;
01216         
01217         // convert "," to " "
01218         G4String s;
01219         for(unsigned i=0; i<fmt.size(); ++i)
01220                 s += (fmt[i]==',' ? ' ' : fmt[i]);
01221 
01222         BLCommand::parseArgs(s,args,namedArgs);
01223         appendVerboseFormat("");
01224         for(unsigned i=0; i<args.size(); ++i)
01225                 appendVerboseFormat(args[i]);
01226         appendVerboseFormat("NEWLINE");
01227 }

G4String BLManager::getFormatHelp (  ) 

getFormatHelp() returns a string with help text about valid format items.

01230 {
01231         G4String s = 
01232                 "        EXT       toggle extended precision (3 more digits)\n"
01233                 "        TAG       print a '>' (useful to grep output)\n"
01234                 "        N         step number\n"
01235                 "        NSTEP     Synonym of N\n"
01236                 "        GLOBAL    X,Y,Z,T in global coords\n"
01237                 "        XYZT      Synonym of GLOBAL\n"
01238                 "        CL        X,Y,Z,dxdz,dydz in CL coords\n"
01239                 "        CLX       extended precision CL\n"
01240                 "        KE        kinetic energy\n"
01241                 "        STEP      step length\n"
01242                 "        STEPLEN   Synonym of STEP\n"
01243                 "        VOL       volume name\n"
01244                 "        VOLNAME   Synonym of VOL\n"
01245                 "        PROCESS   process name\n"
01246                 "        B         magnetic field\n"
01247                 "        E         electric field\n"
01248                 "        P         3-momentum\n"
01249                 "        MAT       material name\n"
01250                 "        ID        event ID, track ID, parent ID\n"
01251                 "        PART      particle name\n"
01252                 "        SEG       centerline coord segment number\n"
01253                 "        WT        weight\n"
01254                 "        NL        <newline>\n"
01255                 "        NEWLINE   Synonym of NL\n"
01256                 "        \\n        Synonym of NL\n";
01257 
01258         return s;
01259 }

void BLManager::appendVerboseFormat ( G4String  fmt  ) 

appendVerboseFormat() appends fmt to the format for printing when param steppingVerbose is nonzero

References B, CL, CLX, E, EXT, GLOBAL, ID, KE, MAT, NEWLINE, NSTEP, P, PART, PROCESS, SEG, STEP, TAG, verboseFormat, VOL, and WT.

Referenced by setSteppingFormat().

01262 {
01263         G4String f;
01264         for(const char *p=fmt.c_str(); *p; ++p)
01265                 f += toupper(*p);
01266         if(f == "") verboseFormat.clear();
01267         else if(f == "TAG") verboseFormat.push_back(TAG);
01268         else if(f == "EXT") verboseFormat.push_back(EXT);
01269         else if(f == "N") verboseFormat.push_back(NSTEP);
01270         else if(f == "NSTEP") verboseFormat.push_back(NSTEP);
01271         else if(f == "GLOBAL") verboseFormat.push_back(GLOBAL);
01272         else if(f == "XYZT") verboseFormat.push_back(GLOBAL);
01273         else if(f == "CL") verboseFormat.push_back(CL);
01274         else if(f == "CLX") verboseFormat.push_back(CLX);
01275         else if(f == "KE") verboseFormat.push_back(KE);
01276         else if(f == "STEP") verboseFormat.push_back(STEP);
01277         else if(f == "STEPLEN") verboseFormat.push_back(STEP);
01278         else if(f == "VOL") verboseFormat.push_back(VOL);
01279         else if(f == "VOLNAME") verboseFormat.push_back(VOL);
01280         else if(f == "PROCESS") verboseFormat.push_back(PROCESS);
01281         else if(f == "B") verboseFormat.push_back(B);
01282         else if(f == "E") verboseFormat.push_back(E);
01283         else if(f == "P") verboseFormat.push_back(P);
01284         else if(f == "MAT") verboseFormat.push_back(MAT);
01285         else if(f == "ID") verboseFormat.push_back(ID);
01286         else if(f == "PART") verboseFormat.push_back(PART);
01287         else if(f == "SEG") verboseFormat.push_back(SEG);
01288         else if(f == "WT") verboseFormat.push_back(WT);
01289         else if(f == "\n") verboseFormat.push_back(NEWLINE);
01290         else if(f == "NL") verboseFormat.push_back(NEWLINE);
01291         else if(f == "NEWLINE") verboseFormat.push_back(NEWLINE);
01292         else printf("BLManager: invalid verbose format '%s'\n",
01293                 fmt.c_str());
01294 }

void BLManager::steppingVerbosePrint ( const G4Step *  step,
const G4Track *  track,
int  header = 0 
)

steppingVerbosePrint() will print this step according to the current verboseFormat. Prints header if header != 0).

References B, BLCOORD_CENTERLINE, CL, CLX, E, EXT, BLCoordinates::getCoords(), getExternalParentID(), getExternalTrackID(), BLGlobalField::GetFieldValue(), getObject(), BLGlobalField::getObject(), BLCoordinates::getRotation(), BLCoordinates::getSegmentCL(), GLOBAL, ID, BLCoordinates::isValid(), KE, MAT, NEWLINE, NSTEP, P, PART, PROCESS, runManager, SEG, STEP, TAG, verboseFormat, VOL, and WT.

Referenced by UserSteppingAction().

01298 {
01299         G4StepPoint *prePoint=0;
01300         G4VPhysicalVolume *preVol=0;
01301         G4StepPoint *postPoint=0;
01302         G4double time=0;
01303         G4ThreeVector position;
01304         G4Material *mat=0;
01305         int eventNum = -999;
01306         int nstep=0;
01307         const char *procName=0;
01308         BLCoordinates *coord = (track!=0 ? (BLCoordinates *)track->GetUserInformation() : 0);
01309         if(coord && !coord->isValid()) coord = 0;
01310 
01311         if(!header) {
01312                 prePoint = step->GetPreStepPoint();
01313                 preVol = prePoint->GetPhysicalVolume();
01314                 postPoint = step->GetPostStepPoint();
01315                 time = track->GetGlobalTime();
01316                 position = track->GetPosition();
01317                 mat = preVol->GetLogicalVolume()->GetMaterial();
01318                 nstep = track->GetCurrentStepNumber();
01319                 procName = (postPoint->GetProcessDefinedStep() ?
01320                   postPoint->GetProcessDefinedStep()->GetProcessName().c_str() :
01321                   "UserLimit"); 
01322                 const G4Event* event = runManager->GetCurrentEvent();
01323                 eventNum = event->GetEventID();
01324         }
01325 
01326         bool extended=false;
01327         for(unsigned i=0; i<verboseFormat.size(); ++i) {
01328                 switch(verboseFormat[i]) {
01329                 case EXT:
01330                         extended = !extended;
01331                         break;
01332                 case TAG:
01333                         if(header) printf(" ");
01334                         else printf(">");
01335                         break;
01336                 case NSTEP:
01337                         if(header) printf(" Step");
01338                         else printf(" %4d",nstep);
01339                         break;
01340                 case GLOBAL:
01341                     if(!extended) {
01342                         if(header) printf("   X(mm)   Y(mm)   Z(mm)    T(ns)");
01343                         else printf(" %7.1f %7.1f %7.1f %8.2f",
01344                                 position[0],position[1],position[2],
01345                                 time);
01346                     } else {
01347                         if(header) printf("      X(mm)      Y(mm)      Z(mm)      T(ns) ");
01348                         else printf(" %10.4f %10.4f %10.4f %11.5f",
01349                                 position[0],position[1],position[2],
01350                                 time);
01351                     }
01352                         break;
01353                 case CL:
01354                     if(extended) goto clx;
01355                     if(header) {
01356                         printf("   CL: X        Y        Z     dxdz    dydz");
01357                     } else {
01358                         G4ThreeVector dir = track->GetMomentumDirection();
01359                         // transform to centerline coordinates
01360                         G4ThreeVector clpos;
01361                         coord->getCoords(BLCOORD_CENTERLINE,clpos);
01362                         dir = coord->getRotation() * dir;
01363                         printf(" %8.1f %8.1f %8.1f %7.4f %7.4f",
01364                                         clpos[0],clpos[1],clpos[2],
01365                                         dir[0]/dir[2],dir[1]/dir[2]);
01366                     }
01367                         break;
01368                 case CLX:
01369 clx:                if(header) {
01370                         printf("   CL: X           Y           Z        dxdz       dydz   ");
01371                     } else {
01372                         G4ThreeVector dir = track->GetMomentumDirection();
01373                         // transform to centerline coordinates
01374                         G4ThreeVector clpos;
01375                         coord->getCoords(BLCOORD_CENTERLINE,clpos);
01376                         dir = coord->getRotation() * dir;
01377                         printf(" %11.4f %11.4f %11.4f %10.7f %10.7f",
01378                                         clpos[0],clpos[1],clpos[2],
01379                                         dir[0]/dir[2],dir[1]/dir[2]);
01380                     }
01381                         break;
01382                 case KE:
01383                     if(!extended) {
01384                         if(header) printf("  KE(MeV)");
01385                         else printf(" %8.1f",track->GetKineticEnergy());
01386                     } else {
01387                         if(header) printf("    KE(MeV) ");
01388                         else printf(" %11.4f",track->GetKineticEnergy());
01389                     }
01390                         break;
01391                 case STEP:
01392                     if(!extended) {
01393                         if(header) printf("  StepLen");
01394                         else printf(" %8.2f",step->GetStepLength());
01395                     } else {
01396                         if(header) printf("    StepLen ");
01397                         else printf(" %11.5f",step->GetStepLength());
01398                     }
01399                         break;
01400                 case VOL:
01401                         if(header) printf(" This Volume     ");
01402                         else printf(" %-16s",preVol->GetName().c_str());
01403                         break;
01404                 case PROCESS:
01405                         if(header) printf(" Process        ");
01406                         else printf(" %-16s",procName);
01407                         break;
01408                 case B:
01409                     if(header) {
01410                         printf("       Bx,By,Bz (Tesla)     ");
01411                     } else {
01412                         G4double point[4], field[6];
01413                         point[0]=position[0],point[1]=position[1],point[2]=position[2];
01414                         point[3]=time;
01415                         BLGlobalField::getObject()->GetFieldValue(point,field);
01416                         G4ThreeVector B(field[0],field[1],field[2]);
01417                         B = coord->getRotation(BLCOORD_CENTERLINE) * B;
01418                         printf(" %8.4f %8.4f %8.4f",
01419                                 B[0]/tesla,B[1]/tesla,B[2]/tesla);
01420                     }
01421                         break;
01422                 case E:
01423                     if(header) {
01424                         printf("      Ex,Ey,Ez (MV/meter)  ");
01425                     } else {
01426                         G4double point[4], field[6];
01427                         point[0]=position[0],point[1]=position[1],point[2]=position[2];
01428                         point[3]=time;
01429                         BLGlobalField::getObject()->GetFieldValue(point,field);
01430                         G4ThreeVector E(field[3],field[4],field[5]);
01431                         E = coord->getRotation(BLCOORD_CENTERLINE) * E;
01432                         printf(" %8.4f %8.4f %8.4f",E[0]/(megavolt/meter),
01433                                 E[1]/(megavolt/meter),E[2]/(megavolt/meter));
01434                     }
01435                         break;
01436                 case P:
01437                   if(!extended) {
01438                     if(header) {
01439                         printf("      Px,Py,Pz (MeV/c)     ");
01440                     } else {
01441                         G4ThreeVector P = track->GetMomentum();
01442                         P = coord->getRotation(BLCOORD_CENTERLINE) * P;
01443                         printf(" %8.1f %8.1f %8.1f",P[0]/(MeV),
01444                                 P[1]/(MeV),P[2]/(MeV));
01445                     }
01446                   } else {
01447                     if(header) {
01448                         printf("          Px,Py,Pz (MeV/c)          ");
01449                     } else {
01450                         G4ThreeVector P = track->GetMomentum();
01451                         P = coord->getRotation(BLCOORD_CENTERLINE) * P;
01452                         printf(" %11.4f %11.4f %11.4f",P[0]/(MeV),
01453                                 P[1]/(MeV),P[2]/(MeV));
01454                     }
01455                   }
01456                         break;
01457                 case MAT:
01458                     if(header) {
01459                         printf(" Material ");
01460                     } else {
01461                         printf(" %-9s",mat->GetName().c_str());
01462                     }
01463                         break;
01464                 case ID:
01465                     if(header) {
01466                         printf(" Event# Trk Prnt");
01467                     } else {
01468                         printf(" %6d %3d %4d",eventNum,
01469                                 BLManager::getObject()->
01470                                     getExternalTrackID(track),
01471                                 BLManager::getObject()->
01472                                     getExternalParentID(track));
01473                     }
01474                         break;
01475                 case PART:
01476                     if(header) {
01477                         printf(" Particle");
01478                     } else {
01479                         printf(" %8s",track->GetDefinition()->GetParticleName().c_str());
01480                     }
01481                         break;
01482                 case SEG:
01483                     if(header) {
01484                         printf(" Seg");
01485                     } else {
01486                         printf(" %3d",coord->getSegmentCL());
01487                     }
01488                         break;
01489                 case WT:
01490                         if(header) {
01491                             printf("  Wt  ");
01492                         } else {
01493                             printf(" %.3f",track->GetWeight());
01494                         }
01495                         break;
01496                 case NEWLINE:
01497                         printf("\n");
01498                         break;
01499                 }
01500         }
01501 
01502         fflush(stdout); // help debugging location of crashes
01503 }

void BLManager::setSteppingManager ( G4SteppingManager *  p  )  [inline]

setSteppingManager() sets the pointer to the current G4SteppingManager.

References fpSteppingManager.

Referenced by BLManager_UserSteppingAction::UserSteppingAction().

00310 { fpSteppingManager = p; }

G4SteppingManager* BLManager::getSteppingManager (  )  [inline]

getSteppingManager() returns a pointer to the current G4SteppingManager.

References fpSteppingManager.

Referenced by PillboxField::UserSteppingAction(), BLCMDtune::UserZSteppingAction(), and BLCMDreference::UserZSteppingAction().

00314 { return fpSteppingManager; }

void BLManager::setTrackingManager ( G4TrackingManager *  p  )  [inline]

setTrackingManager() sets the pointer to the current G4TrackingManager.

References fpTrackingManager.

Referenced by BLManager_UserTrackingAction::PreUserTrackingAction().

00318 { fpTrackingManager = p; }

void BLManager::initialize (  ) 

initialize() will initialize the BLManager object, and the geant4 kernel, thus constructing the geometry in the world group. Note that registerPhysics() must be called before initialize() (normally done by a "physics" command in the input file).

physics->getPhysicsList() is set in registerPhysicsList()

References beamIndex, beamVector, eventTimeLimit, BLParam::getDouble(), BLParam::getInt(), histoUpdate, initialized, Param, physics, referenceVector, runManager, setSteppingFormat(), steppingVerbose, and zTolerance.

Referenced by BLCMDlist::command(), and main().

00315 {
00316         if(initialized) {
00317                 G4Exception("BLManager","Already Initialized",FatalException,
00318                                                                         "");
00319         }
00320         initialized = true;
00321 
00322         G4UImanager* UI = G4UImanager::GetUIpointer();
00323 
00324         if(!physics) {
00325                 G4Exception("BLManager","No physics registered",FatalException,
00326                                                                         "");
00327         }
00328 
00329         setSteppingFormat();
00330 
00331         // set initialization and user classes
00332         runManager->SetUserInitialization((G4VUserDetectorConstruction*)this);
00333         /// physics->getPhysicsList() is set in registerPhysicsList()
00334         runManager->SetUserAction(new BLManager_UserSteppingAction(this));
00335         runManager->SetUserAction(new BLManager_UserTrackingAction(this));
00336         runManager->SetUserAction(new BLManager_UserRunAction(this));
00337         runManager->SetUserAction(new BLManager_UserEventAction(this));
00338         runManager->SetUserAction(new BLManager_PrimaryGeneratorAction(this));
00339         runManager->SetUserAction(new BLManager_UserStackingAction(this));
00340 
00341         // apply initial range cuts
00342         UI->ApplyCommand("/range/cutG 2 mm");
00343         UI->ApplyCommand("/range/cutE 2 mm");
00344 
00345         // initialize the RunManager, construct the apparatus, etc.
00346         runManager->Initialize();
00347 
00348         // set initial verbosities
00349         UI->ApplyCommand("/control/verbose 0");
00350         UI->ApplyCommand("/run/verbose 0");
00351         UI->ApplyCommand("/event/verbose 0");
00352         UI->ApplyCommand("/tracking/verbose 0");
00353         UI->ApplyCommand("/hits/verbose 0");
00354         UI->ApplyCommand("/material/verbose 0");
00355         UI->ApplyCommand("/process/setVerbose 0 all");
00356         UI->ApplyCommand("/process/verbose 0");
00357         UI->ApplyCommand("/process/eLoss/verbose 0");
00358         G4HadronicProcessStore::Instance()->SetVerbose(0);
00359 
00360         steppingVerbose = Param.getInt("steppingVerbose");
00361 
00362         // get various parameter values (may have changed since constructor)
00363         histoUpdate = Param.getInt("histoUpdate");
00364         eventTimeLimit = Param.getInt("eventTimeLimit");
00365         zTolerance = Param.getDouble("zTolerance");
00366 
00367         // initialize the beam
00368         for(beamIndex=0; beamIndex<beamVector.size(); ++beamIndex)
00369                 beamVector[beamIndex]->init();
00370         for(unsigned i=0; i<referenceVector.size(); ++i)
00371                 referenceVector[i]->init();
00372         beamIndex = 0;
00373 
00374         // create the ZStepLimiter
00375         new ZStepLimiter();
00376 }

static bool BLManager::isInitialized (  )  [inline, static]

isInitialized() returns true if the BLManager has been initialized.

References initialized.

Referenced by BLCMDlist::command().

00328 { return initialized; }

void BLManager::trackTuneAndReferenceParticles (  ) 

trackTuneAndReferenceParticles() will generate and track the tune particle and then the reference particle.

References beamIndex, BLRunManager::BeamOn(), FORCE_OFF, BLRunManager::getCollectiveMode(), IDLE, NORMAL, physics, REFERENCE, referenceVector, runManager, BLRunManager::setCollectiveMode(), BLPhysics::setDoStochastics(), setEventID(), state, and TUNE.

Referenced by main().

00406 {
00407         if(referenceVector.size() == 0)
00408                 return;
00409 
00410         if(!physics) {
00411                 G4Exception("BLManager","No physics registered",FatalException,
00412                                                                         "");
00413         }
00414 
00415         // Tune and Reference particles cannot use collective mode
00416         bool collectiveMode = runManager->getCollectiveMode();
00417         runManager->setCollectiveMode(false);
00418 
00419         printf("================= Prepare Tune Particle(s) ===========\n");
00420         physics->setDoStochastics(FORCE_OFF);
00421         runManager->Initialize();
00422 
00423         printf("================= Begin Tune Particle(s) =============\n");
00424         state = TUNE;
00425         setEventID(-2);
00426         beamIndex = 0;
00427         runManager->BeamOn(referenceVector.size());
00428         state = IDLE;
00429 
00430         // now track center particle
00431         printf("================== Begin Reference Particle(s) ===============\n");
00432         state = REFERENCE;
00433         setEventID(-1);
00434         beamIndex = 0;
00435         runManager->BeamOn(referenceVector.size());
00436         state = IDLE;
00437         beamIndex = 0;
00438 
00439         physics->setDoStochastics(NORMAL);
00440         runManager->setCollectiveMode(collectiveMode);
00441 }

void BLManager::trackBeam (  ) 

trackBeam() will generate the defined beam and track each event.

References BEAM, beamIndex, BLRunManager::BeamOn(), beamVector, IDLE, NORMAL, physics, runManager, BLPhysics::setDoStochastics(), setEventID(), and state.

Referenced by main().

00444 {
00445         if(beamVector.size() == 0) {
00446                 G4Exception("BLManager","No beam registered",FatalException,
00447                                                                         "");
00448         }
00449         if(!physics) {
00450                 G4Exception("BLManager","No physics registered",FatalException,
00451                                                                         "");
00452         }
00453 
00454         printf("================== Prepare Tracking Beam ==================\n");
00455         physics->setDoStochastics(NORMAL);
00456         runManager->Initialize();
00457 
00458         printf("================== Begin Tracking Beam ===============\n");
00459         state = BEAM;
00460         setEventID(1);
00461         beamIndex = 0;
00462         runManager->BeamOn(0x7FFFFFFF);
00463         for(beamIndex=0; beamIndex<beamVector.size(); ++beamIndex)
00464                 beamVector[beamIndex]->summary();
00465         beamIndex = 0;
00466         state = IDLE;
00467 }

void BLManager::displayVisual (  ) 

displayVisual() will display the detector visually.

References beamIndex, BLParam::getString(), handleCallbacks(), IDLE, BLVisManager::init(), NORMAL, Param, physics, runManager, BLPhysics::setDoStochastics(), setEventID(), state, and VISUAL.

Referenced by main().

00470 {
00471 #ifdef G4VIS_USE
00472         printf("================== Prepare Visualization ==================\n");
00473         state = VISUAL;
00474         physics->setDoStochastics(NORMAL);
00475         runManager->Initialize();
00476         setEventID(1);
00477         beamIndex = 0;
00478 
00479         BLVisManager *visManager = new BLVisManager(Param.getString("viewer"));
00480         visManager->init();
00481 
00482         handleCallbacks(4);
00483 
00484         printf("To display a run with 10 events in the viewer, type '/run/beamOn 10<cr>'\n");
00485         G4UIterminal *terminal = new G4UIterminal(0,false);
00486         terminal->SessionStart();
00487         
00488         delete terminal;
00489         delete visManager;
00490         state = IDLE;
00491         beamIndex = 0;
00492 #else
00493         G4Exception("BLManager","G4VIS_USE not defined",FatalException,"");
00494 #endif
00495 }

void BLManager::displayGeometry ( G4VPhysicalVolume *  phys = 0,
int  level = 0 
)

displayGeometry() will display the geant4 geometry. This is a hierarchical ASCII listing of all volumes. if phys==0 then use the worldPhysicalVolume.

References BLCommand::dumpRotation(), and worldPhysicalVolume.

Referenced by BLCMDgeometry::callback().

00498 {
00499         if(phys == 0) phys = worldPhysicalVolume;
00500         for(int i=0; i<level; ++i) printf("  ");
00501         G4ThreeVector pos = phys->GetObjectTranslation();
00502         printf("%s x=%.3f y=%.3f z=%.3f\n",phys->GetName().c_str(),
00503                         pos[0],pos[1],pos[2]);
00504         G4RotationMatrix rot = phys->GetObjectRotationValue();
00505         if(!rot.isIdentity()) {
00506                 for(int i=0; i<level; ++i) printf("  ");
00507                 BLCommand::dumpRotation(&rot,"");
00508         }
00509         G4LogicalVolume *log = phys->GetLogicalVolume();
00510         int n = log->GetNoDaughters();
00511         for(int i=0; i<n; ++i) {
00512                 G4VPhysicalVolume *p = log->GetDaughter(i);
00513                 if(p)
00514                         displayGeometry(p,level+1);
00515         }
00516 }

void BLManager::registerPhysics ( BLPhysics _physics  )  [inline]

registerPhysics() registers the BLPhysics object. It also sets the physics list to the BLRunManager, so following commands can find particle by name.

References BLPhysics::getPhysicsList(), physics, and runManager.

Referenced by BLCMDphysics::command().

00349                 { physics = _physics;
00350                   runManager->SetUserInitialization(physics->getPhysicsList());
00351                 }

BLPhysics* BLManager::getPhysics (  )  [inline]

void BLManager::registerBeam ( BLBeam _beam  )  [inline]

registerBeam() registers a BLBeam object for beam generation. Multiple BLBeam objects can be registered, used in order.

References beamVector.

Referenced by BLCMDparticlesource::command(), BLCMDcosmicraybeam::command(), and BLCMDbeam::command().

00358 { beamVector.push_back(_beam); }

void BLManager::registerReference ( BLBeam _beam  )  [inline]

registerReference() registers a BLBeam object for reference particle generation. Multiple BLBeam objects can be registered, used in order.

References referenceVector.

Referenced by BLCMDreference::command().

00364                 { referenceVector.push_back(_beam); }

int BLManager::nReference (  )  [inline]

nReference() returns the number of reference particles registered.

References referenceVector.

Referenced by main().

00367 { return referenceVector.size(); }

void BLManager::registerRunAction ( RunAction a,
G4bool  beamOnly = true 
) [inline]

registerRunAction() registers a UserRunAction. If beamOnly is true (the default), the callback is made only if state==BEAM.

References beamRunActionVector, and runActionVector.

Referenced by BLCMDtune::command(), BLCMDreference::command(), and BLCMDhelicaldipole::construct().

00373                 { if(beamOnly)
00374                         beamRunActionVector.push_back(a);
00375                   else
00376                         runActionVector.push_back(a);
00377                 }

void BLManager::registerEventAction ( EventAction a,
G4bool  beamOnly = true 
) [inline]

registerEventAction() registers a UserEventAction. If beamOnly is true (the default), the callback is made only if state==BEAM.

References beamEventActionVector, and eventActionVector.

Referenced by BLCMDparticlecolor::command(), and BLCMDntuple::command().

00383                 { if(beamOnly)
00384                         beamEventActionVector.push_back(a);
00385                   else
00386                         eventActionVector.push_back(a);
00387                 }

void BLManager::registerTrackingAction ( TrackingAction a  )  [inline]

void BLManager::registerUserCode ( BLUserCode instance  )  [inline]

registerUserCode() registers a BLUserCode instance.

References userCodeVector.

Referenced by BLUserCode::registerUserCode().

00395                 { userCodeVector.push_back(instance); }

std::vector<BLUserCode*> BLManager::getUserCodeInstances ( G4String  type  )  [inline]

getUserCodeInstances() returns a vector of all registered instances of BLUserCode with the specified type.

References userCodeVector.

00400                 { std::vector<BLUserCode*> ret;
00401                   for(unsigned i=0; i<userCodeVector.size(); ++i) {
00402                         if(type == userCodeVector[i]->getType())
00403                                 ret.push_back(userCodeVector[i]);
00404                   }
00405                   return ret;
00406                 }

void BLManager::registerCallback ( BLCallback cb,
int  type 
) [inline]

registerCallback() registers a BLCallback. type=0 for pre-Tune particle, type=1 for post-Reference (pre-beam tracking), type=2 for post-beam tracking. type=3 for replacing the main program loop. type=4 for visualization. NOTE: if there are type=3 callbacks, when the last one returns the closes up by summarizing NTuples and callng handleCallbacks(2), and then the program exits. This prevents the main program loop from executing. type=3 callbacks are called just after the type=1 callbacks (i.e. after the post-Reference callbacks).

References postReferenceCallbackVector, postTrackingCallbackVector, preReferenceCallbackVector, replaceMainLoopCallbackVector, and visualizationCallbackVector.

Referenced by BLCMDzntuple::command(), BLCMDusertrackfilter::command(), BLCMDtrackermode::command(), BLCMDtracker::command(), BLCMDtotalenergy::command(), BLCMDshowmaterial::command(), BLCMDsetdecay::command(), BLCMDprofile::command(), BLCMDprobefield::command(), BLCMDprintfield::command(), BLCMDprintf::command(), BLCMDphysics::command(), BLCMDmovie::command(), BLCMDgeometry::command(), BLCMDg4ui::command(), BLCMDfieldntuple::command(), BLCMDbug1021::command(), SetDecayInstance::fatalError(), ParticleFilterPlacement::ParticleFilterPlacement(), and Surface::Surface().

00420                                                         {
00421                 if(type==0) preReferenceCallbackVector.push_back(cb);
00422                 else if(type==1) postReferenceCallbackVector.push_back(cb);
00423                 else if(type==2) postTrackingCallbackVector.push_back(cb);
00424                 else if(type==3) replaceMainLoopCallbackVector.push_back(cb);
00425                 else if(type==4) visualizationCallbackVector.push_back(cb);
00426         }

void BLManager::handleCallbacks ( int  type  ) 

handleCallbacks() calls all applicable registered callbacks. type=0 for pre-reference particle, type=1 for post-center (pre-beam tracking), type=2 for post-beam tracking (just before program exit). type=3 for replacing the main program loop. type=4 for visualization.

References blManager, BLRunManager::callback(), BLNTuple::closeAll(), exceptionSummary(), g4bl_exit(), BLRunManager::getObject(), postReferenceCallbackVector, postTrackingCallbackVector, preReferenceCallbackVector, replaceMainLoopCallbackVector, BLNTuple::summary(), and visualizationCallbackVector.

Referenced by displayVisual(), main(), and Notify().

00530 {
00531         // startup sequencing prevents RunManager from registering
00532         BLRunManager::getObject()->callback(type);
00533 
00534         if(type == 0) {
00535                 for(unsigned i=0; i<preReferenceCallbackVector.size(); ++i)
00536                         preReferenceCallbackVector[i]->callback(type);
00537         } else if(type == 1) {
00538                 for(unsigned i=0; i<postReferenceCallbackVector.size(); ++i)
00539                         postReferenceCallbackVector[i]->callback(type);
00540         } else if(type == 2) {
00541                 for(unsigned i=0; i<postTrackingCallbackVector.size(); ++i)
00542                         postTrackingCallbackVector[i]->callback(type);
00543                 exceptionSummary();
00544         } else if(type == 3) {
00545                 for(unsigned i=0; i<replaceMainLoopCallbackVector.size(); ++i)
00546                         replaceMainLoopCallbackVector[i]->callback(type);
00547                 if(replaceMainLoopCallbackVector.size() > 0) {
00548                         // we have replaced the main loop, so closeup and exit
00549                         BLNTuple::summary();
00550                         BLNTuple::closeAll();
00551                         // handle post-tracking callbacks
00552                         handleCallbacks(2);
00553                         // Instead of deleting things, just open the geometry
00554                         // to avoid warnings when exit() is called.
00555                         G4GeometryManager::GetInstance()->OpenGeometry();
00556                         blManager = 0;
00557                         extern void g4bl_exit(int);
00558                         g4bl_exit(0);
00559                 }
00560         } else if(type == 4) {
00561                 for(unsigned i=0; i<visualizationCallbackVector.size(); ++i)
00562                         visualizationCallbackVector[i]->callback(type);
00563         }
00564 }

G4VPhysicalVolume* BLManager::getWorldPhysicalVolume (  )  [inline]

getWorldPhysicalVolume() returns a pointer to the world PV. Note that it must already have been consructed, so this function returns NULL before construct() is called (by main).

References worldPhysicalVolume.

Referenced by BLCMDtotalenergy::callback(), BLCMDshowmaterial::callback(), and BLCMDmovie::callback().

00440                 { return worldPhysicalVolume; }

void BLManager::setPRNGSeedMethod ( PRNGSeedMethod  method  )  [inline]

setPRNGSeedMethod() will set the method used to seed the pseudo random number generator at the start of each event.

References seedMethod.

Referenced by BLCMDrandomseed::command(), and BLBeam::setRandomSeedToGenerate().

00445                 { seedMethod = method; }

PRNGSeedMethod BLManager::getPRNGSeedMethod (  )  [inline]

getPRNGSeedMethod() will return the method used to seed the pseudo random number generator at the start of each event.

References seedMethod.

00449 { return seedMethod; }

void BLManager::setEventCutFile ( G4String &  f  )  [inline]

setEventCutFile() will set the file containing event s for cutting. The file is ASCII, with one event number per line; lines beginning with # are ignored.

References eventCutFile.

Referenced by BLCMDeventcuts::command().

00454 { eventCutFile = f; }

void BLManager::UserSteppingAction ( const G4Step *  step  ) 

UserSteppingAction() from G4UserSteppingAction.

References allStepMap, allStepVector, BEAM, beamStepMap, beamStepVector, currentZStep, eventTimeLimit, BLCoordinates::getCLZ(), indexZStep, BLCoordinates::isValid(), nStuck, prevZ, REFERENCE, rpStepMap, rpStepVector, ZStepLimiter::setMaxStep(), startEvent, state, steppingVerbose, steppingVerbosePrint(), BLTime::time(), tpStepMap, tpStepVector, tracksAborted, TUNE, BLCoordinates::update(), and zTolerance.

Referenced by BLManager_UserSteppingAction::UserSteppingAction().

00768 {
00769         G4Track *track = step->GetTrack();
00770         G4StepPoint *prePoint = step->GetPreStepPoint();
00771         G4StepPoint *postPoint = step->GetPostStepPoint();
00772         G4VPhysicalVolume *preVol=0;
00773         G4VPhysicalVolume *postVol=0;
00774         if(prePoint) preVol = prePoint->GetPhysicalVolume();
00775         if(postPoint) postVol = postPoint->GetPhysicalVolume();
00776 
00777         // need to update BLCoordinates before steppingVerbose print.
00778         // So BLCoordinates is not registered with SteppingAction, it
00779         // is handled specially, right here.
00780         BLCoordinates::update(track);
00781         BLCoordinates *coord = (BLCoordinates *)track->GetUserInformation();
00782         if(!coord || !coord->isValid()) {
00783                 G4Exception("BLManager","Coordinates Got Lost",FatalException,
00784                                                                         "");
00785         }
00786 
00787         // steppingVerbose print moved to the top of this routine
00788         if(steppingVerbose > 0) { 
00789                 static bool first=true;
00790                 int nstep = track->GetCurrentStepNumber();
00791                 if(first || (nstep <= 1 && !runManager->getCollectiveMode()))
00792                         steppingVerbosePrint(0,0,1);
00793                 steppingVerbosePrint(step,track,0);
00794                 first = false;
00795         }
00796 
00797         /* check for stuck track */
00798         if(step->GetStepLength() > 0.0001*mm ||
00799            step->GetDeltaTime() > 0.0001*ns) {
00800                 nStuck = 0;
00801         } else {
00802                 if(++nStuck >= 100) {
00803                     G4Exception("BLManager Stepping Action",
00804                     "Stuck Track -- Killed",JustWarning,
00805                     "100 steps in a row, each less than 0.1 micron and 0.1 ps");
00806                     ++tracksAborted;
00807                     track->SetTrackStatus(fStopAndKill);
00808                 }
00809         }
00810 
00811         G4TrackStatus status = track->GetTrackStatus();
00812         if(status == fStopAndKill || status == fKillTrackAndSecondaries)
00813                 goto quit;
00814 
00815         // check event time limit
00816         if(eventTimeLimit > 0 && BLTime::time()-startEvent > eventTimeLimit) {
00817                 G4Exception("BLManager","Event Time Limit",EventMustBeAborted,"");
00818                 goto quit;
00819         }
00820 
00821         // call ZStep actions
00822         // NOTE: indexZstep points to the next entry in the +Z direction.
00823         // Note also the bookends constrain indexZStep to always be valid.
00824         if(currentZStep->size() > 2) { // i.e. not just 2 bookends
00825                 G4double thisZ = coord->getCLZ();
00826                 if(thisZ == prevZ) goto noZstep;
00827                 G4double minZ = (thisZ<prevZ ? thisZ : prevZ);
00828                 G4double maxZ = (thisZ>prevZ ? thisZ : prevZ);
00829                 // find first entry that might be spanned by the step
00830                 while(indexZStep > 0 && minZ <= (*currentZStep)[indexZStep-1].z)
00831                         --indexZStep;
00832                 while(minZ > (*currentZStep)[indexZStep].z)
00833                         ++indexZStep;
00834                 int indexPrev = indexZStep-1;
00835                 // loop over all entries actually spanned by this step
00836                 for( ; indexZStep<currentZStep->size()-1; ++indexZStep) {
00837                         G4double z=(*currentZStep)[indexZStep].z;
00838                         if(maxZ <= z) break;
00839                         if((*currentZStep)[indexZStep].action == 0) continue;
00840                         // save current values from track
00841                         G4ThreeVector pos = track->GetPosition();
00842                         G4ThreeVector mom = track->GetMomentum();
00843                         G4double time = track->GetGlobalTime();
00844                         G4double ke = track->GetKineticEnergy();
00845                         // interpolate linearly to z, global coords
00846                         G4ThreeVector deltaPos=step->GetDeltaPosition();
00847                         G4ThreeVector deltaMom=step->GetDeltaMomentum();
00848                         G4double deltaTime=step->GetDeltaTime();
00849                         G4double deltaE=step->GetDeltaEnergy();
00850                         G4double f=(z-thisZ)/(prevZ-thisZ);
00851                         track->SetPosition(pos-f*deltaPos);
00852                         track->SetMomentumDirection((mom-f*deltaMom).unit());
00853                         track->SetGlobalTime(time-f*deltaTime);
00854                         track->SetKineticEnergy(ke-f*deltaE);
00855                         BLCoordinates::update(track);
00856                         (*currentZStep)[indexZStep].action->
00857                                                 UserZSteppingAction(track);
00858                         // restore current values to track
00859                         track->SetPosition(pos);
00860                         track->SetMomentumDirection(mom.unit());
00861                         track->SetGlobalTime(time);
00862                         track->SetKineticEnergy(ke);
00863                         BLCoordinates::update(track);
00864                 }
00865                 double dz = (thisZ>prevZ ? (*currentZStep)[indexZStep].z-thisZ : 
00866                                            thisZ-(*currentZStep)[indexPrev].z);
00867                 dz = (dz>zTolerance*2.0 ? dz : zTolerance*2.0);
00868                 ZStepLimiter::setMaxStep(dz);
00869                 prevZ = thisZ;
00870         }
00871 noZstep:
00872 
00873         // call all-state stepping actions before per-state actions.
00874         { std::vector<BLManager::SteppingAction*>::iterator i;
00875           for(i=allStepVector.begin(); i!=allStepVector.end(); ++i) {
00876                   (*i)->UserSteppingAction(step);
00877                   if(track->GetTrackStatus() != fAlive)
00878                           goto quit;
00879           }
00880           if(preVol && allStepMap.count(preVol) > 0) {
00881                 allStepMap[preVol]->UserSteppingAction(step);
00882                 if(track->GetTrackStatus() != fAlive)
00883                         goto quit;
00884           }
00885           if(postVol  && preVol != postVol && allStepMap.count(postVol) > 0) {
00886                 allStepMap[postVol]->UserSteppingAction(step);
00887                 if(track->GetTrackStatus() != fAlive)
00888                         goto quit;
00889           }
00890         }
00891 
00892         // call Tune Particle stepping actions
00893         if(state == TUNE) {
00894                 std::vector<BLManager::SteppingAction*>::iterator i;
00895                 for(i=tpStepVector.begin(); i!=tpStepVector.end(); ++i) {
00896                         (*i)->UserSteppingAction(step);
00897                         if(track->GetTrackStatus() != fAlive)
00898                                 goto quit;
00899                 }
00900                 if(preVol && tpStepMap.count(preVol) > 0) {
00901                         tpStepMap[preVol]->UserSteppingAction(step);
00902                         if(track->GetTrackStatus() != fAlive)
00903                                 goto quit;
00904                 }
00905                 if(postVol && preVol != postVol && 
00906                                                 tpStepMap.count(postVol) > 0) {
00907                         tpStepMap[postVol]->UserSteppingAction(step);
00908                         if(track->GetTrackStatus() != fAlive)
00909                                 goto quit;
00910                 }
00911         }
00912 
00913         // call Reference Particle stepping actions
00914         if(state == REFERENCE) {
00915                 std::vector<BLManager::SteppingAction*>::iterator i;
00916                 for(i=rpStepVector.begin(); i!=rpStepVector.end(); ++i) {
00917                         (*i)->UserSteppingAction(step);
00918                         if(track->GetTrackStatus() != fAlive)
00919                                 goto quit;
00920                 }
00921                 if(preVol && rpStepMap.count(preVol) > 0) {
00922                         rpStepMap[preVol]->UserSteppingAction(step);
00923                         if(track->GetTrackStatus() != fAlive)
00924                                 goto quit;
00925                 }
00926                 if(postVol && preVol != postVol && 
00927                                                 rpStepMap.count(postVol) > 0) {
00928                         rpStepMap[postVol]->UserSteppingAction(step);
00929                         if(track->GetTrackStatus() != fAlive)
00930                                 goto quit;
00931                 }
00932         }
00933 
00934         // call beam stepping actions
00935         if(state == BEAM) {
00936                 std::vector<BLManager::SteppingAction*>::iterator i;
00937                 for(i=beamStepVector.begin(); i!=beamStepVector.end(); ++i) {
00938                         (*i)->UserSteppingAction(step);
00939                         if(track->GetTrackStatus() != fAlive)
00940                                 goto quit;
00941                 }
00942                 if(preVol && beamStepMap.count(preVol) > 0) {
00943                         beamStepMap[preVol]->UserSteppingAction(step);
00944                         if(track->GetTrackStatus() != fAlive)
00945                                 goto quit;
00946                 }
00947                 if(postVol  && preVol != postVol && 
00948                                              beamStepMap.count(postVol) > 0) {
00949                         beamStepMap[postVol]->UserSteppingAction(step);
00950                         if(track->GetTrackStatus() != fAlive)
00951                                 goto quit;
00952                 }
00953         }
00954 
00955 quit:   ;
00956 }

G4VPhysicalVolume * BLManager::Construct (  ) 

Construct() from G4VUserDetectorConstruction.

References BLGroup::constructWorld(), BLGlobalField::getObject(), and worldPhysicalVolume.

00519 {
00520         // ensure the global field is initialized
00521         (void)BLGlobalField::getObject();
00522 
00523         // construct the world
00524         worldPhysicalVolume =  BLGroup::constructWorld();
00525 
00526         return worldPhysicalVolume;
00527 }

void BLManager::PreUserTrackingAction ( const G4Track *  track  ) 

PreUserTrackingAction() from G4UserTrackingAction.

References beamZStep, currentZStep, eventID, fpTrackingManager, BLCoordinates::getCLZ(), BLRunManager::getCollectiveMode(), getExternalParentID(), BLTrackInfo::getExternalTrackID(), indexZStep, BLCoordinates::isValid(), nextSecondaryTrackID, nStuck, prevZ, primaryParentID, primaryTrackID, REFERENCE, referenceZStep, runManager, BLTrackInfo::setExternalParentID(), BLTrackInfo::setExternalTrackID(), BLCoordinates::setGlobal(), ZStepLimiter::setMaxStep(), state, steppingVerbose, trackID, trackIDMap, trackingActionVector, trackParticleDef, TUNE, tuneZStep, and zTolerance.

Referenced by BLManager_UserTrackingAction::PreUserTrackingAction().

00669 {
00670         static bool first=true;
00671         if(first) {
00672                 first=false;
00673                 std::vector<BLManager::TrackingAction*>::iterator i;
00674                 for(i=trackingActionVector.begin(); 
00675                                          i<trackingActionVector.end(); ++i) {
00676                         (*i)->SetTrackingManagerPointer(fpTrackingManager);
00677                 }
00678         }
00679 
00680         nStuck = 0;
00681 
00682         // link a BLTrackInfo into the track (called coord -- historical)
00683         BLTrackInfo *coord = (BLTrackInfo *)track->GetUserInformation();
00684         if(!coord || !coord->isValid()) {
00685                 coord = new BLTrackInfo();
00686                 ((G4Track *)track)->SetUserInformation(coord);
00687         }
00688         coord->setGlobal(track->GetPosition(),track->GetGlobalTime());
00689 
00690         // set the external trackID and parentID
00691         trackID = coord->getExternalTrackID();
00692         if(trackID <= 0) {
00693                 if(primaryTrackID >= 0) {
00694                         trackID = primaryTrackID;
00695                         coord->setExternalTrackID(trackID);
00696                         coord->setExternalParentID(primaryParentID);
00697                         primaryTrackID = primaryParentID = -1;
00698                 } else {
00699                         trackID = nextSecondaryTrackID++;
00700                         coord->setExternalTrackID(trackID);
00701                         coord->setExternalParentID(
00702                                         trackIDMap[track->GetParentID()]);
00703                 }
00704         }
00705         trackIDMap[track->GetTrackID()] = trackID;
00706         trackParticleDef = track->GetDefinition();
00707 
00708         // print header, if appropriate
00709         if(steppingVerbose && !runManager->getCollectiveMode()) {
00710             printf("=========== Event %d Track %d %s    Parent %d",
00711                         eventID,trackID,
00712                         trackParticleDef->GetParticleName().c_str(),
00713                         getExternalParentID(track));
00714             if(getExternalParentID(track) == 0) {
00715                 printf("   CreatorProcess=Beam");
00716             } else {
00717                 const G4VProcess *proc = track->GetCreatorProcess();
00718                 if(proc != 0)
00719                         printf("   CreatorProcess=%s",
00720                                         proc->GetProcessName().c_str());
00721             }
00722             printf(" ===========\n");
00723         }
00724 
00725         // set currentZStep
00726         if(state == TUNE)
00727                 currentZStep = &tuneZStep;
00728         else if(state == REFERENCE)
00729                 currentZStep = &referenceZStep;
00730         else
00731                 currentZStep = &beamZStep;
00732         indexZStep = 1;
00733         prevZ = coord->getCLZ();
00734         while(prevZ > (*currentZStep)[indexZStep].z) {
00735                 ++indexZStep;
00736         }
00737         assert(indexZStep < currentZStep->size());
00738         G4double dz = (*currentZStep)[indexZStep].z - prevZ;
00739         if(dz < zTolerance) dz = zTolerance;
00740         ZStepLimiter::setMaxStep(dz);
00741 
00742         // loop over all registered user tracking actions
00743         std::vector<BLManager::TrackingAction*>::iterator i;
00744         for(i=trackingActionVector.begin(); i<trackingActionVector.end(); ++i) {
00745                 (*i)->PreUserTrackingAction(track);
00746         }
00747 }

void BLManager::PostUserTrackingAction ( const G4Track *  track  ) 

PostUserTrackingAction() from G4UserTrackingAction.

References BLRunManager::getCollectiveMode(), runManager, trackID, trackingActionVector, and trackParticleDef.

Referenced by BLManager_UserTrackingAction::PostUserTrackingAction().

00750 {
00751         std::vector<BLManager::TrackingAction*>::iterator i;
00752         for(i=trackingActionVector.begin(); i<trackingActionVector.end(); ++i) {
00753                 (*i)->PostUserTrackingAction(track);
00754         }
00755 
00756         // delete the BLcoordinates
00757         G4VUserTrackInformation *ti = track->GetUserInformation();
00758         if(ti && !runManager->getCollectiveMode()) {
00759                 delete ti;
00760                 ((G4Track *)track)->SetUserInformation(0);
00761         }
00762 
00763         trackID = -9999;
00764         trackParticleDef = 0;
00765 }

void BLManager::BeginOfRunAction ( const G4Run *  run  ) 

BeginOfRunAction() from G4UserRunAction.

References BEAM, beamRunActionVector, endRun, eventsProcessed, readEventCutFile(), runActionVector, startEvent, startRun, state, and BLTime::time().

Referenced by BLManager_UserRunAction::BeginOfRunAction().

00567 {
00568         readEventCutFile();
00569 
00570         startEvent = startRun = BLTime::time();
00571         eventsProcessed = 0;
00572         endRun = false;
00573 
00574         G4UImanager* UI = G4UImanager::GetUIpointer();
00575         UI->ApplyCommand("/control/verbose 0");
00576         UI->ApplyCommand("/run/verbose 0");
00577         UI->ApplyCommand("/event/verbose 0");
00578         UI->ApplyCommand("/tracking/verbose 0");
00579         UI->ApplyCommand("/hits/verbose 0");
00580         UI->ApplyCommand("/material/verbose 0");
00581         UI->ApplyCommand("/process/setVerbose 0 all");
00582         UI->ApplyCommand("/process/verbose 0");
00583         UI->ApplyCommand("/process/eLoss/verbose 0");
00584 
00585         for(unsigned int i=0; i<runActionVector.size(); ++i)
00586                 runActionVector[i]->BeginOfRunAction(run);
00587 
00588         if(state == BEAM) {
00589                 for(unsigned int i=0; i<beamRunActionVector.size(); ++i)
00590                         beamRunActionVector[i]->BeginOfRunAction(run);
00591         }
00592 }

void BLManager::EndOfRunAction ( const G4Run *  run  ) 

EndOfRunAction() from G4UserRunAction.

References BEAM, beamRunActionVector, eventsProcessed, runActionVector, startRun, state, and BLTime::time().

Referenced by BLManager_UserRunAction::EndOfRunAction().

00595 {
00596         if(state == BEAM) {
00597                 for(unsigned int i=0; i<beamRunActionVector.size(); ++i)
00598                         beamRunActionVector[i]->EndOfRunAction(run);
00599         }
00600 
00601         for(unsigned int i=0; i<runActionVector.size(); ++i)
00602                 runActionVector[i]->EndOfRunAction(run);
00603 
00604         printf("Run complete  %d Events  %ld seconds\n",
00605                         eventsProcessed,(long)(BLTime::time()-startRun));
00606 }

void BLManager::BeginOfEventAction ( const G4Event *  event  ) 

BeginOfEventAction() from G4UserEventAction.

References BEAM, beamEventActionVector, endRun, eventActionVector, eventTimeLimit, BLRunManager::getCollectiveMode(), runManager, BLAlarm::set(), startEvent, state, steppingVerbose, and BLTime::time().

Referenced by BLRunManager::beamOnCollective(), BLRunManager::beginEvent(), and BLManager_UserEventAction::BeginOfEventAction().

00609 {
00610         if(endRun || event->IsAborted()) return;
00611 
00612         if(state == BEAM) {
00613 
00614             int evId = event->GetEventID();
00615             if(steppingVerbose && !runManager->getCollectiveMode()) 
00616                 printf("\n\n=================== Event %d ==================\n",
00617                                                                 evId);
00618 
00619             // call all registered action-s
00620             std::vector<BLManager::EventAction*>::iterator i;
00621             for(i=beamEventActionVector.begin(); i<beamEventActionVector.end(); ++i) {
00622                 (*i)->BeginOfEventAction(event);
00623             }
00624         }
00625 
00626         // call all registered action-s
00627         std::vector<BLManager::EventAction*>::iterator i;
00628         for(i=eventActionVector.begin(); i<eventActionVector.end(); ++i) {
00629                 (*i)->BeginOfEventAction(event);
00630         }
00631 
00632         startEvent = BLTime::time();
00633         // add 2 seconds so the test in UserSteppingAction() can try first
00634         // (it kills 1 event, BLAlarm kills the entire job).
00635         if(eventTimeLimit > 0) BLAlarm::set(eventTimeLimit+2);
00636 }

void BLManager::EndOfEventAction ( const G4Event *  event  ) 

EndOfEventAction() from G4UserEventAction.

References BEAM, beamEventActionVector, BLAlarm::clear(), endRun, eventActionVector, eventsProcessed, BLNTuple::flushAll(), histoUpdate, incrEventsProcessed(), state, and BLNTuple::update().

Referenced by BLRunManager::endEvent(), and BLManager_UserEventAction::EndOfEventAction().

00639 {
00640         BLAlarm::clear();
00641 
00642         if(endRun || event->IsAborted()) return;
00643 
00644         // call all registered action-s
00645         std::vector<BLManager::EventAction*>::iterator i;
00646         for(i=eventActionVector.begin(); i<eventActionVector.end(); ++i) {
00647                 (*i)->EndOfEventAction(event);
00648         }
00649 
00650         if(state != BEAM) return;
00651 
00652         // call all registered action-s
00653         for(i=beamEventActionVector.begin(); i<beamEventActionVector.end(); ++i) {
00654                 (*i)->EndOfEventAction(event);
00655         }
00656         BLNTuple::update();
00657 
00658         if(endRun) return;
00659 
00660         incrEventsProcessed(event->GetEventID());
00661 
00662         // update histogram file, if appropriate
00663         if(histoUpdate > 0 && eventsProcessed%histoUpdate == 0 && 
00664            eventsProcessed > 0)
00665                 BLNTuple::flushAll();
00666 }

void BLManager::GeneratePrimaries ( G4Event *  event  ) 

GeneratePrimaries() from G4VUserPrimaryGeneratorAction.

References BEAM, beamIndex, beamVector, endRun, eventID, IDLE, prevEventID, REFERENCE, referenceVector, runManager, setEventID(), skipEvent(), SPECIAL, state, TUNE, and VISUAL.

Referenced by BLCMDtrackermode::callback(), and BLManager_PrimaryGeneratorAction::GeneratePrimaries().

00959 {
00960         if(beamVector.size() == 0) {
00961                 G4Exception("BLManager","No beam registered",FatalException,
00962                                                                         "");
00963         }
00964 
00965         switch(state) {
00966         case IDLE:
00967         case SPECIAL:
00968                 G4Exception("BLManager","Erroneous call to GeneratePrimaries",
00969                                                         FatalException,"");
00970         case TUNE:
00971                 setEventID(-2);
00972                 event->SetEventID(-2);
00973                 while(beamIndex < referenceVector.size()) {
00974                         if(referenceVector[beamIndex++]->
00975                                         generateReferenceParticle(event))
00976                                 return;
00977                 }
00978                 goto end_run;
00979         case REFERENCE:
00980                 setEventID(-1);
00981                 event->SetEventID(-1);
00982                 while(beamIndex < referenceVector.size()) {
00983                         if(referenceVector[beamIndex++]->
00984                                         generateReferenceParticle(event))
00985                                 return;
00986                 }
00987                 goto end_run;
00988         case VISUAL:
00989         case BEAM:
00990                 for(;;) {
00991                         // use tempEvent so we can skip events when necessary
00992                         G4Event *tempEvent = new G4Event(*event);
00993                         // default event # (nextBeamEvent() can change it)
00994                         eventID = prevEventID + 1;
00995                         event->SetEventID(eventID);
00996                         // generate the tempEvent, indexing through beamVector[]
00997                         if(beamIndex >= beamVector.size()) {
00998                                 goto end_run;
00999                         }
01000                         while(!beamVector[beamIndex]->nextBeamEvent(tempEvent)){
01001                                 if(++beamIndex >= beamVector.size())
01002                                         goto end_run;
01003                         }
01004                         // update eventID (may have changed)
01005                         eventID = tempEvent->GetEventID();
01006                         prevEventID = eventID; // (this increments it, too)
01007                         // skip events when necessary
01008                         if(!skipEvent(eventID)) {
01009                                 *event = *tempEvent;
01010                                 delete tempEvent;
01011                                 break;
01012                         }
01013                         delete tempEvent;
01014                 }
01015                 break;
01016         }
01017         return;
01018 end_run:
01019         // RunManager cannot abort the event from inside
01020         // UserGeneratePrimaries(), so we do a soft abort
01021         // to the RunManager, and abort the event ourself.
01022         // The result is the same as a hard abort.
01023         runManager->AbortRun(true);
01024         event->SetEventAborted();
01025         beamIndex = 0;
01026         endRun = true;
01027 }

G4ClassificationOfNewTrack BLManager::ClassifyNewTrack ( const G4Track *  track  ) 

ClassifyNewTrack() from G4UserStackingAction.

References stackingActionVector.

Referenced by BLManager_UserStackingAction::ClassifyNewTrack().

01030 {
01031         std::vector<BLManager::StackingAction*>::iterator i;
01032         for(i=stackingActionVector.begin(); i<stackingActionVector.end(); ++i) {
01033                 G4ClassificationOfNewTrack c = (*i)->ClassifyNewTrack(track);
01034                 if(c == fKill) return c;
01035         }
01036         return fUrgent;
01037 }

void BLManager::NewStage (  ) 

NewStage() from G4userStackingAction.

References stackingActionVector.

Referenced by BLManager_UserStackingAction::NewStage().

01040 {
01041         std::vector<BLManager::StackingAction*>::iterator i;
01042         for(i=stackingActionVector.begin(); i<stackingActionVector.end(); ++i) {
01043                 (*i)->NewStage();
01044         }
01045 }

void BLManager::PrepareNewEvent (  ) 

PrepareNewEvent() from G4StackingAction.

References stackingActionVector.

Referenced by BLManager_UserStackingAction::PrepareNewEvent().

01048 {
01049         std::vector<BLManager::StackingAction*>::iterator i;
01050         for(i=stackingActionVector.begin(); i<stackingActionVector.end(); ++i) {
01051                 (*i)->PrepareNewEvent();
01052         }
01053 }

void BLManager::clearTrackIDMap (  )  [inline]

void BLManager::setNextSecondaryTrackID ( int  next  )  [inline]

setNextSecondaryTrackID() sets the external TrackID for the next secondary track. Automatically incremented for subsequent secondaries.

References nextSecondaryTrackID.

Referenced by BLCMDreference::generateReferenceParticle(), BLCMDparticlesource::nextBeamEvent(), and BLCMDbeam::nextBeamEvent().

00504 { nextSecondaryTrackID = next; }

int BLManager::getNextSecondaryTrackID (  )  [inline]

References nextSecondaryTrackID.

Referenced by BLRunManager::beamOnCollective().

00505 { return nextSecondaryTrackID; }

int BLManager::getExternalTrackID ( const G4Track *  track  ) 

int BLManager::getExternalParentID ( const G4Track *  track  ) 

void BLManager::setExternalTrackID ( G4Track *  track,
int  trackID,
int  parentID 
)

setExternalTrackID() will set the external trackID and parentID. if trackID<0 uses nextSecondarTrackID++.

References BLCoordinates::isValid(), nextSecondaryTrackID, BLTrackInfo::setExternalParentID(), and BLTrackInfo::setExternalTrackID().

01074 {
01075         if(trackID < 0) trackID = nextSecondaryTrackID++;
01076 
01077         BLTrackInfo *ti = (BLTrackInfo *)track->GetUserInformation();
01078         if(!ti || !ti->isValid()) {
01079                 ti = new BLTrackInfo();
01080                 track->SetUserInformation(ti);
01081         }
01082         ti->setExternalTrackID(trackID);
01083         ti->setExternalParentID(parentID);
01084 }

int BLManager::getPrimaryTrackID (  )  [inline]

getPrimaryTrackID() returns the primaryTrackID set by the beam command

References primaryTrackID.

Referenced by BLRunManager::beamOnCollective(), and BLRunManager::getNextBeamEventAndTrack().

00522 { return primaryTrackID; }

int BLManager::getPrimaryParentID (  )  [inline]

getPrimaryParentID() returns the primaryParentID set by the beam command

References primaryParentID.

Referenced by BLRunManager::getNextBeamEventAndTrack().

00525 { return primaryParentID; }

void BLManager::setPrimaryTrackID ( int  t,
int  p 
) [inline]

G4bool BLManager::Notify ( const char *  originOfException,
const char *  exceptionCode,
G4ExceptionSeverity  severity,
const char *  description 
)

Notify() from G4VExceptionHandler.

References BLNTuple::closeAll(), eventID, eventsAborted, exceptionCount, fatalExceptions, g4bl_exit(), handleCallbacks(), BLMPI::isMPI(), runManager, BLNTuple::summary(), trackID, trackParticleDef, and warnings.

01109 {
01110         FILE *out=(BLMPI::isMPI() ? stderr : stdout);   // originally was stderr
01111         fflush(stdout);
01112 
01113         G4bool toBeAborted = false;
01114         const char *p="UNKNOWN";
01115         switch(severity) {
01116         case FatalException:
01117                 p = "Fatal Exception";          
01118                 toBeAborted = true;
01119                 ++fatalExceptions;
01120                 break;
01121         case FatalErrorInArgument:
01122                 p = "Fatal Error in Argument";  
01123                 toBeAborted = true;
01124                 ++fatalExceptions;
01125                 break;
01126         case RunMustBeAborted:
01127                 p = "Run Must Be Aborted";
01128                 runManager->AbortRun(false);
01129                 toBeAborted = false;
01130                 ++fatalExceptions;
01131                 break;
01132         case EventMustBeAborted:
01133                 p = "Event must be Aborted";
01134                 runManager->AbortEvent();
01135                 toBeAborted = false;
01136                 ++eventsAborted;
01137                 break;
01138         default:
01139                 p = "Warning";
01140                 toBeAborted = false;
01141                 ++warnings;
01142                 break;
01143         }
01144 
01145         // thin out the printing of many similar exceptions
01146         int n = ++exceptionCount[exceptionCode];
01147         int interval = 1;
01148         if(n >= 1000)
01149                 interval = 1000;
01150         else if(n >= 100)
01151                 interval = 100;
01152         else if(n >= 10)
01153                 interval = 10;
01154         if(n%interval != 0)
01155                 return toBeAborted;
01156 
01157         fprintf(out,"**************************************************************************\n");
01158         fprintf(out,"*** G4Exception: %s\n",exceptionCode);
01159         fprintf(out,"***    severity: %s\n",p);
01160         fprintf(out,"***   issued by: %s\n",originOfException);
01161         if(strlen(description) > 0)
01162                 fprintf(out,"*** description: %s\n",description);
01163         if(strstr(description,"Missing mandatory data")) {
01164                 // brief description already printed, just URL is missing
01165                 const char *p = getenv("G4BL_DIR");
01166                 if(!p) p = "Installation_Directory";
01167                 fprintf(out,
01168                     "***              The URL is http://geant4.cern.ch\n"
01169                     "***              If you unpack it into %s/data\n"
01170                     "***              then the g4bl and g4blgui scripts will"
01171                     " automatically find it.\n",p);
01172         } else if(G4StateManager::GetStateManager()->GetCurrentState() ==
01173                                 G4State_EventProc && trackParticleDef != 0) {
01174                 fprintf(out,"***     EventID: %d     TrackID: %d   %s\n",
01175                                 eventID,trackID,
01176                                 trackParticleDef->GetParticleName().c_str());
01177         }
01178         if(interval > 1)
01179                 fprintf(out,"***    printing: every %d-th occurrence\n",
01180                         interval);
01181         fprintf(out,"**************************************************************************\n");
01182 
01183         fflush(out); // (could be stdout)
01184 
01185         if(toBeAborted) {
01186                 fprintf(out,"g4beamline: attempting to close up after fatal G4Exception...\n");
01187                 BLNTuple::summary();
01188                 BLNTuple::closeAll();
01189                 handleCallbacks(2);
01190                 fprintf(out,"g4beamline: simulation aborted\n");
01191                 extern void g4bl_exit(int);
01192                 g4bl_exit(99);
01193         }
01194 
01195         return toBeAborted;
01196 }

void BLManager::exceptionSummary (  ) 

exceptionSummary() prints a summary of all exceptions

References eventsAborted, exceptionCount, fatalExceptions, tracksAborted, and warnings.

Referenced by handleCallbacks().

01199 {
01200         printf("\nException Summary: %d Fatal, %d Events Aborted, "
01201                         "%d Tracks Aborted, %d Warnings\n",
01202                         fatalExceptions,eventsAborted,tracksAborted,warnings);
01203         std::map<G4String,int>::iterator it;
01204         for(it=exceptionCount.begin(); it!=exceptionCount.end(); ++it) {
01205                 G4String except = it->first;
01206                 int count = it->second;
01207                 printf(" %6d times: %s\n",count,except.c_str());
01208         }
01209 }


Member Data Documentation

BLManager * BLManager::blManager = 0 [static, private]

bool BLManager::initialized = false [static, private]

Referenced by initialize(), and isInitialized().

G4int BLManager::steppingVerbose [private]

unsigned int BLManager::beamIndex [private]

G4int BLManager::histoUpdate [private]

time_t BLManager::startRun [private]

time_t BLManager::startEvent [private]

G4int BLManager::eventTimeLimit [private]

G4VPhysicalVolume* BLManager::worldPhysicalVolume [private]

G4int BLManager::eventID [private]

G4int BLManager::trackID [private]

G4ParticleDefinition* BLManager::trackParticleDef [private]

G4int BLManager::fatalExceptions [private]

Referenced by BLManager(), exceptionSummary(), and Notify().

G4int BLManager::eventsAborted [private]

Referenced by BLManager(), exceptionSummary(), and Notify().

G4int BLManager::tracksAborted [private]

G4int BLManager::warnings [private]

Referenced by BLManager(), exceptionSummary(), and Notify().

G4int BLManager::prevEventID [private]

G4int BLManager::eventsProcessed [private]

bool BLManager::endRun [private]

G4String BLManager::eventCutFile [private]

std::set<int> BLManager::eventList [private]

Referenced by readEventCutFile(), and skipEvent().

G4SteppingManager* BLManager::fpSteppingManager [private]

G4TrackingManager* BLManager::fpTrackingManager [private]

std::vector<BLBeam*> BLManager::beamVector [private]

std::vector<BLBeam*> BLManager::referenceVector [private]

std::vector<RunAction*> BLManager::runActionVector [private]

std::vector<RunAction*> BLManager::beamRunActionVector [private]

std::vector<EventAction*> BLManager::eventActionVector [private]

std::vector<SteppingAction*> BLManager::allStepVector [private]

std::map<G4VPhysicalVolume*,SteppingAction*> BLManager::allStepMap [private]

std::map<G4VPhysicalVolume*,SteppingAction*> BLManager::tpStepMap [private]

std::map<G4VPhysicalVolume*,SteppingAction*> BLManager::rpStepMap [private]

std::vector<SteppingAction*> BLManager::tpStepVector [private]

std::vector<SteppingAction*> BLManager::rpStepVector [private]

std::map<G4VPhysicalVolume*,SteppingAction*> BLManager::beamStepMap [private]

std::vector<SteppingAction*> BLManager::beamStepVector [private]

std::vector<int> BLManager::verboseFormat [private]

G4double BLManager::zTolerance [private]

std::vector<ZStep> BLManager::tuneZStep [private]

std::vector<ZStep> BLManager::referenceZStep [private]

std::vector<ZStep> BLManager::beamZStep [private]

std::vector<ZStep>* BLManager::currentZStep [private]

unsigned BLManager::indexZStep [private]

G4double BLManager::prevZ [private]

G4int BLManager::nStuck [private]

std::map<G4int,G4int> BLManager::trackIDMap [private]

std::vector<BLUserCode*> BLManager::userCodeVector [private]

std::map<G4String,int> BLManager::exceptionCount [private]

Referenced by exceptionSummary(), and Notify().


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