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.

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

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         if(!physListFactory) physListFactory = new G4PhysListFactory();
00136         if(isLocalPhysicsList(argv[0]))
00137                 physicsList = getLocalPhysicsList(argv[0]);
00138         if(!physicsList) {
00139                 if(!physListFactory->IsReferencePhysList(argv[0])) {
00140                         printError("physics: Unknown physics list '%s'",
00141                                                         argv[0].c_str());
00142                         printError("The choices are:");
00143                         listPhysicsLists();
00144                         return -1;
00145                 }
00146                 physicsList = physListFactory->GetReferencePhysList(argv[0]);
00147         }
00148 
00149         handleNamedArgs(namedArgs);
00150 
00151         disableList = splitString(disable,',',true);
00152 
00153         BLManager::getObject()->registerPhysics(this);
00154 
00155         // now is too early to manage the physics processes; do it in a callback
00156         BLManager::getObject()->registerCallback(this,0);
00157 
00158         print(argv[0]);
00159 
00160         return 0;
00161 }

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.

00164 {
00165         argString(disable,"disable","Comma-separated list of processes to disable (e.g. 'Decay,msc').");
00166         argString(disable,"inactivate","Synonym of disable.");
00167         argString(disable,"deactivate","Synonym of disable.");
00168         argInt(doStochastics,"doStochastics","Set to zero to disable all stochastic processes.");
00169         argDouble(minRangeCut,"minRangeCut","Minimum range cut for particle production (1 mm)");
00170         argInt(list,"list","Nonzero to list the processes (later in output).");
00171 }

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

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

Reimplemented from BLCommand.

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

00174 {
00175         BLCommand::help(detailed);
00176         if(detailed) {
00177                 printf("\n             ----- PHYSICS LISTS -----\n\n");
00178                 printf("             Further guidance in selecting physics lists is available at:\n"
00179                 "http://geant4.web.cern.ch/geant4/support/physicsLists/referencePL/index.shtml\n\n");
00180                 printf("             The default all-around physics list for HEP is called 'QGSP_BERT'.\n\n");
00181                 printf("             LHEP uses exclusively parameterized modeling.\n");
00182                 printf("             FTF lists use the FRITIOF description of string excitation\n"
00183                        "                  and fragmentation.\n");
00184                 printf("             QGSP lists use the quark gluon string model\n");
00185                 printf("             QGSC are as QGSP except applying CHIPS modeling for the nuclear\n"
00186                        "                  de-excitation.\n");
00187                 printf("             QGSP_BERT is like QGSP, but using Geant4 Bertini cascade.\n");
00188                 printf("             QGSP_BIC is like QGSP, but using Geant4 Binary cascade.\n");
00189                 printf("             _EMV suffix indicates a faster but less accurate EM modeling.\n");
00190                 printf("             _HP suffix uses the data driven high precision neutron package\n"
00191                        "                  (thermal to 20 MeV).\n");
00192                 printf("             _NQE suffix indicates a list for comparison with earlier release.\n");
00193                 printf("\n             List of available physics lists:\n");
00194                 listPhysicsLists();
00195         }
00196 }

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

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

Reimplemented from BLCallback.

00259 {
00260         G4ParticleTable::G4PTblDicIterator *theParticleIterator =
00261                         G4ParticleTable::GetParticleTable()->GetIterator();
00262         theParticleIterator -> reset();
00263         while( (*theParticleIterator)() ) {
00264                 G4ParticleDefinition* particle = theParticleIterator->value();
00265                 G4ProcessManager* pmanager = particle -> GetProcessManager();
00266                 if(!pmanager) {
00267                         printf("BLCMDphysics: particle '%s' has no ProcessManager!\n",
00268                                 particle->GetParticleName().c_str());
00269                         continue;
00270                 }
00271                 pmanager -> AddProcess(new G4StepLimiter(),-1,-1,3);
00272         }
00273 }

bool BLCMDphysics::isLocalPhysicsList ( G4String  name  ) 

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

Referenced by command().

00276 {
00277         if(name == "MICEPhysicsList") return true;
00278         if(name == "LISAPhysicsList") return true;
00279         return false;
00280 }

G4VUserPhysicsList * BLCMDphysics::getLocalPhysicsList ( G4String  name  ) 

getLocalPhysicsList() returns a local physics list.

Referenced by command().

00283 {
00284         if(name == "MICEPhysicsList") return new MICEPhysicsList();
00285         if(name == "LISAPhysicsList") return new LISAPhysicsList();
00286         return 0;
00287 }

void BLCMDphysics::listPhysicsLists (  ) 

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

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

Referenced by command(), and help().

00291 {
00292         if(!physListFactory) physListFactory = new G4PhysListFactory();
00293         const std::vector<G4String> v = physListFactory->AvailablePhysLists();
00294         for(unsigned i=0; i<v.size(); ++i) {
00295                 G4String n = v[i];
00296                 G4String s = getSynopsis(n);
00297                 printCase(n,s);
00298         }
00299         printCase("MICEPhysicsList",getSynopsis("MICEPhysicsList"));
00300         printCase("LISAPhysicsList",getSynopsis("LISAPhysicsList"));
00301 }

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

printCase() prints the description of a physics list.

References BLCommand::wrapWords().

Referenced by listPhysicsLists().

00304 {
00305         G4String indent1("             "), indent("       ");
00306         indent += indent1;
00307         indent1 += name;
00308         do { indent1 += " "; } while(indent1.size() < indent.size());
00309         printf("%s",wrapWords(description,indent1,indent).c_str());
00310 }

G4String BLCMDphysics::getSynopsis ( G4String  listName  ) 

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

References CASE.

Referenced by listPhysicsLists().

00313 {
00314 #define CASE(NAME,SYNOPSIS)                             \
00315         if(listName == #NAME) return SYNOPSIS;
00316 again:
00317 #ifndef PHYSICS_NO_USE_CASES
00318         CASE(FTFP,"For calorimetry -- is similar to QGSP for the treatment of the fragmentation, but the string excitation/fragmentation is changed from quark-gluon string model to a diffractive string excitation similar to that in FRITJOF, and the Lund fragmentation functions.");
00319         CASE(LBE,"For low background experiments (e.g. underground)");
00320         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.");
00321         CASE(LHEP_BERT,"For 'typical' HEP collider detectors and low-energy penetration shielding (< 5 GeV) -- is as LHEP_GN 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.");
00322         CASE(LHEP_PRECO_HP,"For low-energy penetration shielding (< 5 GeV) -- is like LHEP_PRECO, but uses point-wise evaluated cross-section data to model neutron interactions from thermal neutron energies to ~20MeV. This applies to capture, elastic scattering, fission, and inelastic scattering. Again fission fragments are not available. Please use together with G4NDL3.5 or higher.");
00323         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.");
00324         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.");
00325         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.");
00326         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.");
00327 #endif
00328         CASE(LISAPhysicsList,"Geant4 examples/advanced/cosmicray_charging");
00329         CASE(MICEPhysicsList,"The current MICE physics list.");
00330 
00331         return "No synopsis available.";
00332 }


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