BLCMDphysics Class Reference

Inheritance diagram for BLCMDphysics:

BLPhysics BLCommand BLCallback

List of all members.


Detailed Description

class BLCMDphysics implements a command to select from the geant4 physics lists of physics processes.

Public Member Functions

 BLCMDphysics ()
 Default Constructor.
 ~BLCMDphysics ()
 Destructor.
void setDoStochastics (BLSetValue value, G4bool quiet=false)
 setDoStochastics() sets whether or not stochastic processes are to be enabled.
G4VUserPhysicsList * getPhysicsList ()
 getPhysicsList() returns the G4PhysicsList.
G4String commandName ()
 commandName() returns "physics".
int command (BLArgumentVector &argv, BLArgumentMap &namedArgs)
 command() implements the command.
void defineNamedArgs ()
 defineNamedArgs() defines the named arguments.
void help (bool detailed)
 help() provides special help text for the physics command.
void callback (int type)
 callback() is from BLCallback, and adds a StepLimiter to the physics list for every particle.
bool isLocalPhysicsList (G4String name)
 isLocalPhysicsList() returns true if thsi is a local physics list.
G4VUserPhysicsList * getLocalPhysicsList (G4String name)
 getLocalPhysicsList() returns a local physics list.
void listPhysicsLists ()
 listPhysicsLists() will list all known physics lists, with synopsis.
void printCase (G4String name, G4String description)
 printCase() prints the description of a physics list.
G4String getSynopsis (G4String listName)
 getSynopsis() returns the synopsis for a physics list, by name.

Private Attributes

G4VUserPhysicsList * physicsList
G4int doStochastics
G4String disable
G4int list
std::vector< G4String > disableList
G4PhysListFactory * physListFactory

Constructor & Destructor Documentation

BLCMDphysics::BLCMDphysics (  ) 

Default Constructor.

References BLCMDTYPE_AUX, disable, doStochastics, list, BLPhysics::minRangeCut, physicsList, physListFactory, BLCommand::registerCommand(), BLCommand::setDescription(), and BLCommand::setSynopsis().

00095                            : BLPhysics(), BLCommand(), disableList()
00096 {
00097         registerCommand(BLCMDTYPE_AUX);
00098         setSynopsis("Defines the physics processes and controls them.");
00099         setDescription("Exactly one physics command must be present.\n"
00100                 "This command implements the geant4 physics lists\n"
00101                 "of physics processes.\n"
00102                 "The command is 'physics QGSP' for the QGSP set, and\n"
00103                 "similarly for the other sets. With no argument it prints\n"
00104                 "the available physics lists.\n"
00105                 "Note that stochastic processes are always disabled while\n"
00106                 "tracking the tune and reference particles.\n"
00107                 "The only non-stochastic processes are Transportation and\n"
00108                 "ionization energy loss (with fluctuations disabled).\n"
00109                 "For muon beam studies you may want to disable the Decay "
00110                 "process.\n\n"
00111                 "NOTE: this command defines the particles used throughout "
00112                 "the simulation, so this command must come before others "
00113                 "that use particle names.\n\n"
00114                 "The default all-around physics list for HEP is called "
00115                 "'QGSP_BERT'.");
00116 
00117         physicsList = 0;
00118         doStochastics = 1;
00119         list = 0;
00120         disable = "";
00121         minRangeCut = 1.0 * mm;
00122         physListFactory = 0;
00123 }

BLCMDphysics::~BLCMDphysics (  )  [inline]

Destructor.

00052 { }


Member Function Documentation

void BLCMDphysics::setDoStochastics ( BLSetValue  value,
G4bool  quiet = false 
) [virtual]

setDoStochastics() sets whether or not stochastic processes are to be enabled.

Implements BLPhysics.

References disableList, doStochastics, FORCE_OFF, FORCE_ON, BLPhysics::isStochasticProcess(), list, BLPhysics::minRangeCut, NORMAL, and BLPhysics::stochasticsEnabled.

