BLManager Class Reference

#include <BLManager.hh>

List of all members.


Detailed Description

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

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

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

Public Member Functions

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

Static Public Member Functions

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

Private Member Functions

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

Private Attributes

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

Static Private Attributes

static BLManagerblManager = 0
static bool initialized = false

Classes

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


Constructor & Destructor Documentation

BLManager::BLManager (  )  [private]

private constructor -- immediate construction only.

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

Referenced by getObject().

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

BLManager::~BLManager (  ) 

Destructor.

References blManager.

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


Member Function Documentation

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

Referenced by registerZStep().

00380 {
00381         // the bookends in vector (see BLManager()) ensure this works
00382         std::vector<ZStep>::iterator i;
00383         for(i=vector.begin(); i<vector.end(); ++i) {
00384                 if(z < i->z) {  // preserve order for ==
00385                         ZStep zs(z,action);
00386                         vector.insert(i,zs);
00387                         break;
00388                 }
00389         }
00390 }

void BLManager::readEventCutFile (  )  [private]

References eventCutFile, and eventList.

Referenced by BeginOfRunAction().

01453 {
01454         eventList.clear();
01455         if(eventCutFile.size() == 0) return;
01456 
01457         FILE *in = fopen(eventCutFile.c_str(),"r");
01458         if(!in) {
01459                 printf("WARNING: Cannot read eventCutFile '%s' -- IGNORED\n",
01460                                                         eventCutFile.c_str());
01461                 return;
01462         }
01463 
01464         char line[1025];
01465         while(fgets(line,sizeof(line),in)) {
01466                 assert(strlen(line) < sizeof(line)-1);
01467                 if(line[0] == '#') continue;
01468                 long ev = strtol(line,0,0);
01469                 eventList.insert(ev);
01470         }
01471 
01472         fclose(in);
01473 
01474         // eventList is ignored if empty, so insert a crazy entry
01475         if(eventList.size() == 0) eventList.insert(0x80000000);
01476 }

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

00373 {
00374         if(!blManager) new BLManager();
00375         return blManager;
00376 }

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

00287 {
00288         runManager = new BLRunManager();
00289         G4StateManager * stateManager = G4StateManager::GetStateManager();
00290         stateManager->SetExceptionHandler(this);
00291 }

BLManagerState BLManager::getState (  )  [inline]

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

int BLManager::getSteppingVerbose (  )  [inline]

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

References steppingVerbose.

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

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

00194 { steppingVerbose = v; }

int BLManager::getEventTimeLimit (  )  [inline]

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

References eventTimeLimit.

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

00202 { eventTimeLimit = sec; }

G4int BLManager::getEventID (  )  const [inline]

getEventID() gets the current eventID;

References eventID.

Referenced by BLCMDbeam::nextBeamEvent(), and BLCMDusertrackfilter::UserSteppingAction().

00205 { return eventID; }

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

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

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

References eventID, and eventList.

Referenced by GeneratePrimaries().

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

void BLManager::incrEventsProcessed ( int  eventID  ) 

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

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

Referenced by EndOfEventAction().

01069 {
01070         // print event number, if appropriate
01071         ++eventsProcessed;
01072         if(!runManager->getCollectiveMode()) {
01073                 if(eventsProcessed <= 10 ||
01074                    (eventsProcessed < 100 && eventsProcessed%10 == 0) ||
01075                    (eventsProcessed < 1000 && eventsProcessed%100 == 0) ||
01076                    eventsProcessed%1000 == 0) {
01077                         printf("Event %d Completed",eventID);
01078                         int t = BLTime::time() - startRun;
01079                         if(t <= 0) t = 1;
01080                         printf("  %d events  realTime=%d sec  %.1f ev/sec",
01081                                 eventsProcessed,t,(double)eventsProcessed/t);
01082                         printf("\n");
01083                         fflush(stdout);
01084                 }
01085         }
01086 }

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.

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

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

