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() determines if this EventID should be skipped. It checks keepEventList and skipEventList.
void setKeepEvent (int EventID)
 setKeepEvent() puts EventID into the keepEventList. (if keepEventList is empty, all events are kept.)
void setSkipEvent (int EventID)
 setSkipEvent() puts EventID into the skipEventList.
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, bool front=false)
 registerTrackingAction() registers a UserTackingAction. By default puts new entry at the back end of the vector; if front is true, puts new entry at the front.
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 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.
static void gsl_error_handler (const char *reason, const char *file, int lineno, int gsl_errno)
 GSL error handler.

Private Member Functions

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

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
std::set< int > keepEventList
std::set< int > skipEventList
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, eventID, eventsAborted, eventsProcessed, eventTimeLimit, fatalExceptions, BLParam::getDouble(), GLOBAL, gsl_error_handler(), histoUpdate, IDLE, indexZStep, BLAlarm::init(), BLSignal::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().

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

BLManager::~BLManager (  ) 

Destructor.

References blManager.

00308 {
00309         // cannot delete runManager, as it deletes the detector, 
00310         // steppingaction, eventaction, ..., and that becomes multiple
00311         // deletes of this object. But we are exiting, so we don't really
00312         // need to delete anything.
00313         //delete runManager;
00314 
00315         blManager = 0;
00316 }


Member Function Documentation

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

Referenced by registerZStep().

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

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(), BLTrackNTuple::appendTrack(), BLRunManager::BeamOn(), BLRunManager::beamOnCollective(), BLCMDmovie::beamStep(), BLRunManager::beginEvent(), BLRunManager::beginRun(), BLVirtualDetectorNTuple::BLVirtualDetectorNTuple(), BLCMDzntuple::callback(), BLCMDtrackermode::callback(), BLCMDtotalenergy::callback(), BLCMDshowmaterial::callback(), BLCMDgeometry::callback(), BLCMDtrackcuts::ClassifyNewTrack(), BLCMDzntuple::command(), BLCMDusertrackfilter::command(), BLCMDtune::command(), BLCMDtrackermode::command(), BLCMDtracker::command(), BLCMDtrackcuts::command(), BLCMDtrace::command(), BLCMDtotalenergy::command(), BLCMDtimentuple::command(), BLCMDspacechargelw::command(), BLCMDspacecharge::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(), BLCMDmuminuscapturefix::command(), BLCMDmovie::command(), BLCMDlist::command(), BLCMDhelp::command(), BLCMDgeometry::command(), BLCMDg4ui::command(), BLCMDfieldntuple::command(), BLCMDfieldlines::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(), BLCMDhelicalharmonic::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(), BLCMDeventcuts::readFile(), 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().

00383 {
00384         if(!blManager) new BLManager();
00385         return blManager;
00386 }

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

00301 {
00302         runManager = new BLRunManager();
00303         G4StateManager * stateManager = G4StateManager::GetStateManager();
00304         stateManager->SetExceptionHandler(this);
00305 }

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 BLCMDspacechargelw::beginTrack(), BLVirtualDetectorNTuple::BLVirtualDetectorNTuple(), BLCoordinates::init(), and BLCoordinates::start().

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

00197 { steppingVerbose = v; }

int BLManager::getEventTimeLimit (  )  [inline]

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

References eventTimeLimit.

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

00205 { eventTimeLimit = sec; }

G4int BLManager::getEventID (  )  const [inline]

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

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

skipEvent() determines if this EventID should be skipped. It checks keepEventList and skipEventList.

References keepEventList, and skipEventList.

Referenced by GeneratePrimaries().

00215                                     {
00216           if(skipEventList.count(EventID) > 0) return true;
00217           return keepEventList.size() > 0 && keepEventList.count(EventID) == 0;
00218         }

void BLManager::setKeepEvent ( int  EventID  )  [inline]

setKeepEvent() puts EventID into the keepEventList. (if keepEventList is empty, all events are kept.)

References keepEventList.

Referenced by BLCMDeventcuts::readFile().

00222 { keepEventList.insert(EventID); }

void BLManager::setSkipEvent ( int  EventID  )  [inline]

setSkipEvent() puts EventID into the skipEventList.

References skipEventList.

Referenced by BLCMDeventcuts::readFile().

00225 { skipEventList.insert(EventID); }

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