00210 {
00211         if(value == NORMAL) {
00212                 if(doStochastics)
00213                         value = FORCE_ON;
00214                 else
00215                         value = FORCE_OFF;
00216         }
00217 
00218         G4UImanager* UI = G4UImanager::GetUIpointer();
00219         G4String cmd;
00220         char rangeCmd[64];
00221         switch(value) {
00222         case FORCE_OFF:
00223                 // turn off delta rays (and other particle production) by
00224                 // setting production cut very high in energy by setting
00225                 // the required range large.
00226                 UI->ApplyCommand("/run/setCut 1 km");
00227                 // Turn off enegry-loss straggling
00228                 UI->ApplyCommand("/process/eLoss/fluct false");
00229                 cmd = "/process/inactivate ";
00230                 if(!quiet) printf("*** NOTE: Stochastic processes are disabled!\n");
00231                 stochasticsEnabled = false;
00232                 break;
00233         case FORCE_ON:
00234                 // turn on delta rays (and other particle production)
00235                 sprintf(rangeCmd,"/run/setCut %.6g mm",minRangeCut/mm);
00236                 UI->ApplyCommand(rangeCmd);
00237                 // Turn on enegry-loss straggling
00238                 UI->ApplyCommand("/process/eLoss/fluct true");
00239                 cmd = "/process/activate ";
00240                 if(!quiet) printf("Stochastic processes are enabled.\n");
00241                 stochasticsEnabled = true;
00242                 break;
00243         case NORMAL:
00244                 return;
00245         }
00246 
00247         // run through list of processes, applying cmd to stochastic ones
00248         G4ProcessTable *pt = G4ProcessTable::GetProcessTable();
00249         G4ProcessTable::G4ProcNameVector *nv = pt->GetNameList();
00250         if(list) printf("List of physics processes:\n");
00251         for(unsigned int i=0; i<nv->size(); ++i) {
00252                 G4String name = (*nv)[i];
00253                 if(list) printf("\t%s\n",name.c_str());
00254                 if(isStochasticProcess(name))
00255                         UI->ApplyCommand(cmd+name);
00256         }
00257         if(list) printf("\n");
00258         list = 0;       // print list only once
00259 
00260         // now disable all proceses in the disableList
00261         cmd = "/process/inactivate ";
00262         for(unsigned i=0; i<disableList.size(); ++i) {
00263                 if(!quiet) printf("*** NOTE: Disabling physics process '%s'\n",
00264                                                 disableList[i].c_str());
00265                 UI->ApplyCommand(cmd+disableList[i]);
00266         }
00267 }

G4VUserPhysicsList* BLCMDphysics::getPhysicsList (  )  [inline, virtual]

getPhysicsList() returns the G4PhysicsList.

Implements BLPhysics.

References physicsList.

00059 { return physicsList; }

G4String BLCMDphysics::commandName (  )  [inline, virtual]

commandName() returns "physics".

Implements BLCommand.

00062 { return "physics"; }

int BLCMDphysics::command ( BLArgumentVector argv,
BLArgumentMap namedArgs 
) [virtual]

command() implements the command.

Implements BLCommand.

References disable, disableList, getLocalPhysicsList(), BLManager::getObject(), BLCommand::handleNamedArgs(), isLocalPhysicsList(), listPhysicsLists(), physicsList, physListFactory, BLCommand::print(), BLCommand::printError(), BLManager::registerCallback(), BLManager::registerPhysics(), and BLCommand::splitString().