00259                 { if(physicalVol != 0) rpStepMap[physicalVol] = sa; 
00260                   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 BLCMDtimentuple::command(), BLCMDmovie::command(), BLCMDvirtualdetector::construct(), and BLCMDtrackerplane::construct().

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

00393 {
00394         if(when & 1) insertZStep(tuneZStep,z,sa);
00395         if(when & 2) insertZStep(referenceZStep,z,sa);
00396         if(when & 4) insertZStep(beamZStep,z,sa);
00397 }

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

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

References stackingActionVector.

Referenced by BLCMDtrackcuts::command().

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

01159 {
01160         G4String fmt = Param.getString("steppingFormat");
01161         BLArgumentVector args;
01162         BLArgumentMap namedArgs;
01163         
01164         // convert "," to " "
01165         G4String s;
01166         for(unsigned i=0; i<fmt.size(); ++i)
01167                 s += (fmt[i]==',' ? ' ' : fmt[i]);
01168 
01169         BLCommand::parseArgs(s,args,namedArgs);
01170         appendVerboseFormat("");
01171         for(unsigned i=0; i<args.size(); ++i)
01172                 appendVerboseFormat(args[i]);
01173         appendVerboseFormat("NEWLINE");
01174 }

G4String BLManager::getFormatHelp (  ) 

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

01177 {
01178         G4String s = 
01179                 "        EXT       toggle extended precision (3 more digits)\n"
01180                 "        TAG       print a '>' (useful to grep output)\n"
01181                 "        N         step number\n"
01182                 "        NSTEP     Synonym of N\n"
01183                 "        GLOBAL    X,Y,Z,T in global coords\n"
01184                 "        XYZT      Synonym of GLOBAL\n"
01185                 "        CL        X,Y,Z,dxdz,dydz in CL coords\n"
01186                 "        CLX       extended precision CL\n"
01187                 "        KE        kinetic energy\n"
01188                 "        STEP      step length\n"
01189                 "        STEPLEN   Synonym of STEP\n"
01190                 "        VOL       volume name\n"
01191                 "        VOLNAME   Synonym of VOL\n"
01192                 "        PROCESS   process name\n"
01193                 "        B         magnetic field\n"
01194                 "        E         electric field\n"
01195                 "        P         3-momentum\n"
01196                 "        MAT       material name\n"
01197                 "        ID        event ID, track ID, parent ID\n"
01198                 "        PART      particle name\n"
01199                 "        SEG       centerline coord segment number\n"
01200                 "        WT        weight\n"
01201                 "        NL        <newline>\n"
01202                 "        NEWLINE   Synonym of NL\n"
01203                 "        \\n        Synonym of NL\n";
01204 
01205         return s;
01206 }

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

01209 {
01210         G4String f;
01211         for(const char *p=fmt.c_str(); *p; ++p)
01212                 f += toupper(*p);
01213         if(f == "") verboseFormat.clear();
01214         else if(f == "TAG") verboseFormat.push_back(TAG);
01215         else if(f == "EXT") verboseFormat.push_back(EXT);
01216         else if(f == "N") verboseFormat.push_back(NSTEP);
01217         else if(f == "NSTEP") verboseFormat.push_back(NSTEP);
01218         else if(f == "GLOBAL") verboseFormat.push_back(GLOBAL);
01219         else if(f == "XYZT") verboseFormat.push_back(GLOBAL);
01220         else if(f == "CL") verboseFormat.push_back(CL);
01221         else if(f == "CLX") verboseFormat.push_back(CLX);
01222         else if(f == "KE") verboseFormat.push_back(KE);
01223         else if(f == "STEP") verboseFormat.push_back(STEP);
01224         else if(f == "STEPLEN") verboseFormat.push_back(STEP);
01225         else if(f == "VOL") verboseFormat.push_back(VOL);
01226         else if(f == "VOLNAME") verboseFormat.push_back(VOL);
01227         else if(f == "PROCESS") verboseFormat.push_back(PROCESS);
01228         else if(f == "B") verboseFormat.push_back(B);
01229         else if(f == "E") verboseFormat.push_back(E);
01230         else if(f == "P") verboseFormat.push_back(P);
01231         else if(f == "MAT") verboseFormat.push_back(MAT);
01232         else if(f == "ID") verboseFormat.push_back(ID);
01233         else if(f == "PART") verboseFormat.push_back(PART);
01234         else if(f == "SEG") verboseFormat.push_back(SEG);
01235         else if(f == "WT") verboseFormat.push_back(WT);
01236         else if(f == "\n") verboseFormat.push_back(NEWLINE);
01237         else if(f == "NL") verboseFormat.push_back(NEWLINE);
01238         else if(f == "NEWLINE") verboseFormat.push_back(NEWLINE);
01239         else printf("BLManager: invalid verbose format '%s'\n",
01240                 fmt.c_str());
01241 }

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

01245 {
01246         G4StepPoint *prePoint=0;
01247         G4VPhysicalVolume *preVol=0;
01248         G4StepPoint *postPoint=0;
01249         G4double time=0;
01250         G4ThreeVector position;
01251         G4Material *mat=0;
01252         int eventNum = -999;
01253         int nstep=0;
01254         const char *procName=0;
01255         BLCoordinates *coord = (track!=0 ? (BLCoordinates *)track->GetUserInformation() : 0);
01256         if(coord && !coord->isValid()) coord = 0;
01257 
01258         if(!header) {
01259                 prePoint = step->GetPreStepPoint();
01260                 preVol = prePoint->GetPhysicalVolume();
01261                 postPoint = step->GetPostStepPoint();
01262                 time = track->GetGlobalTime();
01263                 position = track->GetPosition();
01264                 mat = preVol->GetLogicalVolume()->GetMaterial();
01265                 nstep = track->GetCurrentStepNumber();
01266                 procName = (postPoint->GetProcessDefinedStep() ?
01267                   postPoint->GetProcessDefinedStep()->GetProcessName().c_str() :
01268                   "UserLimit"); 
01269                 const G4Event* event = runManager->GetCurrentEvent();
01270                 eventNum = event->GetEventID();
01271         }
01272 
01273         bool extended=false;
01274         for(unsigned i=0; i<verboseFormat.size(); ++i) {
01275                 switch(verboseFormat[i]) {
01276                 case EXT:
01277                         extended = !extended;
01278                         break;
01279                 case TAG:
01280                         if(header) printf(" ");
01281                         else printf(">");
01282                         break;
01283                 case NSTEP:
01284                         if(header) printf(" Step");
01285                         else printf(" %4d",nstep);
01286                         break;
01287                 case GLOBAL:
01288                     if(!extended) {
01289                         if(header) printf("   X(mm)   Y(mm)   Z(mm)    T(ns)");
01290                         else printf(" %7.1f %7.1f %7.1f %8.2f",
01291                                 position[0],position[1],position[2],
01292                                 time);
01293                     } else {
01294                         if(header) printf("      X(mm)      Y(mm)      Z(mm)      T(ns) ");
01295                         else printf(" %10.4f %10.4f %10.4f %11.5f",
01296                                 position[0],position[1],position[2],
01297                                 time);
01298                     }
01299                         break;
01300                 case CL:
01301                     if(extended) goto clx;
01302                     if(header) {
01303                         printf("   CL: X        Y        Z     dxdz    dydz");
01304                     } else {
01305                         G4ThreeVector dir = track->GetMomentumDirection();
01306                         // transform to centerline coordinates
01307                         G4ThreeVector clpos;
01308                         coord->getCoords(BLCOORD_CENTERLINE,clpos);
01309                         dir = coord->getRotation() * dir;
01310                         printf(" %8.1f %8.1f %8.1f %7.4f %7.4f",
01311                                         clpos[0],clpos[1],clpos[2],
01312                                         dir[0]/dir[2],dir[1]/dir[2]);
01313                     }
01314                         break;
01315                 case CLX:
01316 clx:                if(header) {
01317                         printf("   CL: X           Y           Z        dxdz       dydz   ");
01318                     } else {
01319                         G4ThreeVector dir = track->GetMomentumDirection();
01320                         // transform to centerline coordinates
01321                         G4ThreeVector clpos;
01322                         coord->getCoords(BLCOORD_CENTERLINE,clpos);
01323                         dir = coord->getRotation() * dir;
01324                         printf(" %11.4f %11.4f %11.4f %10.7f %10.7f",
01325                                         clpos[0],clpos[1],clpos[2],
01326                                         dir[0]/dir[2],dir[1]/dir[2]);
01327                     }
01328                         break;
01329                 case KE:
01330                     if(!extended) {
01331                         if(header) printf("  KE(MeV)");
01332                         else printf(" %8.1f",track->GetKineticEnergy());
01333                     } else {
01334                         if(header) printf("    KE(MeV) ");
01335                         else printf(" %11.4f",track->GetKineticEnergy());
01336                     }
01337                         break;
01338                 case STEP:
01339                     if(!extended) {
01340                         if(header) printf("  StepLen");
01341                         else printf(" %8.2f",step->GetStepLength());
01342                     } else {
01343                         if(header) printf("    StepLen ");
01344                         else printf(" %11.5f",step->GetStepLength());
01345                     }
01346                         break;
01347                 case VOL:
01348                         if(header) printf(" This Volume     ");
01349                         else printf(" %-16s",preVol->GetName().c_str());
01350                         break;
01351                 case PROCESS:
01352                         if(header) printf(" Process        ");
01353                         else printf(" %-16s",procName);
01354                         break;
01355                 case B:
01356                     if(header) {
01357                         printf("       Bx,By,Bz (Tesla)     ");
01358                     } else {
01359                         G4double point[4], field[6];
01360                         point[0]=position[0],point[1]=position[1],point[2]=position[2];
01361                         point[3]=time;
01362                         BLGlobalField::getObject()->GetFieldValue(point,field);
01363                         G4ThreeVector B(field[0],field[1],field[2]);
01364                         B = coord->getRotation(BLCOORD_CENTERLINE) * B;
01365                         printf(" %8.4f %8.4f %8.4f",
01366                                 B[0]/tesla,B[1]/tesla,B[2]/tesla);
01367                     }
01368                         break;
01369                 case E:
01370                     if(header) {
01371                         printf("      Ex,Ey,Ez (MV/meter)  ");
01372                     } else {
01373                         G4double point[4], field[6];
01374                         point[0]=position[0],point[1]=position[1],point[2]=position[2];
01375                         point[3]=time;
01376                         BLGlobalField::getObject()->GetFieldValue(point,field);
01377                         G4ThreeVector E(field[3],field[4],field[5]);
01378                         E = coord->getRotation(BLCOORD_CENTERLINE) * E;
01379                         printf(" %8.4f %8.4f %8.4f",E[0]/(megavolt/meter),
01380                                 E[1]/(megavolt/meter),E[2]/(megavolt/meter));
01381                     }
01382                         break;
01383                 case P:
01384                   if(!extended) {
01385                     if(header) {
01386                         printf("      Px,Py,Pz (MeV/c)     ");
01387                     } else {
01388                         G4ThreeVector P = track->GetMomentum();
01389                         P = coord->getRotation(BLCOORD_CENTERLINE) * P;
01390                         printf(" %8.1f %8.1f %8.1f",P[0]/(MeV),
01391                                 P[1]/(MeV),P[2]/(MeV));
01392                     }
01393                   } else {
01394                     if(header) {
01395                         printf("          Px,Py,Pz (MeV/c)          ");
01396                     } else {
01397                         G4ThreeVector P = track->GetMomentum();
01398                         P = coord->getRotation(BLCOORD_CENTERLINE) * P;
01399                         printf(" %11.4f %11.4f %11.4f",P[0]/(MeV),
01400                                 P[1]/(MeV),P[2]/(MeV));
01401                     }
01402                   }
01403                         break;
01404                 case MAT:
01405                     if(header) {
01406                         printf(" Material ");
01407                     } else {
01408                         printf(" %-9s",mat->GetName().c_str());
01409                     }
01410                         break;
01411                 case ID:
01412                     if(header) {
01413                         printf(" Event# Trk Prnt");
01414                     } else {
01415                         printf(" %6d %3d %4d",eventNum,
01416                                 BLManager::getObject()->
01417                                     getExternalTrackID(track),
01418                                 BLManager::getObject()->
01419                                     getExternalParentID(track));
01420                     }
01421                         break;
01422                 case PART:
01423                     if(header) {
01424                         printf(" Particle");
01425                     } else {
01426                         printf(" %8s",track->GetDefinition()->GetParticleName().c_str());
01427                     }
01428                         break;
01429                 case SEG:
01430                     if(header) {
01431                         printf(" Seg");
01432                     } else {
01433                         printf(" %3d",coord->getSegmentCL());
01434                     }
01435                         break;
01436                 case WT:
01437                         if(header) {
01438                             printf("  Wt  ");
01439                         } else {
01440                             printf(" %.3f",track->GetWeight());
01441                         }
01442                         break;
01443                 case NEWLINE:
01444                         printf("\n");
01445                         break;
01446                 }
01447         }
01448 
01449         fflush(stdout); // help debugging location of crashes
01450 }

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

setSteppingManager() sets the pointer to the current G4SteppingManager.

References fpSteppingManager.

Referenced by BLManager_UserSteppingAction::UserSteppingAction().

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

00310 { return fpSteppingManager; }

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

setTrackingManager() sets the pointer to the current G4TrackingManager.

References fpTrackingManager.

Referenced by BLManager_UserTrackingAction::PreUserTrackingAction().

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

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

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

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

References initialized.

Referenced by BLCMDlist::command().

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

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

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

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

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

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

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

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

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

00345                 { physics = _physics;
00346                   runManager->SetUserInitialization(physics->getPhysicsList());
00347                 }

BLPhysics* BLManager::getPhysics (  )  [inline]

getPhysics returns the registered BLPhysics object.

References physics.

Referenced by BLCMDtrackermode::callback(), and BLCMDtracker::fitTrack().

00350 { return physics; }

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 BLCMDcosmicraybeam::command(), and BLCMDbeam::command().

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

00360                 { referenceVector.push_back(_beam); }

int BLManager::nReference (  )  [inline]

nReference() returns the number of reference particles registered.

References referenceVector.

Referenced by main().

00363 { return referenceVector.size(); }

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

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

References beamRunActionVector, and runActionVector.

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

00369                 { if(beamOnly)
00370                         beamRunActionVector.push_back(a);
00371                   else
00372                         runActionVector.push_back(a);
00373                 }

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

00379                 { if(beamOnly)
00380                         beamEventActionVector.push_back(a);
00381                   else
00382                         eventActionVector.push_back(a);
00383                 }

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

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

registerUserCode() registers a BLUserCode instance.

References userCodeVector.

Referenced by BLUserCode::registerUserCode().

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

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

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(), BLCMDshowmaterial::command(), BLCMDprofile::command(), BLCMDprobefield::command(), BLCMDprintfield::command(), BLCMDprintf::command(), BLCMDphysics::command(), BLCMDmovie::command(), BLCMDgeometry::command(), BLCMDg4ui::command(), BLCMDbug1021::command(), ParticleFilterPlacement::ParticleFilterPlacement(), and Surface::Surface().

00416                                                         {
00417                 if(type==0) preReferenceCallbackVector.push_back(cb);
00418                 else if(type==1) postReferenceCallbackVector.push_back(cb);
00419                 else if(type==2) postTrackingCallbackVector.push_back(cb);
00420                 else if(type==3) replaceMainLoopCallbackVector.push_back(cb);
00421                 else if(type==4) visualizationCallbackVector.push_back(cb);
00422         }

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(), eventsAborted, g4bl_exit(), BLRunManager::getObject(), postReferenceCallbackVector, postTrackingCallbackVector, preReferenceCallbackVector, replaceMainLoopCallbackVector, BLNTuple::summary(), visualizationCallbackVector, and warnings.

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

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

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 BLCMDshowmaterial::callback(), and BLCMDmovie::callback().

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

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

00445 { return seedMethod; }

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

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

References eventCutFile.

Referenced by BLCMDeventcuts::command().

00450 { eventCutFile = f; }

void BLManager::UserSteppingAction ( const G4Step *  step  ) 

UserSteppingAction() from G4UserSteppingAction.

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

Referenced by BLManager_UserSteppingAction::UserSteppingAction().

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

G4VPhysicalVolume * BLManager::Construct (  ) 

Construct() from G4VUserDetectorConstruction.

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

00513 {
00514         // ensure the global field is initialized
00515         (void)BLGlobalField::getObject();
00516 
00517         // construct the world
00518         worldPhysicalVolume =  BLGroup::constructWorld();
00519 
00520         return worldPhysicalVolume;
00521 }

void BLManager::PreUserTrackingAction ( const G4Track *  track  ) 

PreUserTrackingAction() from G4UserTrackingAction.

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

Referenced by BLManager_UserTrackingAction::PreUserTrackingAction().

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

void BLManager::PostUserTrackingAction ( const G4Track *  track  ) 

PostUserTrackingAction() from G4UserTrackingAction.

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

Referenced by BLManager_UserTrackingAction::PostUserTrackingAction().

00736 {
00737         std::vector<BLManager::TrackingAction*>::iterator i;
00738         for(i=trackingActionVector.begin(); i<trackingActionVector.end(); ++i) {
00739                 (*i)->PostUserTrackingAction(track);
00740         }
00741 
00742         // delete the BLcoordinates
00743         G4VUserTrackInformation *ti = track->GetUserInformation();
00744         if(ti && !runManager->getCollectiveMode()) {
00745                 delete ti;
00746                 ((G4Track *)track)->SetUserInformation(0);
00747         }
00748 
00749         trackID = -9999;
00750 }

void BLManager::BeginOfRunAction ( const G4Run *  run  ) 

BeginOfRunAction() from G4UserRunAction.

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

Referenced by BLManager_UserRunAction::BeginOfRunAction().

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

void BLManager::EndOfRunAction ( const G4Run *  run  ) 

EndOfRunAction() from G4UserRunAction.

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

Referenced by BLManager_UserRunAction::EndOfRunAction().

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

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

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

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

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

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

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

G4ClassificationOfNewTrack BLManager::ClassifyNewTrack ( const G4Track *  track  ) 

ClassifyNewTrack() from G4UserStackingAction.

References stackingActionVector.

Referenced by BLManager_UserStackingAction::ClassifyNewTrack().

01012 {
01013         std::vector<BLManager::StackingAction*>::iterator i;
01014         for(i=stackingActionVector.begin(); i<stackingActionVector.end(); ++i) {
01015                 G4ClassificationOfNewTrack c = (*i)->ClassifyNewTrack(track);
01016                 if(c == fKill) return c;
01017         }
01018         return fUrgent;
01019 }

void BLManager::NewStage (  ) 

NewStage() from G4userStackingAction.

References stackingActionVector.

Referenced by BLManager_UserStackingAction::NewStage().

01022 {
01023         std::vector<BLManager::StackingAction*>::iterator i;
01024         for(i=stackingActionVector.begin(); i<stackingActionVector.end(); ++i) {
01025                 (*i)->NewStage();
01026         }
01027 }

void BLManager::PrepareNewEvent (  ) 

PrepareNewEvent() from G4StackingAction.

References stackingActionVector.

Referenced by BLManager_UserStackingAction::PrepareNewEvent().

01030 {
01031         std::vector<BLManager::StackingAction*>::iterator i;
01032         for(i=stackingActionVector.begin(); i<stackingActionVector.end(); ++i) {
01033                 (*i)->PrepareNewEvent();
01034         }
01035 }

void BLManager::clearTrackIDMap (  )  [inline]

clearTrackIDMap() clears the TrackID map.

References trackIDMap.

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

00495 { trackIDMap.clear(); }

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(), and BLCMDbeam::nextBeamEvent().

00500 { nextSecondaryTrackID = next; }

int BLManager::getNextSecondaryTrackID (  )  [inline]

References nextSecondaryTrackID.

Referenced by BLRunManager::beamOnCollective().

00501 { return nextSecondaryTrackID; }

int BLManager::getExternalTrackID ( const G4Track *  track  ) 

int BLManager::getExternalParentID ( const G4Track *  track  ) 

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

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

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

01056 {
01057         if(trackID < 0) trackID = nextSecondaryTrackID++;
01058 
01059         BLTrackInfo *ti = (BLTrackInfo *)track->GetUserInformation();
01060         if(!ti || !ti->isValid()) {
01061                 ti = new BLTrackInfo();
01062                 track->SetUserInformation(ti);
01063         }
01064         ti->setExternalTrackID(trackID);
01065         ti->setExternalParentID(parentID);
01066 }

int BLManager::getPrimaryTrackID (  )  [inline]

getPrimaryTrackID() returns the primaryTrackID set by the beam command

References primaryTrackID.

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

00518 { return primaryTrackID; }

int BLManager::getPrimaryParentID (  )  [inline]

getPrimaryParentID() returns the primaryParentID set by the beam command

References primaryParentID.

Referenced by BLRunManager::getNextBeamEventAndTrack().

00521 { return primaryParentID; }

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

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

Notify() from G4VExceptionHandler.

References BLNTuple::closeAll(), eventID, eventsAborted, g4bl_exit(), handleCallbacks(), runManager, BLNTuple::summary(), trackID, and warnings.

01091 {
01092         fflush(stdout);
01093 
01094         G4bool toBeAborted = false;
01095         const char *p="UNKNOWN";
01096         switch(severity) {
01097         case FatalException:
01098                 p = "Fatal Exception";          
01099                 toBeAborted = true;
01100                 break;
01101         case FatalErrorInArgument:
01102                 p = "Fatal Error in Argument";  
01103                 toBeAborted = true;
01104                 break;
01105         case RunMustBeAborted:
01106                 p = "Run Must Be Aborted";
01107                 runManager->AbortRun(false);
01108                 toBeAborted = false;
01109                 break;
01110         case EventMustBeAborted:
01111                 p = "Event must be Aborted";
01112                 runManager->AbortEvent();
01113                 toBeAborted = false;
01114                 ++eventsAborted;
01115                 break;
01116         default:
01117                 p = "Warning";
01118                 toBeAborted = false;
01119                 ++warnings;
01120                 break;
01121         }
01122 
01123         fprintf(stderr,"**************************************************************************\n");
01124         fprintf(stderr,"*** G4Exception: %s\n",exceptionCode);
01125         fprintf(stderr,"***    severity: %s\n",p);
01126         fprintf(stderr,"***   issued by: %s\n",originOfException);
01127         if(strlen(description) > 0)
01128                 fprintf(stderr,"*** description: %s\n",description);
01129         if(strstr(description,"Missing mandatory data")) {
01130                 // brief description already printed, just URL is missing
01131                 const char *p = getenv("G4BL_DIR");
01132                 if(!p) p = "Installation_Directory";
01133                 fprintf(stderr,
01134                     "***              The URL is http://geant4.cern.ch\n"
01135                     "***              If you unpack it into %s/data\n"
01136                     "***              then the g4bl and g4blgui scripts will"
01137                     " automatically find it.\n",p);
01138         } else if(G4StateManager::GetStateManager()->GetCurrentState() ==
01139                                                             G4State_EventProc) {
01140                 fprintf(stderr,"***     EventID: %d     TrackID: %d\n",
01141                                                         eventID,trackID);
01142         }
01143         fprintf(stderr,"**************************************************************************\n");
01144 
01145         if(toBeAborted) {
01146                 fprintf(stderr,"g4beamline: attempting to close up after G4Exception...\n");
01147                 BLNTuple::summary();
01148                 BLNTuple::closeAll();
01149                 handleCallbacks(2);
01150                 fprintf(stderr,"g4beamline: simulation aborted\n");
01151                 extern void g4bl_exit(int);
01152                 g4bl_exit(99);
01153         }
01154 
01155         return toBeAborted;
01156 }


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]

G4int BLManager::eventsAborted [private]

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

G4int BLManager::warnings [private]

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

G4int BLManager::prevEventID [private]

G4int BLManager::eventsProcessed [private]

bool BLManager::endRun [private]

G4String BLManager::eventCutFile [private]

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

Referenced by readEventCutFile(), and skipEvent().

G4SteppingManager* BLManager::fpSteppingManager [private]

G4TrackingManager* BLManager::fpTrackingManager [private]

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

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

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

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

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

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

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

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

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

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

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

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

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

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

G4double BLManager::zTolerance [private]

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

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

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

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

unsigned BLManager::indexZStep [private]

G4double BLManager::prevZ [private]

G4int BLManager::nStuck [private]

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

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


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