01095 {
01096         // print event number, if appropriate
01097         ++eventsProcessed;
01098         if(!runManager->getCollectiveMode()) {
01099                 if(eventsProcessed <= 10 ||
01100                    (eventsProcessed < 100 && eventsProcessed%10 == 0) ||
01101                    (eventsProcessed < 1000 && eventsProcessed%100 == 0) ||
01102                    eventsProcessed%1000 == 0) {
01103                         printf("Event %d Completed",eventID);
01104                         int t = BLTime::time() - startRun;
01105                         if(t <= 0) t = 1;
01106                         printf("  %d events  realTime=%d sec  %.1f ev/sec",
01107                                 eventsProcessed,t,(double)eventsProcessed/t);
01108                         printf("\n");
01109                         fflush(stdout);
01110                 }
01111         }
01112 }

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.

00246                 { if(physicalVol != 0) allStepMap[physicalVol] = sa;
00247                   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().

00258                 { if(physicalVol != 0) tpStepMap[physicalVol] = sa; 
00259                   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().

00271                 { if(physicalVol != 0) rpStepMap[physicalVol] = sa; 
00272                   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().

00284                 { if(physicalVol != 0) beamStepMap[physicalVol] = sa; 
00285                   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().

00403 {
00404         if(when & 1) insertZStep(tuneZStep,z,sa);
00405         if(when & 2) insertZStep(referenceZStep,z,sa);
00406         if(when & 4) insertZStep(beamZStep,z,sa);
00407 }

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

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

References stackingActionVector.

Referenced by BLCMDtrackcuts::command().

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

01220 {
01221         G4String fmt = Param.getString("steppingFormat");
01222         BLArgumentVector args;
01223         BLArgumentMap namedArgs;
01224         
01225         // convert "," to " "
01226         G4String s;
01227         for(unsigned i=0; i<fmt.size(); ++i)
01228                 s += (fmt[i]==',' ? ' ' : fmt[i]);
01229 
01230         BLCommand::parseArgs(s,args,namedArgs);
01231         appendVerboseFormat("");
01232         for(unsigned i=0; i<args.size(); ++i)
01233                 appendVerboseFormat(args[i]);
01234         appendVerboseFormat("NEWLINE");
01235 }

G4String BLManager::getFormatHelp (  ) 

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

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

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

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

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

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

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

setSteppingManager() sets the pointer to the current G4SteppingManager.

References fpSteppingManager.

Referenced by BLManager_UserSteppingAction::UserSteppingAction().

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

00322 { return fpSteppingManager; }

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

setTrackingManager() sets the pointer to the current G4TrackingManager.

References fpTrackingManager.

Referenced by BLManager_UserTrackingAction::PreUserTrackingAction().

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

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

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

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

References initialized.

Referenced by BLCMDlist::command().

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

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

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

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

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

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

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

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

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

00357                 { physics = _physics;
00358                   runManager->SetUserInitialization(physics->getPhysicsList());
00359                 }

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

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

00372                 { referenceVector.push_back(_beam); }

int BLManager::nReference (  )  [inline]

nReference() returns the number of reference particles registered.

References referenceVector.

Referenced by main().

00375 { 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(), BLCMDhelicalharmonic::construct(), and BLCMDhelicaldipole::construct().

00381                 { if(beamOnly)
00382                         beamRunActionVector.push_back(a);
00383                   else
00384                         runActionVector.push_back(a);
00385                 }

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

00391                 { if(beamOnly)
00392                         beamEventActionVector.push_back(a);
00393                   else
00394                         eventActionVector.push_back(a);
00395                 }

void BLManager::registerTrackingAction ( TrackingAction a,
bool  front = false 
) [inline]

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

registerUserCode() registers a BLUserCode instance.

References userCodeVector.

Referenced by BLUserCode::registerUserCode().

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

00414                 { std::vector<BLUserCode*> ret;
00415                   for(unsigned i=0; i<userCodeVector.size(); ++i) {
00416                         if(type == userCodeVector[i]->getType())
00417                                 ret.push_back(userCodeVector[i]);
00418                   }
00419                   return ret;
00420                 }

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(), BLCMDmuminuscapturefix::command(), BLCMDmovie::command(), BLCMDgeometry::command(), BLCMDg4ui::command(), BLCMDfieldntuple::command(), BLCMDfieldlines::command(), BLCMDbug1021::command(), SetDecayInstance::fatalError(), ParticleFilterPlacement::ParticleFilterPlacement(), and Surface::Surface().

00434                                                         {
00435                 if(type==0) preReferenceCallbackVector.push_back(cb);
00436                 else if(type==1) postReferenceCallbackVector.push_back(cb);
00437                 else if(type==2) postTrackingCallbackVector.push_back(cb);
00438                 else if(type==3) replaceMainLoopCallbackVector.push_back(cb);
00439                 else if(type==4) visualizationCallbackVector.push_back(cb);
00440         }

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 BLCMDlist::command(), displayVisual(), main(), and Notify().

00537 {
00538         // startup sequencing prevents RunManager from registering
00539         // so just do it
00540         BLRunManager::getObject()->callback(type);
00541 
00542         if(type == 0) {
00543                 for(unsigned i=0; i<preReferenceCallbackVector.size(); ++i)
00544                         preReferenceCallbackVector[i]->callback(type);
00545         } else if(type == 1) {
00546                 for(unsigned i=0; i<postReferenceCallbackVector.size(); ++i)
00547                         postReferenceCallbackVector[i]->callback(type);
00548         } else if(type == 2) {
00549                 for(unsigned i=0; i<postTrackingCallbackVector.size(); ++i)
00550                         postTrackingCallbackVector[i]->callback(type);
00551                 exceptionSummary();
00552         } else if(type == 3) {
00553                 for(unsigned i=0; i<replaceMainLoopCallbackVector.size(); ++i)
00554                         replaceMainLoopCallbackVector[i]->callback(type);
00555                 if(replaceMainLoopCallbackVector.size() > 0) {
00556                         // we have replaced the main loop, so closeup and exit
00557                         BLNTuple::summary();
00558                         BLNTuple::closeAll();
00559                         // handle post-tracking callbacks
00560                         handleCallbacks(2);
00561                         blManager = 0;
00562                         g4bl_exit(0);
00563                 }
00564         } else if(type == 4) {
00565                 for(unsigned i=0; i<visualizationCallbackVector.size(); ++i)
00566                         visualizationCallbackVector[i]->callback(type);
00567         }
00568 }

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

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

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

00463 { return seedMethod; }

void BLManager::UserSteppingAction ( const G4Step *  step  ) 

UserSteppingAction() from G4UserSteppingAction.

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

Referenced by BLManager_UserSteppingAction::UserSteppingAction().

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

G4VPhysicalVolume * BLManager::Construct (  ) 

Construct() from G4VUserDetectorConstruction.

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

00523 {
00524         // ensure the global field is initialized
00525         (void)BLGlobalField::getObject();
00526 
00527         // construct the world
00528         worldPhysicalVolume =  BLGroup::constructWorld();
00529 
00530         // make it invisible
00531         worldPhysicalVolume->GetLogicalVolume()->SetVisAttributes(G4VisAttributes::GetInvisible());
00532 
00533         return worldPhysicalVolume;
00534 }

void BLManager::PreUserTrackingAction ( const G4Track *  track  ) 

PreUserTrackingAction() from G4UserTrackingAction.

References beamZStep, BLAssert, 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().

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

void BLManager::PostUserTrackingAction ( const G4Track *  track  ) 

PostUserTrackingAction() from G4UserTrackingAction.

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

Referenced by BLManager_UserTrackingAction::PostUserTrackingAction().

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

void BLManager::BeginOfRunAction ( const G4Run *  run  ) 

BeginOfRunAction() from G4UserRunAction.

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

Referenced by BLManager_UserRunAction::BeginOfRunAction().

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

void BLManager::EndOfRunAction ( const G4Run *  run  ) 

EndOfRunAction() from G4UserRunAction.

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

Referenced by BLManager_UserRunAction::EndOfRunAction().

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

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

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

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

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

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

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

G4ClassificationOfNewTrack BLManager::ClassifyNewTrack ( const G4Track *  track  ) 

ClassifyNewTrack() from G4UserStackingAction.

References stackingActionVector.

Referenced by BLManager_UserStackingAction::ClassifyNewTrack().

01038 {
01039         std::vector<BLManager::StackingAction*>::iterator i;
01040         for(i=stackingActionVector.begin(); i<stackingActionVector.end(); ++i) {
01041                 G4ClassificationOfNewTrack c = (*i)->ClassifyNewTrack(track);
01042                 if(c == fKill) return c;
01043         }
01044         return fUrgent;
01045 }

void BLManager::NewStage (  ) 

NewStage() from G4userStackingAction.

References stackingActionVector.

Referenced by BLManager_UserStackingAction::NewStage().

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

void BLManager::PrepareNewEvent (  ) 

PrepareNewEvent() from G4StackingAction.

References stackingActionVector.

Referenced by BLManager_UserStackingAction::PrepareNewEvent().

01056 {
01057         std::vector<BLManager::StackingAction*>::iterator i;
01058         for(i=stackingActionVector.begin(); i<stackingActionVector.end(); ++i) {
01059                 (*i)->PrepareNewEvent();
01060         }
01061 }

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

00513 { nextSecondaryTrackID = next; }

int BLManager::getNextSecondaryTrackID (  )  [inline]

References nextSecondaryTrackID.

Referenced by BLRunManager::beamOnCollective().

00514 { return nextSecondaryTrackID; }

int BLManager::getExternalTrackID ( const G4Track *  track  ) 

getExternalTrackID() returns the external TrackID for the given track. In collective mode, internal and external trackID-s are the same.

References BLRunManager::getCollectiveMode(), BLTrackInfo::getExternalTrackID(), BLCoordinates::isValid(), runManager, and trackIDMap.

Referenced by BLTrackNTuple::appendTrack(), BLCMDmovie::beamStep(), BLCMDtrace::newTrace(), steppingVerbosePrint(), and BLCMDtracker::UserZSteppingAction().

01064 {
01065         if(runManager->getCollectiveMode()) return track->GetTrackID();
01066         BLTrackInfo *p = (BLTrackInfo *)track->GetUserInformation();
01067         if(p && p->isValid()) 
01068                 return p->getExternalTrackID();
01069         return trackIDMap[track->GetTrackID()];
01070 }

int BLManager::getExternalParentID ( const G4Track *  track  ) 

getExternalParentID() returns the external ParentID for the given track. In collective mode, internal and external trackID-s are the same.

References BLRunManager::getCollectiveMode(), BLTrackInfo::getExternalParentID(), BLCoordinates::isValid(), runManager, and trackIDMap.

Referenced by BLTrackNTuple::appendTrack(), BLCMDmovie::beamStep(), PreUserTrackingAction(), steppingVerbosePrint(), and BLCMDtracker::UserZSteppingAction().

01073 {
01074         if(runManager->getCollectiveMode()) return track->GetParentID();
01075         BLTrackInfo *p = (BLTrackInfo *)track->GetUserInformation();
01076         if(p && p->isValid()) 
01077                 return p->getExternalParentID();
01078         return trackIDMap[track->GetParentID()];
01079 }

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

01082 {
01083         if(trackID < 0) trackID = nextSecondaryTrackID++;
01084 
01085         BLTrackInfo *ti = (BLTrackInfo *)track->GetUserInformation();
01086         if(!ti || !ti->isValid()) {
01087                 ti = new BLTrackInfo();
01088                 track->SetUserInformation(ti);
01089         }
01090         ti->setExternalTrackID(trackID);
01091         ti->setExternalParentID(parentID);
01092 }

int BLManager::getPrimaryTrackID (  )  [inline]

getPrimaryTrackID() returns the primaryTrackID set by the beam command

References primaryTrackID.

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

00531 { return primaryTrackID; }

int BLManager::getPrimaryParentID (  )  [inline]

getPrimaryParentID() returns the primaryParentID set by the beam command

References primaryParentID.

Referenced by BLRunManager::getNextBeamEventAndTrack().

00534 { 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 BLAlarm::clear(), BLNTuple::closeAll(), eventID, eventsAborted, exceptionCount, fatalExceptions, g4bl_exit(), handleCallbacks(), BLMPI::isMPI(), runManager, BLNTuple::summary(), trackID, trackParticleDef, and warnings.

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

void BLManager::exceptionSummary (  ) 

exceptionSummary() prints a summary of all exceptions

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

Referenced by handleCallbacks().

01207 {
01208         printf("\nException Summary: %d Fatal, %d Events Aborted, "
01209                         "%d Tracks Aborted, %d Warnings\n",
01210                         fatalExceptions,eventsAborted,tracksAborted,warnings);
01211         std::map<G4String,int>::iterator it;
01212         for(it=exceptionCount.begin(); it!=exceptionCount.end(); ++it) {
01213                 G4String except = it->first;
01214                 int count = it->second;
01215                 printf(" %6d times: %s\n",count,except.c_str());
01216         }
01217 }

void BLManager::gsl_error_handler ( const char *  reason,
const char *  file,
int  lineno,
int  gsl_errno 
) [static]

GSL error handler.

References snprintf.

Referenced by BLManager().

01515 {
01516         // ignore underflows
01517         if(gsl_errno == GSL_EUNDRFLW) return;
01518 
01519         char tmp[256];
01520         snprintf(tmp,sizeof(tmp)-1,"file=%s line=%d gsl_errno=%d",file,lineno,
01521                                                                 gsl_errno);
01522         G4Exception("BLManager",reason,EventMustBeAborted,tmp);
01523 }


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]

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

Referenced by setKeepEvent(), and skipEvent().

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

Referenced by setSkipEvent(), 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