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
G4int synchrotronRadiation
G4int synchrotronRadiationMuon
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(), BLCommand::setSynopsis(), synchrotronRadiation, and synchrotronRadiationMuon.

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

BLCMDphysics::~BLCMDphysics (  )  [inline]

Destructor.

00060 { }


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.

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

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

getPhysicsList() returns the G4PhysicsList.

Implements BLPhysics.

References physicsList.

00067 { return physicsList; }

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

commandName() returns "physics".

Implements BLCommand.

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

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

void BLCMDphysics::defineNamedArgs (  )  [virtual]

defineNamedArgs() defines the named arguments.

Reimplemented from BLCommand.

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

00185 {
00186         argString(disable,"disable","Comma-separated list of processes to disable (e.g. 'Decay,msc').");
00187         argString(disable,"inactivate","Synonym of disable.");
00188         argString(disable,"deactivate","Synonym of disable.");
00189         argInt(doStochastics,"doStochastics","Set to zero to disable all stochastic processes.");
00190         argDouble(minRangeCut,"minRangeCut","Minimum range cut for particle production (1 mm)");
00191         argInt(list,"list","Nonzero to list the processes (later in output).");
00192         argInt(synchrotronRadiation,"synchrotronRadiation","Nonzero to add synchrotron radiation to the physics list for e- and e+.");
00193         argInt(synchrotronRadiationMuon,"synchrotronRadiationMuon","Nonzero to add synchrotron radiation to the physics list for mu- and mu+ NOTE: This is experimental, and may not work correctly.");
00194 }

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

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

Reimplemented from BLCommand.

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

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

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

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

ensure only mu-nu is present

ensure only mu-nu is present

Reimplemented from BLCallback.

References BLAssert, synchrotronRadiation, and synchrotronRadiationMuon.

00283 {
00284         G4ParticleTable::G4PTblDicIterator *theParticleIterator =
00285                         G4ParticleTable::GetParticleTable()->GetIterator();
00286         theParticleIterator -> reset();
00287         while( (*theParticleIterator)() ) {
00288                 G4ParticleDefinition* particle = theParticleIterator->value();
00289                 G4ProcessManager* pmanager = particle -> GetProcessManager();
00290                 if(!pmanager) {
00291                         printf("BLCMDphysics: particle '%s' has no ProcessManager!\n",
00292                                 particle->GetParticleName().c_str());
00293                         continue;
00294                 }
00295                 pmanager -> AddProcess(new G4StepLimiter(),-1,-1,3);
00296         }
00297 
00298         // add pi+/pi- decay to e nu
00299         double br=1.230e-4;
00300         G4ParticleTable *table = G4ParticleTable::GetParticleTable();
00301         { G4ParticleDefinition *pi=table->FindParticle("pi+");
00302           G4DecayTable* table = pi->GetDecayTable();
00303           BLAssert(table->entries() == 1); /// ensure only mu-nu is present
00304           G4VDecayChannel *munu = table->GetDecayChannel(0);
00305           munu->SetBR(1.0-br);  
00306           table->Insert(new G4PhaseSpaceDecayChannel("pi+",br,2,"e+","nu_e"));
00307         }
00308         { G4ParticleDefinition *pi=table->FindParticle("pi-");
00309           G4DecayTable* table = pi->GetDecayTable();
00310           BLAssert(table->entries() == 1); /// ensure only mu-nu is present
00311           G4VDecayChannel *munu = table->GetDecayChannel(0);
00312           munu->SetBR(1.0-br);  
00313           table->Insert(new G4PhaseSpaceDecayChannel("pi-",br,2,"e-","anti_nu_e"));
00314         }
00315 
00316         if(synchrotronRadiation) {
00317                 table->FindParticle("e-")->GetProcessManager()->
00318                         AddDiscreteProcess(new G4SynchrotronRadiation);
00319                 table->FindParticle("e+")->GetProcessManager()->
00320                         AddDiscreteProcess(new G4SynchrotronRadiation);
00321         }
00322 
00323         if(synchrotronRadiationMuon) {
00324                 table->FindParticle("mu-")->GetProcessManager()->
00325                         AddDiscreteProcess(new G4SynchrotronRadiation);
00326                 table->FindParticle("mu+")->GetProcessManager()->
00327                         AddDiscreteProcess(new G4SynchrotronRadiation);
00328         }
00329 }