00126 {
00127         if(physicsList != 0 && argv.size() > 0) {
00128                 printError("physics: Multiple physics commands not allowed!");
00129                 return -1;
00130         }
00131 
00132         if(argv.size() == 0)
00133                 return -1;
00134 
00135         // handle obsolete physics lists
00136         if(argv[0] == "QGSC") {
00137                 argv[0] = "QGSC_BERT";
00138                 printf("Physics list QGSC is obsolete, replaced by QGSC_BERT\n");
00139         }
00140         if(argv[0] == "QGSP") {
00141                 argv[0] = "QGSP_BERT";
00142                 printf("Physics list QGSP is obsolete, replaced by QGSP_BERT\n");
00143         }
00144 
00145         if(!physListFactory) physListFactory = new G4PhysListFactory();
00146         if(isLocalPhysicsList(argv[0]))
00147                 physicsList = getLocalPhysicsList(argv[0]);
00148         if(!physicsList) {
00149                 if(!physListFactory->IsReferencePhysList(argv[0])) {
00150                         printError("physics: Unknown physics list '%s'",
00151                                                         argv[0].c_str());
00152                         printError("The choices are:");
00153                         listPhysicsLists();
00154                         return -1;
00155                 }
00156                 physicsList = physListFactory->GetReferencePhysList(argv[0]);
00157         }
00158 
00159         handleNamedArgs(namedArgs);
00160 
00161         disableList = splitString(disable,',',true);
00162 
00163         BLManager::getObject()->registerPhysics(this);
00164 
00165         // now is too early to manage the physics processes; do it in a callback
00166         BLManager::getObject()->registerCallback(this,0);
00167 
00168         print(argv[0]);
00169 
00170         return 0;
00171 }

void BLCMDphysics::defineNamedArgs (  )  [virtual]

defineNamedArgs() defines the named arguments.

Reimplemented from BLCommand.

References BLCommand::argDouble(), BLCommand::argInt(), BLCommand::argString(), disable, doStochastics, list, and BLPhysics::minRangeCut.

00174 {
00175         argString(disable,"disable","Comma-separated list of processes to disable (e.g. 'Decay,msc').");
00176         argString(disable,"inactivate","Synonym of disable.");
00177         argString(disable,"deactivate","Synonym of disable.");
00178         argInt(doStochastics,"doStochastics","Set to zero to disable all stochastic processes.");
00179         argDouble(minRangeCut,"minRangeCut","Minimum range cut for particle production (1 mm)");
00180         argInt(list,"list","Nonzero to list the processes (later in output).");
00181 }

void BLCMDphysics::help ( bool  detailed  )  [virtual]

help() provides special help text for the physics command.

Reimplemented from BLCommand.

References BLCommand::help(), and listPhysicsLists().

00184 {
00185         BLCommand::help(detailed);
00186         if(detailed) {
00187                 printf("\n             ----- PHYSICS LISTS -----\n\n");
00188                 printf("             Further guidance in selecting physics lists is available at:\n"
00189                 "http://geant4.web.cern.ch/geant4/support/physicsLists/referencePL/index.shtml\n\n");
00190                 printf("             The default all-around physics list for HEP is called 'QGSP_BERT'.\n\n");
00191                 printf("             LHEP uses exclusively parameterized modeling.\n");
00192                 printf("             FTF lists use the FRITIOF description of string excitation\n"
00193                        "                  and fragmentation.\n");
00194                 printf("             QGSP lists use the quark gluon string model\n");
00195                 printf("             QGSC are as QGSP except applying CHIPS modeling for the nuclear\n"
00196                        "                  de-excitation.\n");
00197                 printf("             _BERT uses Geant4 Bertini cascade below ~ 10 GeV.\n");
00198                 printf("             _BIC uses Geant4 Binary cascade below ~ 10 GeV.\n");
00199                 printf("             _EMV suffix indicates a faster but less accurate EM modeling.\n");
00200                 printf("             _EMX suffix indicates the low-energy EM processes.\n");
00201                 printf("             _HP suffix uses the data driven high precision neutron package\n"
00202                        "                  (thermal to 20 MeV).\n");
00203                 printf("             _NQE suffix indicates a list for comparison with earlier release.\n");
00204                 printf("\n             List of available physics lists:\n");
00205                 listPhysicsLists();
00206         }
00207 }

void BLCMDphysics::callback ( int  type  )  [virtual]

