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
G4int gammaToMuPair
std::vector< G4String > disableList
G4PhysListFactory * physListFactory

Constructor & Destructor Documentation

BLCMDphysics::BLCMDphysics (  ) 

Default Constructor.

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

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

BLCMDphysics::~BLCMDphysics (  )  [inline]

Destructor.

00064 { }


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.

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

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

getPhysicsList() returns the G4PhysicsList.

Implements BLPhysics.

References physicsList.

00071 { return physicsList; }

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

commandName() returns "physics".

Implements BLCommand.

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

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

void BLCMDphysics::defineNamedArgs (  )  [virtual]

defineNamedArgs() defines the named arguments.

Reimplemented from BLCommand.

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

00190 {
00191         argString(disable,"disable","Comma-separated list of processes to disable (e.g. 'Decay,msc').");
00192         argString(disable,"inactivate","Synonym of disable.");
00193         argString(disable,"deactivate","Synonym of disable.");
00194         argInt(doStochastics,"doStochastics","Set to zero to disable all stochastic processes.");
00195         argDouble(minRangeCut,"minRangeCut","Minimum range cut for particle production (1 mm)");
00196         argInt(list,"list","Nonzero to list the processes (later in output).");
00197         argInt(gammaToMuPair,"gammaToMuPair","Nonzero to add gamma->mu+mu- (0).");
00198         argInt(synchrotronRadiation,"synchrotronRadiation","Nonzero to add synchrotron radiation to the physics list for e- and e+.");
00199         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.");
00200 }

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

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

Reimplemented from BLCommand.

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

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

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, gammaToMuPair, synchrotronRadiation, and synchrotronRadiationMuon.

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

bool BLCMDphysics::isLocalPhysicsList ( G4String  name  ) 

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

Referenced by command().

00343 {
00344 #ifdef MICEANDLISALISTS
00345         if(name == "MICEPhysicsList") return true;
00346         if(name == "LISAPhysicsList") return true;
00347 #endif
00348         return false;
00349 }

G4VUserPhysicsList * BLCMDphysics::getLocalPhysicsList ( G4String  name  ) 

getLocalPhysicsList() returns a local physics list.

Referenced by command().

00352 {
00353 #ifdef MICEANDLISALISTS
00354         if(name == "MICEPhysicsList") return new MICEPhysicsList();
00355         if(name == "LISAPhysicsList") return new LISAPhysicsList();
00356 #endif
00357         return 0;
00358 }

void BLCMDphysics::listPhysicsLists (  ) 

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

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

Referenced by command(), and help().

00362 {
00363         if(!physListFactory) physListFactory = new G4PhysListFactory();
00364         const std::vector<G4String> v = physListFactory->AvailablePhysLists();
00365         for(unsigned i=0; i<v.size(); ++i) {
00366                 G4String n = v[i];
00367                 G4String s = getSynopsis(n);
00368                 printCase(n,s);
00369         }
00370 #ifdef MICEANDLISALISTS
00371         printCase("MICEPhysicsList",getSynopsis("MICEPhysicsList"));
00372         printCase("LISAPhysicsList",getSynopsis("LISAPhysicsList"));
00373 #endif
00374 }

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

printCase() prints the description of a physics list.

References BLCommand::wrapWords().

Referenced by listPhysicsLists().

00377 {
00378         G4String indent1("             "), indent("       ");
00379         indent += indent1;
00380         indent1 += name;
00381         do { indent1 += " "; } while(indent1.size() < indent.size());
00382         printf("%s",wrapWords(description,indent1,indent).c_str());
00383 }

G4String BLCMDphysics::getSynopsis ( G4String  listName  ) 

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

References CASE.

Referenced by listPhysicsLists().

00386 {
00387 #define CASE(NAME,SYNOPSIS)                             \
00388         if(listName == #NAME) return SYNOPSIS;
00389 
00390         CASE(CHIPS,"No synopsis available.");
00391         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.");
00392         CASE(FTFP_BERT_EMV,"Like FTFP_BERT but with faster EM modeling.");
00393         CASE(FTFP_BERT_EMX,"Like FTFP_BERT but with low-energy EM modeling.");
00394         CASE(FTFP_BERT_TRV,"A variant of QGSP_BERT where the Geant4 Bertini cascade is only used for particles below ~5.5 GeV.");
00395         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.");
00396         CASE(LBE,"For low background experiments (e.g. underground)");
00397         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.");
00398         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.");
00399         CASE(LHEP_EMV,"Like LHEP but with faster EM modeling.");
00400         CASE(QBBC,"No synopsis available.");
00401         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.");
00402         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.");
00403         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.");
00404         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.");
00405         CASE(QGSP_BERT_EMV,"Like QGSP_BERT but with faster EM modeling.");
00406         CASE(QGSP_BERT_EMX,"Like QGSP_BERT but with low-energy EM modeling.");
00407         CASE(QGSP_BERT_HP,"Like QGSP_BERT but with _HP modeling for neutrons.");
00408         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.");
00409         CASE(QGSP_BIC_EMX,"Like QGSP_BIC but with low-energy EM modeling.");
00410         CASE(QGSP_BIC_HP,"Like QGSP_BIC but with _HP modeling for neutrons.");
00411 
00412 #ifdef MICEANDLISALISTS
00413         CASE(MICEPhysicsList,"The current MICE physics list.");
00414         CASE(LISAPhysicsList,"Geant4 examples/advanced/cosmicray_charging");
00415 #endif
00416 
00417         return "No synopsis available.";
00418 }


Member Data Documentation

G4VUserPhysicsList* BLCMDphysics::physicsList [private]

G4int BLCMDphysics::doStochastics [private]

G4String BLCMDphysics::disable [private]

G4int BLCMDphysics::list [private]

G4int BLCMDphysics::gammaToMuPair [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