00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050 #include "G4VUserPhysicsList.hh"
00051
00052 #include "globals.hh"
00053 #include "G4SystemOfUnits.hh"
00054 #include "G4PhysicsListHelper.hh"
00055 #include "G4ParticleDefinition.hh"
00056 #include "G4ProcessManager.hh"
00057 #include "G4ParticleTable.hh"
00058 #include "G4ProductionCutsTable.hh"
00059 #include "G4Material.hh"
00060 #include "G4UserPhysicsListMessenger.hh"
00061 #include "G4UImanager.hh"
00062 #include "G4UnitsTable.hh"
00063 #include "G4RegionStore.hh"
00064 #include "G4Region.hh"
00065 #include "G4ProductionCutsTable.hh"
00066 #include "G4ProductionCuts.hh"
00067 #include "G4MaterialCutsCouple.hh"
00068
00069 #include "G4ios.hh"
00070 #include <iomanip>
00071 #include <fstream>
00072
00074 G4VUserPhysicsList::G4VUserPhysicsList()
00075 :verboseLevel(1),
00076 defaultCutValue(1.0 * mm),
00077 isSetDefaultCutValue(false),
00078 fRetrievePhysicsTable(false),
00079 fStoredInAscii(true),
00080 fIsCheckedForRetrievePhysicsTable(false),
00081 fIsRestoredCutValues(false),
00082 directoryPhysicsTable("."),
00083 fDisplayThreshold(0),
00084 fIsPhysicsTableBuilt(false),
00085 fDisableCheckParticleList(false)
00086 {
00087
00088 defaultCutValue = 1.0*mm;
00089
00090
00091 theParticleTable = G4ParticleTable::GetParticleTable();
00092 theParticleIterator = theParticleTable->GetIterator();
00093
00094
00095 fCutsTable = G4ProductionCutsTable::GetProductionCutsTable();
00096
00097
00098 fCutsTable->SetEnergyRange(0.99*keV, 100*TeV);
00099
00100
00101 theMessenger = new G4UserPhysicsListMessenger(this);
00102
00103
00104 thePLHelper = G4PhysicsListHelper::GetPhysicsListHelper();
00105 thePLHelper->SetVerboseLevel(verboseLevel);
00106
00107 }
00108
00110 G4VUserPhysicsList::~G4VUserPhysicsList()
00111 {
00112 if (theMessenger != 0) {
00113 delete theMessenger;
00114 theMessenger = 0;
00115 }
00116 RemoveProcessManager();
00117
00118
00119 theParticleTable->DeleteAllParticles();
00120
00121 }
00122
00124 G4VUserPhysicsList::G4VUserPhysicsList(const G4VUserPhysicsList& right)
00125 :verboseLevel(right.verboseLevel),
00126 defaultCutValue(right.defaultCutValue),
00127 isSetDefaultCutValue(right.isSetDefaultCutValue),
00128 fRetrievePhysicsTable(right.fRetrievePhysicsTable),
00129 fStoredInAscii(right.fStoredInAscii),
00130 fIsCheckedForRetrievePhysicsTable(right.fIsCheckedForRetrievePhysicsTable),
00131 fIsRestoredCutValues(right.fIsRestoredCutValues),
00132 directoryPhysicsTable(right.directoryPhysicsTable),
00133 fDisplayThreshold(right.fDisplayThreshold),
00134 fIsPhysicsTableBuilt(right.fIsPhysicsTableBuilt),
00135 fDisableCheckParticleList(right.fDisableCheckParticleList)
00136 {
00137
00138 theParticleTable = G4ParticleTable::GetParticleTable();
00139 theParticleIterator = theParticleTable->GetIterator();
00140
00141
00142 fCutsTable = G4ProductionCutsTable::GetProductionCutsTable();
00143
00144
00145 theMessenger = new G4UserPhysicsListMessenger(this);
00146
00147
00148 thePLHelper = G4PhysicsListHelper::GetPhysicsListHelper();
00149 thePLHelper->SetVerboseLevel(verboseLevel);
00150
00151 }
00152
00153
00155 G4VUserPhysicsList & G4VUserPhysicsList::operator=(const G4VUserPhysicsList & right)
00156 {
00157 if (this != &right) {
00158 verboseLevel = right.verboseLevel;
00159 defaultCutValue = right.defaultCutValue;
00160 isSetDefaultCutValue = right.isSetDefaultCutValue;
00161 fRetrievePhysicsTable = right.fRetrievePhysicsTable;
00162 fStoredInAscii = right.fStoredInAscii;
00163 fIsCheckedForRetrievePhysicsTable = right.fIsCheckedForRetrievePhysicsTable;
00164 fIsRestoredCutValues = right.fIsRestoredCutValues;
00165 directoryPhysicsTable = right.directoryPhysicsTable;
00166 fDisplayThreshold = right.fDisplayThreshold;
00167 fIsPhysicsTableBuilt = right.fIsPhysicsTableBuilt;
00168 fDisableCheckParticleList = right.fDisableCheckParticleList;
00169 }
00170 return *this;
00171 }
00172
00174 void G4VUserPhysicsList::AddProcessManager(G4ParticleDefinition* newParticle,
00175 G4ProcessManager* newManager)
00176 {
00177 if (newParticle == 0) return;
00178 if (newParticle->GetProcessManager() != 0) {
00179 #ifdef G4VERBOSE
00180 if (verboseLevel >1){
00181 G4cout << "G4VUserPhysicsList::AddProcessManager: "
00182 << newParticle->GetParticleName()
00183 << " already has ProcessManager " << G4endl;
00184 }
00185 #endif
00186 return;
00187 }
00188
00189
00190 if (newManager == 0){
00191
00192 if (newParticle->GetParticleType() == "nucleus") {
00193
00194 G4ParticleDefinition* genericIon =
00195 (G4ParticleTable::GetParticleTable())->FindParticle("GenericIon");
00196
00197 if (genericIon != 0) {
00198 G4ProcessManager* ionMan = genericIon->GetProcessManager();
00199 if (ionMan != 0) {
00200 newManager = new G4ProcessManager(*ionMan);
00201 } else {
00202
00203 newManager = new G4ProcessManager(newParticle);
00204 G4Exception("G4VUserPhysicsList::AddProcessManager",
00205 "Run0251", RunMustBeAborted,
00206 "GenericIon has no ProcessMamanger");
00207 }
00208 } else {
00209
00210 newManager = new G4ProcessManager(newParticle);
00211 G4Exception("G4VUserPhysicsList::AddProcessManager",
00212 "Run0252", RunMustBeAborted,
00213 "GenericIon does not exist");
00214 }
00215
00216 } else {
00217
00218 newManager = new G4ProcessManager(newParticle);
00219 }
00220 }
00221
00222
00223 newManager->SetParticleType(newParticle);
00224
00225
00226 newParticle->SetProcessManager(newManager);
00227
00228 #ifdef G4VERBOSE
00229 if (verboseLevel >2){
00230 G4cout << "G4VUserPhysicsList::AddProcessManager: "
00231 << "adds ProcessManager to "
00232 << newParticle->GetParticleName() << G4endl;
00233 newManager->DumpInfo();
00234 }
00235 #endif
00236 if ( fIsPhysicsTableBuilt
00237 && (newParticle->GetParticleType() == "nucleus")) {
00238 PreparePhysicsTable(newParticle);
00239 BuildPhysicsTable(newParticle);
00240 }
00241 }
00242
00243
00245 void G4VUserPhysicsList::InitializeProcessManager()
00246 {
00247
00248 theParticleIterator->reset();
00249 while( (*theParticleIterator)() ){
00250 G4ParticleDefinition* particle = theParticleIterator->value();
00251 G4ProcessManager* pmanager = particle->GetProcessManager();
00252 if (pmanager==0) {
00253
00254 pmanager = new G4ProcessManager(particle);
00255 particle->SetProcessManager(pmanager);
00256 }
00257 }
00258 }
00259
00261 void G4VUserPhysicsList::RemoveProcessManager()
00262 {
00263
00264 theParticleIterator->reset();
00265 while( (*theParticleIterator)() ){
00266 G4ParticleDefinition* particle = theParticleIterator->value();
00267 G4ProcessManager* pmanager = particle->GetProcessManager();
00268 if (pmanager!=0) delete pmanager;
00269 particle->SetProcessManager(0);
00270 #ifdef G4VERBOSE
00271 if (verboseLevel >2){
00272 G4cout << "G4VUserPhysicsList::RemoveProcessManager: "
00273 << "remove ProcessManager from "
00274 << particle->GetParticleName() << G4endl;
00275 }
00276 #endif
00277 }
00278 }
00279
00280
00282 void G4VUserPhysicsList::SetCuts()
00283 {
00284 if ( !isSetDefaultCutValue ){
00285 SetDefaultCutValue(defaultCutValue);
00286 }
00287
00288 #ifdef G4VERBOSE
00289 if (verboseLevel >1){
00290 G4cout << "G4VUserPhysicsList::SetCuts: " << G4endl;
00291 G4cout << "Cut for gamma: " << GetCutValue("gamma")/mm
00292 << "[mm]" << G4endl;
00293 G4cout << "Cut for e-: " << GetCutValue("e-")/mm
00294 << "[mm]" << G4endl;
00295 G4cout << "Cut for e+: " << GetCutValue("e+")/mm
00296 << "[mm]" << G4endl;
00297 G4cout << "Cut for proton: " << GetCutValue("proton")/mm
00298 << "[mm]" << G4endl;
00299 }
00300 #endif
00301
00302
00303 if (verboseLevel>2) {
00304 DumpCutValuesTable();
00305 }
00306 }
00307
00308
00310 void G4VUserPhysicsList::SetDefaultCutValue(G4double value)
00311 {
00312 if (value<0.0) {
00313 #ifdef G4VERBOSE
00314 if (verboseLevel >0){
00315 G4cout << "G4VUserPhysicsList::SetDefaultCutValue: negative cut values"
00316 << " :" << value/mm << "[mm]" << G4endl;
00317 }
00318 #endif
00319 return;
00320 }
00321
00322 defaultCutValue = value;
00323 isSetDefaultCutValue = true;
00324
00325
00326 SetCutValue(defaultCutValue, "gamma");
00327 SetCutValue(defaultCutValue, "e-");
00328 SetCutValue(defaultCutValue, "e+");
00329 SetCutValue(defaultCutValue, "proton");
00330
00331 #ifdef G4VERBOSE
00332 if (verboseLevel >1){
00333 G4cout << "G4VUserPhysicsList::SetDefaultCutValue:"
00334 << "default cut value is changed to :"
00335 << defaultCutValue/mm << "[mm]" << G4endl;
00336 }
00337 #endif
00338 }
00339
00340
00342 G4double G4VUserPhysicsList::GetCutValue(const G4String& name) const
00343 {
00344 size_t nReg = (G4RegionStore::GetInstance())->size();
00345 if (nReg==0) {
00346 #ifdef G4VERBOSE
00347 if (verboseLevel>0){
00348 G4cout << "G4VUserPhysicsList::GetCutValue "
00349 <<" : No Default Region " <<G4endl;
00350 }
00351 #endif
00352 G4Exception("G4VUserPhysicsList::GetCutValue",
00353 "Run0253", FatalException,
00354 "No Default Region");
00355 return -1.*mm;
00356 }
00357 G4Region* region = (*(G4RegionStore::GetInstance()))[0];
00358 return region->GetProductionCuts()->GetProductionCut(name);
00359 }
00360
00362 void G4VUserPhysicsList::SetCutValue(G4double aCut, const G4String& name)
00363 {
00364 SetParticleCuts( aCut ,name );
00365 }
00366
00368 void G4VUserPhysicsList::SetCutValue
00369 (G4double aCut, const G4String& pname, const G4String& rname)
00370 {
00371 G4Region* region = G4RegionStore::GetInstance()->GetRegion(rname);
00372 if (region != 0){
00373
00374 SetParticleCuts( aCut ,pname, region );
00375 } else {
00376 #ifdef G4VERBOSE
00377 if (verboseLevel>0){
00378 G4cout << "G4VUserPhysicsList::SetCutValue "
00379 <<" : No Region of " << rname << G4endl;
00380 }
00381 #endif
00382 }
00383 }
00384
00385
00387 void G4VUserPhysicsList::SetCutsWithDefault()
00388 {
00389 SetDefaultCutValue(defaultCutValue);
00390 G4VUserPhysicsList::SetCuts();
00391 }
00392
00394 void G4VUserPhysicsList::SetCutsForRegion(G4double aCut, const G4String& rname)
00395 {
00396
00397 SetCutValue(aCut, "gamma", rname);
00398 SetCutValue(aCut, "e-", rname);
00399 SetCutValue(aCut, "e+", rname);
00400 SetCutValue(aCut, "proton", rname);
00401 }
00402
00403
00404
00406 void G4VUserPhysicsList::SetParticleCuts( G4double cut, G4ParticleDefinition* particle, G4Region* region)
00407 {
00408 SetParticleCuts(cut, particle->GetParticleName(), region);
00409 }
00410
00412 void G4VUserPhysicsList::SetParticleCuts( G4double cut, const G4String& particleName, G4Region* region)
00413 {
00414 if (cut<0.0) {
00415 #ifdef G4VERBOSE
00416 if (verboseLevel >0){
00417 G4cout << "G4VUserPhysicsList::SetParticleCuts: negative cut values"
00418 << " :" << cut/mm << "[mm]"
00419 << " for "<< particleName << G4endl;
00420 }
00421 #endif
00422 return;
00423 }
00424
00425 if(!region){
00426 size_t nReg = (G4RegionStore::GetInstance())->size();
00427 if (nReg==0) {
00428 #ifdef G4VERBOSE
00429 if (verboseLevel>0){
00430 G4cout << "G4VUserPhysicsList::SetParticleCuts "
00431 <<" : No Default Region " <<G4endl;
00432 }
00433 #endif
00434 G4Exception("G4VUserPhysicsList::SetParticleCuts ",
00435 "Run0254", FatalException,
00436 "No Default Region");
00437 return;
00438 }
00439 region = (*(G4RegionStore::GetInstance()))[0];
00440 }
00441
00442 if ( !isSetDefaultCutValue ){
00443 SetDefaultCutValue(defaultCutValue);
00444 }
00445
00446 G4ProductionCuts* pcuts = region->GetProductionCuts();
00447 pcuts->SetProductionCut(cut,particleName);
00448 #ifdef G4VERBOSE
00449 if (verboseLevel>2){
00450 G4cout << "G4VUserPhysicsList::SetParticleCuts: "
00451 << " :" << cut/mm << "[mm]"
00452 << " for "<< particleName << G4endl;
00453 }
00454 #endif
00455 }
00456
00458 void G4VUserPhysicsList::BuildPhysicsTable()
00459 {
00460
00461 theParticleIterator->reset();
00462 while( (*theParticleIterator)() ){
00463 G4ParticleDefinition* particle = theParticleIterator->value();
00464 PreparePhysicsTable(particle);
00465 }
00466
00467
00468 if (fRetrievePhysicsTable) {
00469 fIsRestoredCutValues = fCutsTable->RetrieveCutsTable(directoryPhysicsTable, fStoredInAscii);
00470
00471 if (!fIsRestoredCutValues) {
00472 #ifdef G4VERBOSE
00473 if (verboseLevel>0){
00474 G4cout << "G4VUserPhysicsList::BuildPhysicsTable"
00475 << " Retrieve Cut Table failed !!" << G4endl;
00476 }
00477 #endif
00478 G4Exception("G4VUserPhysicsList::BuildPhysicsTable",
00479 "Run0255", RunMustBeAborted,
00480 "Fail to retrieve Production Cut Table");
00481 } else {
00482 #ifdef G4VERBOSE
00483 if (verboseLevel>2){
00484 G4cout << "G4VUserPhysicsList::BuildPhysicsTable"
00485 << " Retrieve Cut Table successfully " << G4endl;
00486 }
00487 #endif
00488 }
00489 } else {
00490 #ifdef G4VERBOSE
00491 if (verboseLevel>2){
00492 G4cout << "G4VUserPhysicsList::BuildPhysicsTable"
00493 << " does not retrieve Cut Table but calculate " << G4endl;
00494 }
00495 #endif
00496 }
00497
00498
00499
00500 G4String particleName;
00501 G4ParticleDefinition* GammaP = theParticleTable->FindParticle("gamma");
00502 if(GammaP) BuildPhysicsTable(GammaP);
00503 G4ParticleDefinition* EMinusP = theParticleTable->FindParticle("e-");
00504 if(EMinusP) BuildPhysicsTable(EMinusP);
00505 G4ParticleDefinition* EPlusP = theParticleTable->FindParticle("e+");
00506 if(EPlusP) BuildPhysicsTable(EPlusP);
00507 G4ParticleDefinition* ProtonP = theParticleTable->FindParticle("proton");
00508 if(ProtonP) BuildPhysicsTable(ProtonP);
00509
00510
00511 theParticleIterator->reset();
00512 while( (*theParticleIterator)() ){
00513 G4ParticleDefinition* particle = theParticleIterator->value();
00514 if( particle!=GammaP &&
00515 particle!=EMinusP &&
00516 particle!=EPlusP &&
00517 particle!=ProtonP ){
00518 BuildPhysicsTable(particle);
00519 }
00520 }
00521
00522
00523 fIsPhysicsTableBuilt = true;
00524
00525 }
00527 void G4VUserPhysicsList::BuildPhysicsTable(G4ParticleDefinition* particle)
00528 {
00529 if (fRetrievePhysicsTable) {
00530 if ( !fIsRestoredCutValues){
00531
00532 #ifdef G4VERBOSE
00533 if (verboseLevel>0){
00534 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
00535 << "Physics table can not be retreived and will be calculated "
00536 << G4endl;
00537 }
00538 #endif
00539 fRetrievePhysicsTable = false;
00540
00541 } else {
00542 #ifdef G4VERBOSE
00543 if (verboseLevel>2){
00544 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
00545 << " Retrieve Physics Table for "
00546 << particle->GetParticleName() << G4endl;
00547 }
00548 #endif
00549
00550 RetrievePhysicsTable(particle, directoryPhysicsTable, fStoredInAscii);
00551 }
00552 }
00553
00554 #ifdef G4VERBOSE
00555 if (verboseLevel>2){
00556 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
00557 << "Calculate Physics Table for "
00558 << particle->GetParticleName() << G4endl;
00559 }
00560 #endif
00561
00562
00563 if(!particle->IsShortLived()) {
00564 G4ProcessManager* pManager = particle->GetProcessManager();
00565 if (!pManager) {
00566 #ifdef G4VERBOSE
00567 if (verboseLevel>0){
00568 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
00569 <<" : No Process Manager for "
00570 << particle->GetParticleName() << G4endl;
00571 G4cout << particle->GetParticleName()
00572 << " should be created in your PhysicsList" <<G4endl;
00573 }
00574 #endif
00575 G4Exception("G4VUserPhysicsList::BuildPhysicsTable",
00576 "Run0271", FatalException,
00577 "No process manager");
00578 return;
00579 }
00580 G4ProcessVector* pVector = pManager->GetProcessList();
00581 if (!pVector) {
00582 #ifdef G4VERBOSE
00583 if (verboseLevel>0){
00584 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
00585 <<" : No Process Vector for "
00586 << particle->GetParticleName() <<G4endl;
00587 }
00588 #endif
00589 G4Exception("G4VUserPhysicsList::BuildPhysicsTable",
00590 "Run0272", FatalException,
00591 "No process Vector");
00592 return;
00593 }
00594 for (G4int j=0; j < pVector->size(); ++j) {
00595 (*pVector)[j]->BuildPhysicsTable(*particle);
00596 }
00597 }
00598 }
00599
00601 void G4VUserPhysicsList::PreparePhysicsTable(G4ParticleDefinition* particle)
00602 {
00603
00604
00605 if(!particle->IsShortLived()) {
00606 G4ProcessManager* pManager = particle->GetProcessManager();
00607 if (!pManager) {
00608 #ifdef G4VERBOSE
00609 if (verboseLevel>0) {
00610 G4cout<< "G4VUserPhysicsList::PreparePhysicsTable "
00611 << ": No Process Manager for "
00612 << particle->GetParticleName() <<G4endl;
00613 G4cout << particle->GetParticleName()
00614 << " should be created in your PhysicsList" <<G4endl;
00615 }
00616 #endif
00617 G4Exception("G4VUserPhysicsList::PreparePhysicsTable",
00618 "Run0273", FatalException,
00619 "No process manager");
00620 return;
00621 }
00622
00623 G4ProcessVector* pVector = pManager->GetProcessList();
00624 if (!pVector) {
00625 #ifdef G4VERBOSE
00626 if (verboseLevel>0) {
00627 G4cout << "G4VUserPhysicsList::PreparePhysicsTable "
00628 << ": No Process Vector for "
00629 << particle->GetParticleName() <<G4endl;
00630 }
00631 #endif
00632 G4Exception("G4VUserPhysicsList::PreparePhysicsTable",
00633 "Run0274", FatalException,
00634 "No process Vector");
00635 return;
00636 }
00637 for (G4int j=0; j < pVector->size(); ++j) {
00638 (*pVector)[j]->PreparePhysicsTable(*particle);
00639 }
00640 }
00641 }
00642
00644 void G4VUserPhysicsList::BuildIntegralPhysicsTable(G4VProcess* process,
00645 G4ParticleDefinition* particle)
00646 {
00647
00648
00649
00650
00651
00652 if ( (process->GetProcessName() == "Imsc") ||
00653 (process->GetProcessName() == "IeIoni") ||
00654 (process->GetProcessName() == "IeBrems") ||
00655 (process->GetProcessName() == "Iannihil") ||
00656 (process->GetProcessName() == "IhIoni") ||
00657 (process->GetProcessName() == "IMuIoni") ||
00658 (process->GetProcessName() == "IMuBrems") ||
00659 (process->GetProcessName() == "IMuPairProd") ) {
00660 #ifdef G4VERBOSE
00661 if (verboseLevel>2){
00662 G4cout << "G4VUserPhysicsList::BuildIntegralPhysicsTable "
00663 << " BuildPhysicsTable is invoked for "
00664 << process->GetProcessName()
00665 << "(" << particle->GetParticleName() << ")" << G4endl;
00666 }
00667 #endif
00668 process->BuildPhysicsTable(*particle);
00669 }
00670 }
00671
00673 void G4VUserPhysicsList::DumpList() const
00674 {
00675 theParticleIterator->reset();
00676 G4int idx = 0;
00677 while( (*theParticleIterator)() ){
00678 G4ParticleDefinition* particle = theParticleIterator->value();
00679 G4cout << particle->GetParticleName();
00680 if ((idx++ % 4) == 3) {
00681 G4cout << G4endl;
00682 } else {
00683 G4cout << ", ";
00684 }
00685 }
00686 G4cout << G4endl;
00687 }
00688
00689
00691 void G4VUserPhysicsList::DumpCutValuesTable(G4int flag)
00692 {
00693 fDisplayThreshold = flag;
00694 }
00695
00697 void G4VUserPhysicsList::DumpCutValuesTableIfRequested()
00698 {
00699 if(fDisplayThreshold==0) return;
00700 G4ProductionCutsTable::GetProductionCutsTable()->DumpCouples();
00701 fDisplayThreshold = 0;
00702 }
00703
00704
00706 G4bool G4VUserPhysicsList::StorePhysicsTable(const G4String& directory)
00707 {
00708 G4bool ascii = fStoredInAscii;
00709 G4String dir = directory;
00710 if (dir.isNull()) dir = directoryPhysicsTable;
00711 else directoryPhysicsTable = dir;
00712
00713
00714 if (!fCutsTable->StoreCutsTable(dir, ascii)) {
00715 G4Exception("G4VUserPhysicsList::StorePhysicsTable",
00716 "Run0281", JustWarning,
00717 "Fail to store Cut Table");
00718 return false;
00719 }
00720 #ifdef G4VERBOSE
00721 if (verboseLevel>2){
00722 G4cout << "G4VUserPhysicsList::StorePhysicsTable "
00723 << " Store material and cut values successfully" << G4endl;
00724 }
00725 #endif
00726
00727 G4bool success= true;
00728
00729
00730 theParticleIterator->reset();
00731 while( (*theParticleIterator)() ){
00732 G4ParticleDefinition* particle = theParticleIterator->value();
00733
00734 G4ProcessVector* pVector = (particle->GetProcessManager())->GetProcessList();
00735 G4int j;
00736 for ( j=0; j < pVector->size(); ++j) {
00737 if (!(*pVector)[j]->StorePhysicsTable(particle,dir,ascii)){
00738 G4String comment = "Fail to store physics table for ";
00739 comment += (*pVector)[j]->GetProcessName();
00740 comment += "(" + particle->GetParticleName() + ")";
00741 G4Exception("G4VUserPhysicsList::StorePhysicsTable",
00742 "Run0282", JustWarning,
00743 comment);
00744 success = false;
00745 }
00746 }
00747
00748 }
00749
00750 return success;
00751 }
00752
00753
00754
00756 void G4VUserPhysicsList::SetPhysicsTableRetrieved(const G4String& directory)
00757 {
00758 fRetrievePhysicsTable = true;
00759 if(!directory.isNull()) {
00760 directoryPhysicsTable = directory;
00761 }
00762 fIsCheckedForRetrievePhysicsTable=false;
00763 fIsRestoredCutValues = false;
00764 }
00765
00767 void G4VUserPhysicsList::RetrievePhysicsTable(G4ParticleDefinition* particle,
00768 const G4String& directory,
00769 G4bool ascii)
00770 {
00771 G4int j;
00772 G4bool success[100];
00773
00774 G4ProcessVector* pVector = (particle->GetProcessManager())->GetProcessList();
00775 for ( j=0; j < pVector->size(); ++j) {
00776 success[j] =
00777 (*pVector)[j]->RetrievePhysicsTable(particle,directory,ascii);
00778
00779 if (!success[j]) {
00780 #ifdef G4VERBOSE
00781 if (verboseLevel>2){
00782 G4cout << "G4VUserPhysicsList::RetrievePhysicsTable "
00783 << " Fail to retrieve Physics Table for "
00784 << (*pVector)[j]->GetProcessName() << G4endl;
00785 G4cout << "Calculate Physics Table for "
00786 << particle->GetParticleName() << G4endl;
00787 }
00788 #endif
00789 (*pVector)[j]->BuildPhysicsTable(*particle);
00790 }
00791 }
00792 for ( j=0; j < pVector->size(); ++j) {
00793
00794 if (!success[j]) BuildIntegralPhysicsTable((*pVector)[j], particle);
00795 }
00796 }
00797
00798
00800 void G4VUserPhysicsList::SetApplyCuts(G4bool value, const G4String& name)
00801 {
00802 #ifdef G4VERBOSE
00803 if (verboseLevel>2){
00804 G4cout << "G4VUserPhysicsList::SetApplyCuts for " << name << G4endl;
00805 }
00806 #endif
00807 if(name=="all") {
00808 theParticleTable->FindParticle("gamma")->SetApplyCutsFlag(value);
00809 theParticleTable->FindParticle("e-")->SetApplyCutsFlag(value);
00810 theParticleTable->FindParticle("e+")->SetApplyCutsFlag(value);
00811 theParticleTable->FindParticle("proton")->SetApplyCutsFlag(value);
00812 } else {
00813 theParticleTable->FindParticle(name)->SetApplyCutsFlag(value);
00814 }
00815 }
00816
00818 G4bool G4VUserPhysicsList::GetApplyCuts(const G4String& name) const
00819 {
00820 return theParticleTable->FindParticle(name)->GetApplyCutsFlag();
00821 }
00822
00823
00825 void G4VUserPhysicsList::CheckParticleList()
00826 {
00827 if (! fDisableCheckParticleList ){
00828 thePLHelper->CheckParticleList();
00829 }
00830 }
00831
00833 void G4VUserPhysicsList::AddTransportation()
00834 {
00835 thePLHelper->AddTransportation();
00836 }
00837
00839 void G4VUserPhysicsList::UseCoupledTransportation(G4bool vl)
00840 {
00841 thePLHelper->UseCoupledTransportation(vl);
00842 }
00843
00845 G4bool G4VUserPhysicsList::RegisterProcess(G4VProcess* process,
00846 G4ParticleDefinition* particle)
00847 {
00848 return thePLHelper->RegisterProcess(process, particle);
00849 }
00850
00852 void G4VUserPhysicsList::SetVerboseLevel(G4int value)
00853 {
00854 verboseLevel = value;
00855
00856 fCutsTable->SetVerboseLevel(verboseLevel);
00857
00858 thePLHelper->SetVerboseLevel(verboseLevel);
00859
00860 #ifdef G4VERBOSE
00861 if (verboseLevel >1){
00862 G4cout << "G4VUserPhysicsList::SetVerboseLevel :"
00863 << " Verbose level is set to " << verboseLevel << G4endl;
00864 }
00865 #endif
00866 }
00867
00868
00871
00873 void G4VUserPhysicsList::ResetCuts()
00874 {
00875 #ifdef G4VERBOSE
00876 if (verboseLevel>0){
00877 G4cout << "G4VUserPhysicsList::ResetCuts() is obsolete."
00878 << " This method gives no effect and you can remove it. "<< G4endl;
00879 }
00880 #endif
00881 }