callback() is from BLCallback, and adds a StepLimiter to the physics list for every particle.

Reimplemented from BLCallback.

00270 {
00271         G4ParticleTable::G4PTblDicIterator *theParticleIterator =
00272                         G4ParticleTable::GetParticleTable()->GetIterator();
00273         theParticleIterator -> reset();
00274         while( (*theParticleIterator)() ) {
00275                 G4ParticleDefinition* particle = theParticleIterator->value();
00276                 G4ProcessManager* pmanager = particle -> GetProcessManager();
00277                 if(!pmanager) {
00278                         printf("BLCMDphysics: particle '%s' has no ProcessManager!\n",
00279                                 particle->GetParticleName().c_str());
00280                         continue;
00281                 }
00282                 pmanager -> AddProcess(new G4StepLimiter(),-1,-1,3);
00283         }
00284 }

bool BLCMDphysics::isLocalPhysicsList ( G4String  name  ) 

isLocalPhysicsList() returns true if thsi is a local physics list.

Referenced by command().

00287 {
00288         if(name == "MICEPhysicsList") return true;
00289         if(name == "LISAPhysicsList") return true;
00290         return false;
00291 }

G4VUserPhysicsList * BLCMDphysics::getLocalPhysicsList ( G4String  name  ) 

getLocalPhysicsList() returns a local physics list.

Referenced by command().

00294 {
00295         if(name == "MICEPhysicsList") return new MICEPhysicsList();
00296         if(name == "LISAPhysicsList") return new LISAPhysicsList();
00297         return 0;
00298 }

void BLCMDphysics::listPhysicsLists (  ) 

listPhysicsLists() will list all known physics lists, with synopsis.

References getSynopsis(), physListFactory, and printCase().

Referenced by command(), and help().

00302 {
00303         if(!physListFactory) physListFactory = new G4PhysListFactory();
00304         const std::vector<G4String> v = physListFactory->AvailablePhysLists();
00305         for(unsigned i=0; i<v.size(); ++i) {
00306                 G4String n = v[i];
00307                 G4String s = getSynopsis(n);
00308                 printCase(n,s);
00309         }
00310         printCase("MICEPhysicsList",getSynopsis("MICEPhysicsList"));
00311         printCase("LISAPhysicsList",getSynopsis("LISAPhysicsList"));
00312 }

void BLCMDphysics::printCase ( G4String  name,
G4String  description 
)

printCase() prints the description of a physics list.

References BLCommand::wrapWords().

Referenced by listPhysicsLists().

00315 {
00316         G4String indent1("             "), indent("       ");
00317         indent += indent1;
00318         indent1 += name;
00319         do { indent1 += " "; } while(indent1.size() < indent.size());
00320         printf("%s",wrapWords(description,indent1,indent).c_str());
00321 }

G4String BLCMDphysics::getSynopsis ( G4String  listName  ) 

getSynopsis() returns the synopsis for a physics list, by name.

References CASE.

Referenced by listPhysicsLists().

