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 |
BLCMDphysics::BLCMDphysics | ( | ) |
Default Constructor.
References BLCMDTYPE_AUX, disable, doStochastics, list, BLPhysics::minRangeCut, physicsList, physListFactory, BLCommand::registerCommand(), BLCommand::setDescription(), and BLCommand::setSynopsis().
00099 : BLPhysics(), BLCommand(), disableList() 00100 { 00101 registerCommand(BLCMDTYPE_AUX); 00102 setSynopsis("Defines the physics processes and controls them."); 00103 setDescription("Exactly one physics command must be present.\n" 00104 "This command implements the geant4 physics lists\n" 00105 "of physics processes.\n" 00106 "The command is 'physics QGSP' for the QGSP set, and\n" 00107 "similarly for the other sets. With no argument it prints\n" 00108 "the available physics lists.\n" 00109 "Note that stochastic processes are always disabled while\n" 00110 "tracking the tune and reference particles.\n" 00111 "The only non-stochastic processes are Transportation and\n" 00112 "ionization energy loss (with fluctuations disabled).\n" 00113 "For muon beam studies you may want to disable the Decay " 00114 "process.\n\n" 00115 "NOTE: this command defines the particles used throughout " 00116 "the simulation, so this command must come before others " 00117 "that use particle names.\n\n" 00118 "NOTE: the rare decy mode for pi+/pi- to e nu is added.\n\n" 00119 "The default all-around physics list for HEP is called " 00120 "'QGSP_BERT'."); 00121 00122 physicsList = 0; 00123 doStochastics = 1; 00124 list = 0; 00125 disable = ""; 00126 minRangeCut = 1.0 * mm; 00127 physListFactory = 0; 00128 }
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.
00215 { 00216 if(value == NORMAL) { 00217 if(doStochastics) 00218 value = FORCE_ON; 00219 else 00220 value = FORCE_OFF; 00221 } 00222 00223 G4UImanager* UI = G4UImanager::GetUIpointer(); 00224 G4String cmd; 00225 char rangeCmd[64]; 00226 switch(value) { 00227 case FORCE_OFF: 00228 // turn off delta rays (and other particle production) by 00229 // setting production cut very high in energy by setting 00230 // the required range large. 00231 UI->ApplyCommand("/run/setCut 1 km"); 00232 // Turn off enegry-loss straggling 00233 UI->ApplyCommand("/process/eLoss/fluct false"); 00234 cmd = "/process/inactivate "; 00235 if(!quiet) printf("*** NOTE: Stochastic processes are disabled!\n"); 00236 stochasticsEnabled = false; 00237 break; 00238 case FORCE_ON: 00239 // turn on delta rays (and other particle production) 00240 sprintf(rangeCmd,"/run/setCut %.6g mm",minRangeCut/mm); 00241 UI->ApplyCommand(rangeCmd); 00242 // Turn on enegry-loss straggling 00243 UI->ApplyCommand("/process/eLoss/fluct true"); 00244 cmd = "/process/activate "; 00245 if(!quiet) printf("Stochastic processes are enabled.\n"); 00246 stochasticsEnabled = true; 00247 break; 00248 case NORMAL: 00249 return; 00250 } 00251 00252 // run through list of processes, applying cmd to stochastic ones 00253 G4ProcessTable *pt = G4ProcessTable::GetProcessTable(); 00254 G4ProcessTable::G4ProcNameVector *nv = pt->GetNameList(); 00255 if(list) printf("List of physics processes:\n"); 00256 for(unsigned int i=0; i<nv->size(); ++i) { 00257 G4String name = (*nv)[i]; 00258 if(list) printf("\t%s\n",name.c_str()); 00259 if(isStochasticProcess(name)) 00260 UI->ApplyCommand(cmd+name); 00261 } 00262 if(list) printf("\n"); 00263 list = 0; // print list only once 00264 00265 // now disable all proceses in the disableList 00266 cmd = "/process/inactivate "; 00267 for(unsigned i=0; i<disableList.size(); ++i) { 00268 if(!quiet) printf("*** NOTE: Disabling physics process '%s'\n", 00269 disableList[i].c_str()); 00270 UI->ApplyCommand(cmd+disableList[i]); 00271 } 00272 }
G4VUserPhysicsList* BLCMDphysics::getPhysicsList | ( | ) | [inline, virtual] |
getPhysicsList() returns the G4PhysicsList.
Implements BLPhysics.
References physicsList.
00063 { return physicsList; }
G4String BLCMDphysics::commandName | ( | ) | [inline, virtual] |
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().
00131 { 00132 if(physicsList != 0 && argv.size() > 0) { 00133 printError("physics: Multiple physics commands not allowed!"); 00134 return -1; 00135 } 00136 00137 if(argv.size() == 0) 00138 return -1; 00139 00140 // handle obsolete physics lists 00141 if(argv[0] == "QGSC") { 00142 argv[0] = "QGSC_BERT"; 00143 printf("Physics list QGSC is obsolete, replaced by QGSC_BERT\n"); 00144 } 00145 if(argv[0] == "QGSP") { 00146 argv[0] = "QGSP_BERT"; 00147 printf("Physics list QGSP is obsolete, replaced by QGSP_BERT\n"); 00148 } 00149 00150 if(!physListFactory) physListFactory = new G4PhysListFactory(); 00151 if(isLocalPhysicsList(argv[0])) 00152 physicsList = getLocalPhysicsList(argv[0]); 00153 if(!physicsList) { 00154 if(!physListFactory->IsReferencePhysList(argv[0])) { 00155 printError("physics: Unknown physics list '%s'", 00156 argv[0].c_str()); 00157 printError("The choices are:"); 00158 listPhysicsLists(); 00159 return -1; 00160 } 00161 physicsList = physListFactory->GetReferencePhysList(argv[0]); 00162 } 00163 00164 handleNamedArgs(namedArgs); 00165 00166 disableList = splitString(disable,',',true); 00167 00168 BLManager::getObject()->registerPhysics(this); 00169 00170 // now is too early to manage the physics processes; do it in a callback 00171 BLManager::getObject()->registerCallback(this,0); 00172 00173 print(argv[0]); 00174 00175 return 0; 00176 }
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.
00179 { 00180 argString(disable,"disable","Comma-separated list of processes to disable (e.g. 'Decay,msc')."); 00181 argString(disable,"inactivate","Synonym of disable."); 00182 argString(disable,"deactivate","Synonym of disable."); 00183 argInt(doStochastics,"doStochastics","Set to zero to disable all stochastic processes."); 00184 argDouble(minRangeCut,"minRangeCut","Minimum range cut for particle production (1 mm)"); 00185 argInt(list,"list","Nonzero to list the processes (later in output)."); 00186 }
void BLCMDphysics::help | ( | bool | detailed | ) | [virtual] |
help() provides special help text for the physics command.
Reimplemented from BLCommand.
References BLCommand::help(), and listPhysicsLists().
00189 { 00190 BLCommand::help(detailed); 00191 if(detailed) { 00192 printf("\n ----- PHYSICS LISTS -----\n\n"); 00193 printf(" Further guidance in selecting physics lists is available at:\n" 00194 "http://geant4.web.cern.ch/geant4/support/physicsLists/referencePL/index.shtml\n\n"); 00195 printf(" The default all-around physics list for HEP is called 'QGSP_BERT'.\n\n"); 00196 printf(" LHEP uses exclusively parameterized modeling.\n"); 00197 printf(" FTF lists use the FRITIOF description of string excitation\n" 00198 " and fragmentation.\n"); 00199 printf(" QGSP lists use the quark gluon string model\n"); 00200 printf(" QGSC are as QGSP except applying CHIPS modeling for the nuclear\n" 00201 " de-excitation.\n"); 00202 printf(" _BERT uses Geant4 Bertini cascade below ~ 10 GeV.\n"); 00203 printf(" _BIC uses Geant4 Binary cascade below ~ 10 GeV.\n"); 00204 printf(" _EMV suffix indicates a faster but less accurate EM modeling.\n"); 00205 printf(" _EMX suffix indicates the low-energy EM processes.\n"); 00206 printf(" _HP suffix uses the data driven high precision neutron package\n" 00207 " (thermal to 20 MeV).\n"); 00208 printf(" _NQE suffix indicates a list for comparison with earlier release.\n"); 00209 printf("\n List of available physics lists:\n"); 00210 listPhysicsLists(); 00211 } 00212 }
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.
00275 { 00276 G4ParticleTable::G4PTblDicIterator *theParticleIterator = 00277 G4ParticleTable::GetParticleTable()->GetIterator(); 00278 theParticleIterator -> reset(); 00279 while( (*theParticleIterator)() ) { 00280 G4ParticleDefinition* particle = theParticleIterator->value(); 00281 G4ProcessManager* pmanager = particle -> GetProcessManager(); 00282 if(!pmanager) { 00283 printf("BLCMDphysics: particle '%s' has no ProcessManager!\n", 00284 particle->GetParticleName().c_str()); 00285 continue; 00286 } 00287 pmanager -> AddProcess(new G4StepLimiter(),-1,-1,3); 00288 } 00289 00290 // add pi+/pi- decay to e nu 00291 double br=1.230e-4; 00292 G4ParticleTable *table = G4ParticleTable::GetParticleTable(); 00293 { G4ParticleDefinition *pi=table->FindParticle("pi+"); 00294 G4DecayTable* table = pi->GetDecayTable(); 00295 assert(table->entries() == 1); /// ensure only mu-nu is present 00296 G4VDecayChannel *munu = table->GetDecayChannel(0); 00297 munu->SetBR(1.0-br); 00298 table->Insert(new G4PhaseSpaceDecayChannel("pi+",br,2,"e+","nu_e")); 00299 } 00300 { G4ParticleDefinition *pi=table->FindParticle("pi-"); 00301 G4DecayTable* table = pi->GetDecayTable(); 00302 assert(table->entries() == 1); /// ensure only mu-nu is present 00303 G4VDecayChannel *munu = table->GetDecayChannel(0); 00304 munu->SetBR(1.0-br); 00305 table->Insert(new G4PhaseSpaceDecayChannel("pi-",br,2,"e-","anti_nu_e")); 00306 } 00307 }
bool BLCMDphysics::isLocalPhysicsList | ( | G4String | name | ) |
isLocalPhysicsList() returns true if thsi is a local physics list.
Referenced by command().
00310 { 00311 if(name == "MICEPhysicsList") return true; 00312 if(name == "LISAPhysicsList") return true; 00313 return false; 00314 }
G4VUserPhysicsList * BLCMDphysics::getLocalPhysicsList | ( | G4String | name | ) |
getLocalPhysicsList() returns a local physics list.
Referenced by command().
00317 { 00318 if(name == "MICEPhysicsList") return new MICEPhysicsList(); 00319 if(name == "LISAPhysicsList") return new LISAPhysicsList(); 00320 return 0; 00321 }
void BLCMDphysics::listPhysicsLists | ( | ) |
listPhysicsLists() will list all known physics lists, with synopsis.
References getSynopsis(), physListFactory, and printCase().
Referenced by command(), and help().
00325 { 00326 if(!physListFactory) physListFactory = new G4PhysListFactory(); 00327 const std::vector<G4String> v = physListFactory->AvailablePhysLists(); 00328 for(unsigned i=0; i<v.size(); ++i) { 00329 G4String n = v[i]; 00330 G4String s = getSynopsis(n); 00331 printCase(n,s); 00332 } 00333 printCase("MICEPhysicsList",getSynopsis("MICEPhysicsList")); 00334 printCase("LISAPhysicsList",getSynopsis("LISAPhysicsList")); 00335 }
void BLCMDphysics::printCase | ( | G4String | name, | |
G4String | description | |||
) |
printCase() prints the description of a physics list.
References BLCommand::wrapWords().
Referenced by listPhysicsLists().
00338 { 00339 G4String indent1(" "), indent(" "); 00340 indent += indent1; 00341 indent1 += name; 00342 do { indent1 += " "; } while(indent1.size() < indent.size()); 00343 printf("%s",wrapWords(description,indent1,indent).c_str()); 00344 }
G4String BLCMDphysics::getSynopsis | ( | G4String | listName | ) |
getSynopsis() returns the synopsis for a physics list, by name.
References CASE.
Referenced by listPhysicsLists().
00347 { 00348 #define CASE(NAME,SYNOPSIS) \ 00349 if(listName == #NAME) return SYNOPSIS; 00350 00351 CASE(CHIPS,"No synopsis available."); 00352 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."); 00353 CASE(FTFP_BERT_EMV,"Like FTFP_BERT but with faster EM modeling."); 00354 CASE(FTFP_BERT_EMX,"Like FTFP_BERT but with low-energy EM modeling."); 00355 CASE(FTFP_BERT_TRV,"A variant of QGSP_BERT where the Geant4 Bertini cascade is only used for particles below ~5.5 GeV."); 00356 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."); 00357 CASE(LBE,"For low background experiments (e.g. underground)"); 00358 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."); 00359 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."); 00360 CASE(LHEP_EMV,"Like LHEP but with faster EM modeling."); 00361 CASE(QBBC,"No synopsis available."); 00362 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."); 00363 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."); 00364 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."); 00365 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."); 00366 CASE(QGSP_BERT_EMV,"Like QGSP_BERT but with faster EM modeling."); 00367 CASE(QGSP_BERT_EMX,"Like QGSP_BERT but with low-energy EM modeling."); 00368 CASE(QGSP_BERT_HP,"Like QGSP_BERT but with _HP modeling for neutrons."); 00369 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."); 00370 CASE(QGSP_BIC_EMX,"Like QGSP_BIC but with low-energy EM modeling."); 00371 CASE(QGSP_BIC_HP,"Like QGSP_BIC but with _HP modeling for neutrons."); 00372 00373 CASE(LISAPhysicsList,"Geant4 examples/advanced/cosmicray_charging"); 00374 CASE(MICEPhysicsList,"The current MICE physics list."); 00375 00376 return "No synopsis available."; 00377 }
G4VUserPhysicsList* BLCMDphysics::physicsList [private] |
Referenced by BLCMDphysics(), command(), and getPhysicsList().
G4int BLCMDphysics::doStochastics [private] |
Referenced by BLCMDphysics(), defineNamedArgs(), and setDoStochastics().
G4String BLCMDphysics::disable [private] |
Referenced by BLCMDphysics(), command(), and defineNamedArgs().
G4int BLCMDphysics::list [private] |
Referenced by BLCMDphysics(), defineNamedArgs(), and setDoStochastics().
std::vector<G4String> BLCMDphysics::disableList [private] |
Referenced by command(), and setDoStochastics().
G4PhysListFactory* BLCMDphysics::physListFactory [private] |
Referenced by BLCMDphysics(), command(), and listPhysicsLists().