bool BLCMDphysics::isLocalPhysicsList ( G4String  name  ) 

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

Referenced by command().

00332 {
00333         if(name == "MICEPhysicsList") return true;
00334         if(name == "LISAPhysicsList") return true;
00335         return false;
00336 }

G4VUserPhysicsList * BLCMDphysics::getLocalPhysicsList ( G4String  name  ) 

getLocalPhysicsList() returns a local physics list.

Referenced by command().

00339 {
00340         if(name == "MICEPhysicsList") return new MICEPhysicsList();
00341         if(name == "LISAPhysicsList") return new LISAPhysicsList();
00342         return 0;
00343 }

void BLCMDphysics::listPhysicsLists (  ) 

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

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

Referenced by command(), and help().

00347 {
00348         if(!physListFactory) physListFactory = new G4PhysListFactory();
00349         const std::vector<G4String> v = physListFactory->AvailablePhysLists();
00350         for(unsigned i=0; i<v.size(); ++i) {
00351                 G4String n = v[i];
00352                 G4String s = getSynopsis(n);
00353                 printCase(n,s);
00354         }
00355         printCase("MICEPhysicsList",getSynopsis("MICEPhysicsList"));
00356         printCase("LISAPhysicsList",getSynopsis("LISAPhysicsList"));
00357 }

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

printCase() prints the description of a physics list.

References BLCommand::wrapWords().

Referenced by listPhysicsLists().

00360 {
00361         G4String indent1("             "), indent("       ");
00362         indent += indent1;
00363         indent1 += name;
00364         do { indent1 += " "; } while(indent1.size() < indent.size());
00365         printf("%s",wrapWords(description,indent1,indent).c_str());
00366 }

G4String BLCMDphysics::getSynopsis ( G4String  listName  ) 

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

References CASE.

Referenced by listPhysicsLists().

00369 {
00370 #define CASE(NAME,SYNOPSIS)                             \
00371         if(listName == #NAME) return SYNOPSIS;
00372 
00373         CASE(CHIPS,"No synopsis available.");
00374         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.");
00375         CASE(FTFP_BERT_EMV,"Like FTFP_BERT but with faster EM modeling.");
00376         CASE(FTFP_BERT_EMX,"Like FTFP_BERT but with low-energy EM modeling.");
00377         CASE(FTFP_BERT_TRV,"A variant of QGSP_BERT where the Geant4 Bertini cascade is only used for particles below ~5.5 GeV.");
00378         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.");
00379         CASE(LBE,"For low background experiments (e.g. underground)");
00380         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.");
00381         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.");
00382         CASE(LHEP_EMV,"Like LHEP but with faster EM modeling.");
00383         CASE(QBBC,"No synopsis available.");
00384         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.");
00385         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.");
00386         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.");
00387         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.");
00388         CASE(QGSP_BERT_EMV,"Like QGSP_BERT but with faster EM modeling.");
00389         CASE(QGSP_BERT_EMX,"Like QGSP_BERT but with low-energy EM modeling.");
00390         CASE(QGSP_BERT_HP,"Like QGSP_BERT but with _HP modeling for neutrons.");
00391         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.");
00392         CASE(QGSP_BIC_EMX,"Like QGSP_BIC but with low-energy EM modeling.");
00393         CASE(QGSP_BIC_HP,"Like QGSP_BIC but with _HP modeling for neutrons.");
00394 
00395         CASE(LISAPhysicsList,"Geant4 examples/advanced/cosmicray_charging");
00396         CASE(MICEPhysicsList,"The current MICE physics list.");
00397 
00398         return "No synopsis available.";
00399 }


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