00324 {
00325 #define CASE(NAME,SYNOPSIS)                             \
00326         if(listName == #NAME) return SYNOPSIS;
00327 
00328         CASE(CHIPS,"No synopsis available.");
00329         CASE(FTFP_BERT,"For calorimetry. The FTF model is based on the FRITIOF description of string excitation and fragmentation. Uses Geant4 Bertini cascade for primary protons, neutrons, pions and Kaons below ~10GeV.");
00330         CASE(FTFP_BERT_EMV,"Like FTFP_BERT but with faster EM modeling.");
00331         CASE(FTFP_BERT_EMX,"Like FTFP_BERT but with low-energy EM modeling.");
00332         CASE(FTFP_BERT_TRV,"A variant of QGSP_BERT where the Geant4 Bertini cascade is only used for particles below ~5.5 GeV.");
00333         CASE(FTFP_BIC,"For calorimetry. The FTF model is based on the FRITIOF description of string excitation and fragmentation. Uses Geant4 binary cascade for primary protons, neutrons, pions and Kaons below ~10GeV.");
00334         CASE(LBE,"For low background experiments (e.g. underground)");
00335         CASE(LHEP,"For calorimetry -- is the fastest, when it comes to CPU. It uses the LEP and HEP parametrized models for inelastic scattering. The modeling parametrizes the final states individual inelastic reactions, so you will not see resonances, and the detailed secondary angular distributions for O(100MeV) reactions may not be described perfectly. The average quantities will be well described.");
00336         CASE(LHEP_BERT,"For 'typical' HEP collider detectors and low-energy penetration shielding (< 5 GeV) -- is as LHEP for energies above 3GeV. Below 3 GeV the Bertini cascade is used for nucleon and pion induced reactions. The price to pay is reduced CPU performance. It does not include gamma and electro nuclear physics.");
00337         CASE(LHEP_EMV,"Like LHEP but with faster EM modeling.");
00338         CASE(QBBC,"No synopsis available.");
00339         CASE(QGSC,"For calorimetry and high energy physics trackers -- is as QGSP for the initial reaction, but uses chiral invariant phase-space decay (multi-quasmon fragmentation) to model the behavior of the system's fragmentation.");
00340         CASE(QGSC_BERT,"For calorimetry and high energy physics trackers -- is as QGSP for the initial reaction, but uses chiral invariant phase-space decay (multi-quasmon fragmentation) to model the behavior of the system's fragmentation. Uses Geant4 Bertini cascade for nucleon and pion induced reactions.");
00341         CASE(QGSP,"For calorimetry and high energy physics trackers and high-energy and medium-energy production targets -- uses theory driven modeling for the reactions of energetic pions, kaons, and nucleons. It employs quark gluon string model for the 'punch-through' interactions of the projectile with a nucleus, the string excitation cross-sections being calculated in quasi-eikonal approximation. A pre-equilibrium decay model with an extensive evaporation phase to model the behavior of the nucleus 'after the punch'. It uses current best pion cross-section.");
00342         CASE(QGSP_BERT,"Like QGSP, but using Geant4 Bertini cascade for primary protons, neutrons, pions and Kaons below ~10GeV. In comparison to experimental data we find improved agreement to data compared to QGSP which uses the low energy parameterised (LEP) model for all particles at these energies. The Bertini model produces more secondary neutrons and protons than the LEP model, yielding a better agreement to experimental data.");
00343         CASE(QGSP_BERT_EMV,"Like QGSP_BERT but with faster EM modeling.");
00344         CASE(QGSP_BERT_EMX,"Like QGSP_BERT but with low-energy EM modeling.");
00345         CASE(QGSP_BERT_HP,"Like QGSP_BERT but with _HP modeling for neutrons.");
00346         CASE(QGSP_BIC,"Like QGSP, but using Geant4 Binary cascade for primary protons and neutrons with energies below ~10GeV, thus replacing the use of the LEP model for protons and neutrons In comparison to teh LEP model, Binary cascade better describes production of secondary particles produced in interactions of protons and neutrons with nuclei.");
00347         CASE(QGSP_BIC_EMX,"Like QGSP_BIC but with low-energy EM modeling.");
00348         CASE(QGSP_BIC_HP,"Like QGSP_BIC but with _HP modeling for neutrons.");
00349 
00350         CASE(LISAPhysicsList,"Geant4 examples/advanced/cosmicray_charging");
00351         CASE(MICEPhysicsList,"The current MICE physics list.");
00352 
00353         return "No synopsis available.";
00354 }


Member Data Documentation

G4VUserPhysicsList* BLCMDphysics::physicsList [private]

G4int BLCMDphysics::doStochastics [private]

G4String BLCMDphysics::disable [private]

G4int BLCMDphysics::list [private]

std::vector<G4String> BLCMDphysics::disableList [private]

Referenced by command(), and setDoStochastics().

G4PhysListFactory* BLCMDphysics::physListFactory [private]


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