TGeant4/TG4G3PhysicsManager.cxx

   1 // $Id$
   2 // Category: physics
   3 //
   4 // See the class description in the header file.
   5
   6 #include "TG4G3PhysicsManager.h"
   7 #include "TG4G3CutVector.h"
   8 #include "TG4G3ControlVector.h"
   9 #include "TG4G3Defaults.h"
  10
  11 #include <G4ParticleDefinition.hh>
  12 #include <G4VProcess.hh>
  13 #include <G4UImessenger.hh>
  14 #include <G4ProcessTable.hh>
  15
  16 TG4G3PhysicsManager* TG4G3PhysicsManager::fgInstance = 0;
  17
  18 //_____________________________________________________________________________
  19 TG4G3PhysicsManager::TG4G3PhysicsManager()
  20   : fLock(false),
  21     fCutVector(0),
  22     fControlVector(0)
  23 {
  24 //
  25   if (fgInstance) {
  26     TG4Globals::Exception(
  27       "TG4G3PhysicsManager: attempt to create two instances of singleton.");
  28   }
  29
  30   fgInstance = this;
  31
  32   // initialize fIsCutVector
  33   fIsCutVector = new TG4boolVector;
  34   G4int i;
  35   //for (i=0; i<kNofParticlesWSP; i++) fIsCutVector->insert(false);
  36   for (i=0; i<kNofParticlesWSP; i++) fIsCutVector->push_back(false);
  37
  38   // initialize fIsControlVector
  39   fIsControlVector = new TG4boolVector;
  40   //for (i=0; i<kNofParticlesWSP; i++) fIsControlVector->insert(false);
  41   for (i=0; i<kNofParticlesWSP; i++) fIsControlVector->push_back(false);
  42
  43   // define fCutNameVector, fControlNameVector
  44   FillG3CutNameVector();
  45   FillG3ControlNameVector();
  46 }
  47
  48 //_____________________________________________________________________________
  49 TG4G3PhysicsManager::TG4G3PhysicsManager(const TG4G3PhysicsManager& right) {
  50 //
  51   TG4Globals::Exception(
  52     "Attempt to copy TG4G3PhysicsManager singleton.");
  53 }
  54
  55 //_____________________________________________________________________________
  56 TG4G3PhysicsManager::~TG4G3PhysicsManager() {
  57 //
  58   delete fIsCutVector;
  59   delete fIsControlVector;
  60 }
  61
  62 // operators
  63
  64 //_____________________________________________________________________________
  65 TG4G3PhysicsManager&
  66 TG4G3PhysicsManager::operator=(const TG4G3PhysicsManager& right)
  67 {
  68   // check assignement to self
  69   if (this == &right) return *this;
  70
  71   TG4Globals::Exception(
  72     "Attempt to assign TG4G3PhysicsManager singleton.");
  73
  74   return *this;
  75 }
  76
  77 // private methods
  78
  79 //_____________________________________________________________________________
  80 void TG4G3PhysicsManager::FillG3CutNameVector()
  81 {
  82 // Defines fCutNameVector.
  83 // ---
  84
  85   fG3CutNameVector.insert("CUTGAM");
  86   fG3CutNameVector.insert("CUTELE");
  87   fG3CutNameVector.insert("CUTNEU");
  88   fG3CutNameVector.insert("CUTHAD");
  89   fG3CutNameVector.insert("CUTMUO");
  90   fG3CutNameVector.insert("BCUTE");
  91   fG3CutNameVector.insert("BCUTM");
  92   fG3CutNameVector.insert("DCUTE");
  93   fG3CutNameVector.insert("DCUTM");
  94   fG3CutNameVector.insert("PPCUTM");
  95 }
  96
  97 //_____________________________________________________________________________
  98 void TG4G3PhysicsManager::FillG3ControlNameVector()
  99 {
 100 // Defines fControlNameVector.
 101 // ---
 102
 103   fG3ControlNameVector.insert("PAIR");
 104   fG3ControlNameVector.insert("COMP");
 105   fG3ControlNameVector.insert("PHOT");
 106   fG3ControlNameVector.insert("PFIS");
 107   fG3ControlNameVector.insert("DRAY");
 108   fG3ControlNameVector.insert("ANNI");
 109   fG3ControlNameVector.insert("BREM");
 110   fG3ControlNameVector.insert("HADR");
 111   fG3ControlNameVector.insert("MUNU");
 112   fG3ControlNameVector.insert("DCAY");
 113   fG3ControlNameVector.insert("LOSS");
 114   fG3ControlNameVector.insert("MULS");
 115 }
 116
 117 //_____________________________________________________________________________
 118 void TG4G3PhysicsManager::SetCut(TG4G3Cut cut, G4double cutValue)
 119 {
 120 // Sets kinetic energy cut (in a G3-like way).
 121 // ---
 122
 123   if (!fCutVector) {
 124     // create vector of kinetic energy cut values
 125     fCutVector = new TG4G3CutVector();
 126   }
 127   fCutVector->SetG3Cut(cut, cutValue);
 128   SwitchIsCutVector(cut);
 129 }
 130
 131 //_____________________________________________________________________________
 132 void TG4G3PhysicsManager::SetProcess(TG4G3Control control, G4int controlValue)
 133 {
 134 // Sets control process control (in a G3-like way).
 135 // ---
 136
 137   if (!fControlVector) {
 138     // create vector of control process control values
 139     fControlVector = new TG4G3ControlVector;
 140   }
 141   fControlVector->SetG3Control(control, controlValue);
 142 }
 143
 144 //_____________________________________________________________________________
 145 void TG4G3PhysicsManager::SwitchIsCutVector(TG4G3Cut cut)
 146 {
 147 // Updates the vector of booleans (fIsCutVector) for the specified cut.
 148 // ---
 149
 150   switch (cut) {
 151     case kCUTGAM:
 152            (*fIsCutVector)[kGamma] = true;
 153            break;
 154     case kBCUTE:
 155            (*fIsCutVector)[kGamma] = true;
 156            break;
 157     case kBCUTM:
 158            (*fIsCutVector)[kGamma] = true;
 159            break;
 160     case kCUTELE:
 161            (*fIsCutVector)[kElectron] = true;
 162            break;
 163     case kDCUTE:
 164            (*fIsCutVector)[kElectron] = true;
 165            break;
 166     case kDCUTM:
 167            (*fIsCutVector)[kElectron] = true;
 168            break;
 169     case kCUTNEU:
 170            (*fIsCutVector)[kNeutralHadron] = true;
 171            break;
 172     case kCUTHAD:
 173            (*fIsCutVector)[kChargedHadron] = true;
 174            break;
 175     case kCUTMUO:
 176            (*fIsCutVector)[kMuon] = true;
 177            break;
 178     default:
 179            break;
 180   }
 181 }
 182
 183 //_____________________________________________________________________________
 184 void TG4G3PhysicsManager::SwitchIsControlVector(TG4G3Control control)
 185 {
 186 // Updates the vector of booleans (fIsControlVector) for the specified control.
 187 // ---
 188
 189   switch (control) {
 190     case kPAIR:
 191            // gamma
 192            (*fIsControlVector)[kGamma] = true;
 193            break;
 194     case kCOMP:
 195            // gamma
 196            (*fIsControlVector)[kGamma] = true;
 197            break;
 198     case kPHOT:
 199            // gamma
 200            (*fIsControlVector)[kGamma] = true;
 201            break;
 202     case kPFIS:
 203            // gamma
 204            (*fIsControlVector)[kGamma] = true;
 205            break;
 206     case kDRAY:
 207            // all charged particles
 208            (*fIsControlVector)[kElectron] = true;
 209            (*fIsControlVector)[kEplus] = true;
 210            (*fIsControlVector)[kChargedHadron] = true;
 211            (*fIsControlVector)[kMuon] = true;
 212            break;
 213     case kANNI:
 214            // e+ only
 215            (*fIsControlVector)[kEplus] = true;
 216            break;
 217     case kBREM:
 218            // e-/e+, muons
 219            (*fIsControlVector)[kElectron] = true;
 220            (*fIsControlVector)[kEplus] = true;
 221            (*fIsControlVector)[kMuon] = true;
 222            break;
 223     case kHADR:
 224            // hadrons
 225            (*fIsControlVector)[kNeutralHadron] = true;
 226            (*fIsControlVector)[kChargedHadron] = true;
 227            break;
 228     case kMUNU:
 229            // muons
 230            (*fIsControlVector)[kMuon] = true;
 231            break;
 232     case kDCAY:
 233            // any
 234            (*fIsControlVector)[kAny] = true;
 235            break;
 236     case kLOSS:
 237            // all charged particles
 238            (*fIsControlVector)[kElectron] = true;
 239            (*fIsControlVector)[kEplus] = true;
 240            (*fIsControlVector)[kChargedHadron] = true;
 241            (*fIsControlVector)[kMuon] = true;
 242            break;
 243     case kMULS:
 244            // all charged particles
 245            (*fIsControlVector)[kElectron] = true;
 246            (*fIsControlVector)[kEplus] = true;
 247            (*fIsControlVector)[kChargedHadron] = true;
 248            (*fIsControlVector)[kMuon] = true;
 249            break;
 250     default:
 251           break;
 252   }
 253 }
 254
 255 //_____________________________________________________________________________
 256 TG4G3Cut TG4G3PhysicsManager::GetG3Cut(G4String cutName)
 257 {
 258 // Retrieves corresponding TG4G3Cut constant from the cutName.
 259 // ---
 260
 261   if      (cutName == fG3CutNameVector[kCUTGAM]) return kCUTGAM;
 262   else if (cutName == fG3CutNameVector[kBCUTE])  return kBCUTE;
 263   else if (cutName == fG3CutNameVector[kBCUTM])  return kBCUTM;
 264   else if (cutName == fG3CutNameVector[kCUTELE]) return kCUTELE;
 265   else if (cutName == fG3CutNameVector[kDCUTE])  return kDCUTE;
 266   else if (cutName == fG3CutNameVector[kDCUTM])  return kDCUTM;
 267   else if (cutName == fG3CutNameVector[kCUTNEU]) return kCUTNEU;
 268   else if (cutName == fG3CutNameVector[kCUTHAD]) return kCUTHAD;
 269   else if (cutName == fG3CutNameVector[kCUTMUO]) return kCUTMUO;
 270   else return kNoG3Cuts;
 271 }
 272
 273 //_____________________________________________________________________________
 274 TG4G3Control TG4G3PhysicsManager::GetG3Control(G4String controlName)
 275 {
 276 // Retrieves corresponding TG4G3Control constant from the controlName.
 277 // ---
 278
 279   if      (controlName == fG3ControlNameVector[kPAIR]) return kPAIR;
 280   else if (controlName == fG3ControlNameVector[kCOMP]) return kCOMP;
 281   else if (controlName == fG3ControlNameVector[kPHOT]) return kPHOT;
 282   else if (controlName == fG3ControlNameVector[kPFIS]) return kPFIS;
 283   else if (controlName == fG3ControlNameVector[kDRAY]) return kDRAY;
 284   else if (controlName == fG3ControlNameVector[kANNI]) return kANNI;
 285   else if (controlName == fG3ControlNameVector[kBREM]) return kBREM;
 286   else if (controlName == fG3ControlNameVector[kHADR]) return kHADR;
 287   else if (controlName == fG3ControlNameVector[kMUNU]) return kMUNU;
 288   else if (controlName == fG3ControlNameVector[kDCAY]) return kDCAY;
 289   else if (controlName == fG3ControlNameVector[kLOSS]) return kLOSS;
 290   else if (controlName == fG3ControlNameVector[kMULS]) return kMULS;
 291   else return kNoG3Controls;
 292 }
 293
 294 // public methods
 295
 296 //_____________________________________________________________________________
 297 void TG4G3PhysicsManager::CheckLock()
 298 {
 299 // Gives exception in case the physics manager is locked.
 300 // Prevents from modifying physics setup after the physics manager is locked.
 301 // ---
 302
 303   if (fLock) {
 304     G4String text = "TG4PhysicsManager: \n";
 305     text = text + "    It is too late to change physics setup. \n";
 306     text = text + "    PhysicsManager has been already locked.";
 307     TG4Globals::Exception(text);
 308   }
 309 }
 310
 311 //_____________________________________________________________________________
 312 G4VProcess* TG4G3PhysicsManager::FindProcess(G4String processName) const
 313 {
 314 // Finds G4VProcess with specified name.
 315 // ---
 316
 317   G4ProcessTable* processTable = G4ProcessTable::GetProcessTable();
 318
 319   G4ProcessVector* processVector
 320     = processTable->FindProcesses(processName);
 321   G4VProcess* firstFoundProcess = 0;
 322   if (processVector->entries()>0) firstFoundProcess= (*processVector)[0];
 323
 324   processVector->clear();
 325   delete processVector;
 326
 327   return firstFoundProcess;
 328 }
 329
 330 //_____________________________________________________________________________
 331 G4bool TG4G3PhysicsManager::CheckCutWithTheVector(G4String name,
 332                                  G4double value, TG4G3Cut& cut)
 333 {
 334 // Retrieves corresponding TG4G3Cut from the name and
 335 // in case the value is different from the value in cutVector
 336 // sets true the value of the fIsCutVector element
 337 // corresponding to this cut and returns true;
 338 // returns false otherwise.
 339 // ---
 340
 341   // convert cut name -> TG4G3Cut
 342   cut = GetG3Cut(name);
 343
 344   // set switch vector element only if the value
 345   // is different from the value in cutVector
 346   if (cut !=kNoG3Cuts) {
 347     // get tolerance from TG4G3Defaults in GeV
 348     G4double tolerance = TG4G3Defaults::CutTolerance()/GeV;
 349     if (!(fCutVector) || (abs(value - (*fCutVector)[cut]) > tolerance)) {
 350       SwitchIsCutVector(cut);
 351       return true;
 352     }
 353     else return false;
 354   }
 355   return false;
 356 }
 357
 358 //_____________________________________________________________________________
 359 G4bool TG4G3PhysicsManager::CheckControlWithTheVector(G4String name,
 360                                  G4double value, TG4G3Control& control)
 361 {
 362 // Retrieves corresponding TG4G3Control from the name and
 363 // in case the value is different from the value in controlVector
 364 // sets true the value of the fIsControlVector element
 365 // corresponding to this control and returns true;
 366 // returns false otherwise.
 367 // ---
 368
 369   // convert control name -> TG4G3Control
 370   control = GetG3Control(name);
 371
 372   // set switch vector element only if the value
 373   // is different from the value in controlVector
 374   if (control !=kNoG3Controls) {
 375     if (!(fControlVector) || (abs(value - (*fControlVector)[control]) > 0.01)) {
 376       SwitchIsControlVector(control);
 377       return true;
 378     }
 379     else return false;
 380   }
 381   return false;
 382 }
 383
 384 //_____________________________________________________________________________
 385 G4bool TG4G3PhysicsManager::CheckCutWithG3Defaults(G4String name,
 386                                  G4double value, TG4G3Cut& cut)
 387 {
 388 // Retrieves corresponding TG4G3Cut from the name and
 389 // in case the value is different from the G3 default value
 390 // sets true the value of the SwitchCutVector element
 391 // corresponding to this cut and returns true;
 392 // returns false otherwise.
 393 // ---
 394
 395   // convert cut name -> TG4G3Cut
 396   cut = GetG3Cut(name);
 397
 398   // set switch vector element only if the value
 399   // is different from G3 default
 400   if (cut !=kNoG3Cuts) {
 401     if (!TG4G3Defaults::IsDefaultCut(cut, value)) {
 402       SwitchIsCutVector(cut);
 403       return true;
 404     }
 405     else return false;
 406   }
 407   return false;
 408 }
 409
 410 //_____________________________________________________________________________
 411 G4bool TG4G3PhysicsManager::CheckControlWithG3Defaults(G4String name,
 412                                  G4double value, TG4G3Control& control)
 413 {
 414 // Retrieves corresponding TG4G3Control from the name and
 415 // in case the value is different from the G3 default value
 416 // sets true the value of the SwitchControlVector element
 417 // corresponding to this control and returns true;
 418 // returns false otherwise.
 419 // ---
 420
 421   // convert control name -> TG4G3Control
 422   control = GetG3Control(name);
 423
 424   // set switch vector element only if the value
 425   // is different from G3 default
 426   if (control !=kNoG3Controls) {
 427     if (!TG4G3Defaults::IsDefaultControl(control, value)) {
 428       SwitchIsControlVector(control);
 429       return true;
 430     }
 431     else return false;
 432   }
 433   return false;
 434 }
 435
 436 //_____________________________________________________________________________
 437 void TG4G3PhysicsManager::SetG3DefaultCuts()
 438 {
 439 // Sets G3 default values of kinetic energy cuts.
 440 // ---
 441
 442   CheckLock();
 443   if (!fCutVector) {
 444     // create vector of kinetic energy cut values
 445     fCutVector = new TG4G3CutVector();
 446   }
 447   fCutVector->SetG3Defaults();
 448 }
 449
 450 //_____________________________________________________________________________
 451 void TG4G3PhysicsManager::SetG3DefaultControls()
 452 {
 453 // Sets G3 default values of control process controls.
 454 // ---
 455
 456   CheckLock();
 457   if (!fControlVector) {
 458     // create vector of control process control values
 459     fControlVector = new TG4G3ControlVector;
 460   }
 461   fControlVector->SetG3Defaults();
 462 }
 463
 464 //_____________________________________________________________________________
 465 G4bool TG4G3PhysicsManager::IsSpecialCuts() const
 466 {
 467 // Returns true if any special cut value is set.
 468 // ---
 469
 470   for (G4int i=0; i<kNofParticlesWSP; i++)
 471   {  if ((*fIsCutVector)[i]) return true; }
 472
 473   return false;
 474 }
 475
 476 //_____________________________________________________________________________
 477 G4bool TG4G3PhysicsManager::IsSpecialControls() const
 478 {
 479 // Returns true if any special control value is set.
 480 // ---
 481
 482   for (G4int i=0; i<kNofParticlesWSP; i++)
 483   {  if ((*fIsControlVector)[i]) return true; }
 484
 485   return false;
 486 }
 487
 488 //_____________________________________________________________________________
 489 TG4G3ParticleWSP TG4G3PhysicsManager::GetG3ParticleWSP(
 490                                       G4ParticleDefinition* particle) const
 491 {
 492 // Returns TG4G3ParticleWSP constant for the specified particle.
 493 // (See TG4G3ParticleWSP.h, too.)
 494 // ---
 495
 496   G4String name = particle->GetParticleName();
 497   G4String pType = particle->GetParticleType();
 498
 499   if (name == "gamma") {
 500     return kGamma;
 501   }
 502   else if (name == "e-") {
 503     return kElectron;
 504   }
 505   else if (name == "e+") {
 506     return kEplus;
 507   }
 508   else if (( pType == "baryon" || pType == "meson" || pType == "nucleus" )) {
 509     if (particle->GetPDGCharge() == 0) {
 510       return kNeutralHadron;
 511     }
 512     else
 513       return kChargedHadron;
 514   }
 515   else if ( name == "mu-" || name == "mu+" ) {
 516     return kMuon;
 517   }
 518   else {
 519     return kNofParticlesWSP;
 520   }
 521 }
 522
 523 //_____________________________________________________________________________
 524 void TG4G3PhysicsManager::GetG3ParticleWSPName(G4int particleWSP,
 525                                                G4String& name) const
 526 {
 527 // Fills the passed name with the name/type of particle specified
 528 // by TG4G3ParticleWSP constant.
 529 // (See TG4G3ParticleWSP.h, too.)
 530 // ---
 531
 532   switch (particleWSP) {
 533     case kGamma:
 534       name = "Gamma";
 535       break;
 536     case kElectron:
 537       name = "Electron";
 538       break;
 539     case kEplus:
 540       name = "Eplus";
 541       break;
 542     case kNeutralHadron:
 543       name = "NeutralHadron";
 544       break;
 545     case kChargedHadron:
 546       name = "ChargedHadron";
 547       break;
 548     case kMuon:
 549       name = "Muon";
 550       break;
 551     case kAny:
 552       name = "Any";
 553       break;
 554     case kNofParticlesWSP:
 555       name = "NoSP";
 556       break;
 557     default:
 558       G4String text = "TG4G3PhysicsManager::GetG3ParticleWSPName:\n";
 559       text = text + "   Wrong particleWSP.";
 560       TG4Globals::Exception(text);
 561       break;
 562   }
 563 }
 564