TGeant4/TG4PhysicsManager.cxx

   1 // $Id$
   2 // Category: physics
   3 //
   4 // See the class description in the header file.
   5
   6 #include "TG4PhysicsManager.h"
   7 #include "TG4PhysicsList.h"
   8 #include "TG4CutVector.h"
   9 #include "TG4FlagVector.h"
  10 #include "TG4G3Defaults.h"
  11
  12 #include <G4ParticleDefinition.hh>
  13 #include <G4VProcess.hh>
  14
  15 #include <TDatabasePDG.h>
  16
  17 TG4PhysicsManager* TG4PhysicsManager::fgInstance = 0;
  18
  19 TG4PhysicsManager::TG4PhysicsManager()
  20   : fLock(false),
  21     fPhysicsList(0),
  22     fCutVector(0),
  23     fFlagVector(0)
  24 {
  25 //
  26   if (fgInstance) {
  27     TG4Globals::Exception(
  28       "TG4PhysicsManager: attempt to create two instances of singleton.");
  29   }
  30
  31   fgInstance = this;
  32
  33   // initialize fIsCutVector
  34   fIsCutVector = new TG4boolVector;
  35   G4int i;
  36   //for (i=0; i<kNofParticlesWSP; i++) fIsCutVector->insert(false);
  37   for (i=0; i<kNofParticlesWSP; i++) fIsCutVector->push_back(false);
  38
  39   // initialize fIsFlagVector
  40   fIsFlagVector = new TG4boolVector;
  41   //for (i=0; i<kNofParticlesWSP; i++) fIsFlagVector->insert(false);
  42   for (i=0; i<kNofParticlesWSP; i++) fIsFlagVector->push_back(false);
  43
  44   // define fCutNameVector, fFlagNameVector
  45   FillG3CutNameVector();
  46   FillG3FlagNameVector();
  47
  48   // fill process name map
  49   FillProcessMap();
  50 }
  51
  52 TG4PhysicsManager::TG4PhysicsManager(const TG4PhysicsManager& right) {
  53 //
  54   TG4Globals::Exception(
  55     "Attempt to copy TG4PhysicsManager singleton.");
  56 }
  57
  58 TG4PhysicsManager::~TG4PhysicsManager() {
  59 //
  60   delete fIsCutVector;
  61   delete fIsFlagVector;
  62 }
  63
  64 // operators
  65
  66 TG4PhysicsManager&
  67 TG4PhysicsManager::operator=(const TG4PhysicsManager& right)
  68 {
  69   // check assignement to self
  70   if (this == &right) return *this;
  71
  72   TG4Globals::Exception(
  73     "Attempt to assign TG4PhysicsManager singleton.");
  74
  75   return *this;
  76 }
  77
  78 // private methods
  79
  80 void TG4PhysicsManager::LockException() const
  81 {
  82 // Gives exception in case of attempt to modified physics
  83 // setup after physics manager was locked.
  84 // ---
  85
  86   G4String text = "TG4PhysicsManager: \n";
  87   text = text + "    It is too late to change physics setup. \n";
  88   text = text + "    PhysicsManager has been already locked.";
  89   TG4Globals::Exception(text);
  90 }
  91
  92 void TG4PhysicsManager::FillG3CutNameVector()
  93 {
  94 // Defines fCutNameVector.
  95 // ---
  96
  97   fG3CutNameVector.insert("CUTGAM");
  98   fG3CutNameVector.insert("CUTELE");
  99   fG3CutNameVector.insert("CUTNEU");
 100   fG3CutNameVector.insert("CUTHAD");
 101   fG3CutNameVector.insert("CUTMUO");
 102   fG3CutNameVector.insert("BCUTE");
 103   fG3CutNameVector.insert("BCUTM");
 104   fG3CutNameVector.insert("DCUTE");
 105   fG3CutNameVector.insert("DCUTM");
 106   fG3CutNameVector.insert("PPCUTM");
 107 }
 108
 109 void TG4PhysicsManager::FillG3FlagNameVector()
 110 {
 111 // Defines fFlagNameVector.
 112 // ---
 113
 114   fG3FlagNameVector.insert("PAIR");
 115   fG3FlagNameVector.insert("COMP");
 116   fG3FlagNameVector.insert("PHOT");
 117   fG3FlagNameVector.insert("PFIS");
 118   fG3FlagNameVector.insert("DRAY");
 119   fG3FlagNameVector.insert("ANNI");
 120   fG3FlagNameVector.insert("BREM");
 121   fG3FlagNameVector.insert("HADR");
 122   fG3FlagNameVector.insert("MUNU");
 123   fG3FlagNameVector.insert("DCAY");
 124   fG3FlagNameVector.insert("LOSS");
 125   fG3FlagNameVector.insert("MULS");
 126 }
 127
 128 void TG4PhysicsManager::FillProcessMap()
 129 {
 130 // Fills fProcessMap.
 131 // The default G4 process names are used in the map.
 132 // ---
 133
 134   // multiple scattering
 135   fProcessMap.Add("msc",  kPMultipleScattering);
 136   fProcessMap.Add("Imsc", kPMultipleScattering);
 137
 138   // continuous energy loss
 139   // !! including delta rays
 140   fProcessMap.Add("eIoni",  kPEnergyLoss);
 141   fProcessMap.Add("IeIoni", kPEnergyLoss);
 142   fProcessMap.Add("LowEnergyIoni", kPEnergyLoss);
 143   fProcessMap.Add("hIoni",  kPEnergyLoss);
 144   fProcessMap.Add("IhIoni", kPEnergyLoss);
 145   fProcessMap.Add("hLowEIoni", kPEnergyLoss);
 146   fProcessMap.Add("MuIoni", kPEnergyLoss);
 147   fProcessMap.Add("IMuIonisation", kPEnergyLoss);
 148   fProcessMap.Add("ionIoni",  kPEnergyLoss);
 149   fProcessMap.Add("ionLowEIoni",  kPEnergyLoss);
 150   fProcessMap.Add("PAIonisation",  kPEnergyLoss);
 151
 152   // bending in mag. field
 153   // kPMagneticFieldL
 154
 155   // particle decay
 156   fProcessMap.Add("Decay", kPDecay);
 157
 158   // photon pair production or
 159   // muon direct pair production
 160   fProcessMap.Add("conv", kPPair);
 161   fProcessMap.Add("LowEnConversion", kPPair);
 162   fProcessMap.Add("MuPairProd", kPPair);
 163   fProcessMap.Add("IMuPairProduction", kPPair);
 164
 165   // Compton scattering
 166   fProcessMap.Add("compt", kPCompton);
 167   fProcessMap.Add("LowEnCompton", kPCompton);
 168   fProcessMap.Add("polarCompt", kPCompton);
 169
 170   // photoelectric effect
 171   fProcessMap.Add("phot", kPPhotoelectric);
 172   fProcessMap.Add("LowEnPhotoElec", kPPhotoelectric);
 173
 174   // bremsstrahlung
 175   fProcessMap.Add("eBrem", kPBrem);
 176   fProcessMap.Add("IeBrem", kPBrem);
 177   fProcessMap.Add("MuBrem", kPBrem);
 178   fProcessMap.Add("IMuBremsstrahlung", kPBrem);
 179   fProcessMap.Add("LowEnBrem", kPBrem);
 180
 181   // delta-ray production
 182   // kPDeltaRay
 183   // has to be distinguished from kPEnergyLoss on flight
 184
 185   // positron annihilation
 186   fProcessMap.Add("annihil", kPAnnihilation);
 187   fProcessMap.Add("Iannihil", kPAnnihilation);
 188
 189   // hadronic interaction
 190   // kPHadronic
 191
 192   // nuclear evaporation
 193   // kPEvaporation
 194
 195   // nuclear fission
 196   // kPNuclearFission
 197
 198   // nuclear absorption
 199   fProcessMap.Add("PionMinusAbsorptionAtRest", kPNuclearAbsorption);
 200   fProcessMap.Add("PiMinusAbsorptionAtRest", kPNuclearAbsorption);
 201   fProcessMap.Add("KaonMinusAbsorption", kPNuclearAbsorption);
 202   fProcessMap.Add("KaonMinusAbsorptionAtRest", kPNuclearAbsorption);
 203
 204   // antiproton annihilation
 205   fProcessMap.Add("AntiProtonAnnihilationAtRest", kPPbarAnnihilation);
 206   // fProcessMap.Add("AntiNeutronAnnihilationAtRest", not defined);
 207
 208   // neutron capture
 209   fProcessMap.Add("NeutronCaptureAtRest", kPNCapture);
 210   // fProcessMap.Add("LCapture", hadron capture not defined);
 211
 212   // hadronic elastic incoherent scattering
 213   fProcessMap.Add("LElastic", kPHElastic);
 214
 215   // hadronic inelastic scattering
 216   fProcessMap.Add("inelastic", kPHInhelastic);
 217
 218   // muon nuclear interaction
 219   fProcessMap.Add("MuNucl", kPMuonNuclear);
 220
 221   // exceeded time of flight cut
 222   // kPTOFlimit
 223
 224   // nuclear photofission
 225   // kPPhotoFission
 226
 227   // Rayleigh scattering
 228   fProcessMap.Add("Rayleigh Scattering", kPRayleigh);
 229
 230   // no mechanism is active, usually at the entrance of a new volume
 231   fProcessMap.Add("Transportation", kPNull);
 232
 233   // particle has fallen below energy threshold and tracking stops
 234   // kPStop
 235
 236   // Cerenkov photon absorption
 237   fProcessMap.Add("Absorption", kPLightAbsorption);
 238
 239   // Cerenkov photon reflection/refraction
 240   // kPLightScattering, kPLightReflection, kPLightRefraction
 241   // has to be inquired from the G4OpBoundary process
 242
 243   // synchrotron radiation
 244   fProcessMap.Add("SynchrotronRadiation", kPSynchrotron);
 245 }
 246
 247
 248 G4int TG4PhysicsManager::GetPDGEncoding(G4ParticleDefinition* particle)
 249 {
 250 // Returns the PDG code of particle;
 251 // if standard PDG code is not defined the TDatabasePDG
 252 // is used.
 253 // ---
 254
 255   // get PDG encoding from G4 particle definition
 256   G4int pdgEncoding = particle->GetPDGEncoding();
 257
 258   if (pdgEncoding == 0) {
 259     // get PDG encoding from TDatabasePDG
 260
 261     // get particle name from the name map
 262     G4String g4name = particle->GetParticleName();
 263     G4String tname = fParticleNameMap.GetSecond(g4name);
 264     if (tname == "Undefined") {
 265       G4String text = "TG4PhysicsManager::GetPDGEncoding: \n";
 266       text = text + "    Particle " + g4name;
 267       text = text + " was not found in the name map.";
 268       TG4Globals::Exception(text);
 269     }
 270
 271     // get particle from TDatabasePDG
 272     TDatabasePDG* pdgDB = TDatabasePDG::Instance();
 273     TParticlePDG* particle = pdgDB->GetParticle(tname);
 274     if (!particle) {
 275       G4String text = "TG4PhysicsManager::GetPDGEncoding: \n";
 276       text = text + "    Particle " + tname;
 277       text = text + " was not found in TDatabasePDG.";
 278       TG4Globals::Exception(text);
 279     }
 280
 281     // get PDG encoding
 282     pdgEncoding = particle->PdgCode();
 283   }
 284
 285   return pdgEncoding;
 286 }
 287
 288 AliMCProcess TG4PhysicsManager::GetMCProcess(const G4VProcess* process)
 289 {
 290 // Returns the AliMCProcess code of the specified G4 process.
 291 // ---
 292
 293   if (!process) return kPNoProcess;
 294
 295   G4String name = process->GetProcessName();
 296   G4int code = fProcessMap.GetSecond(name);
 297
 298   if (code == 0) return kPNoProcess;
 299
 300   return (AliMCProcess)code;
 301 }
 302
 303 G4int TG4PhysicsManager::GetPDGEncoding(G4String particleName)
 304 {
 305 // Returns the PDG code of particle specified by name.
 306 // ---
 307
 308   G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
 309
 310   G4ParticleDefinition* particle = 0;
 311   particle = particleTable->FindParticle(particleName);
 312   if (!particle) {
 313     G4String text = "TG4PhysicsManager::GetPDGEncoding:\n";
 314     text = text + "   G4ParticleTable::FindParticle() " + particleName;
 315     text = text + " failed.";
 316     TG4Globals::Exception(text);
 317   }
 318
 319   return GetPDGEncoding(particle);
 320 }
 321
 322 void TG4PhysicsManager::SetCut(TG3Cut cut, G4double cutValue)
 323 {
 324 // Sets kinetic energy cut (in a G3-like way).
 325 // ---
 326
 327   if (!fCutVector) {
 328     // create vector of kinetic energy cut values
 329     fCutVector = new TG4CutVector();
 330   }
 331   fCutVector->SetG3Cut(cut, cutValue);
 332   SwitchIsCutVector(cut);
 333 }
 334
 335 void TG4PhysicsManager::SetProcess(TG3Flag flag, G4int flagValue)
 336 {
 337 // Sets control process flag (in a G3-like way).
 338 // ---
 339
 340   if (!fFlagVector) {
 341     // create vector of control process flag values
 342     fFlagVector = new TG4FlagVector;
 343   }
 344   fFlagVector->SetG3Flag(flag, flagValue);
 345 }
 346
 347
 348 Float_t TG4PhysicsManager::Xsec(char* ch, Float_t p1, Int_t i1, Int_t i2)
 349 {
 350 // Not yet implemented -> gives exception.
 351 // ---
 352
 353   TG4Globals::Exception(
 354     "TG4PhysicsManager::Xsec: not yet implemented.");
 355
 356   return 0.;
 357 }
 358
 359 Int_t TG4PhysicsManager::IdFromPDG(Int_t pdgID) const
 360 {
 361 // G4 does not use the integer particle identifiers
 362 // Id <-> PDG is identity.
 363 // ---
 364
 365   return pdgID;
 366 }
 367
 368 Int_t TG4PhysicsManager::PDGFromId(Int_t mcID) const
 369 {
 370 // G4 does not use integer particle identifiers
 371 // Id <-> PDG is identity.
 372 // ---
 373
 374   return mcID;
 375 }
 376
 377 void  TG4PhysicsManager::DefineParticles()
 378 {
 379   // ======
 380   // Taken from TGeant3
 381   //
 382   // Use ENDF-6 mapping for ions = 10000*z+10*a+iso
 383   // and add 1 000 000
 384   // and numbers above 5 000 000 for special applications
 385   //
 386
 387   const Int_t kion=10000000;
 388   const Int_t kspe=50000000;
 389
 390   const Double_t kGeV=0.9314943228;
 391   const Double_t kHslash = 1.0545726663e-27;
 392   const Double_t kErgGeV = 1/1.6021773349e-3;
 393   const Double_t kHshGeV = kHslash*kErgGeV;
 394   const Double_t kYearsToSec = 3600*24*365.25;
 395
 396   TDatabasePDG *pdgDB = TDatabasePDG::Instance();
 397
 398   pdgDB->AddParticle("Deuteron","Deuteron",2*kGeV+8.071e-3,kTRUE,
 399                      0,1,"Ion",kion+10020);
 400
 401   pdgDB->AddParticle("Triton","Triton",3*kGeV+14.931e-3,kFALSE,
 402                      kHshGeV/(12.33*kYearsToSec),1,"Ion",kion+10030);
 403
 404   pdgDB->AddParticle("Alpha","Alpha",4*kGeV+2.424e-3,kTRUE,
 405                      kHshGeV/(12.33*kYearsToSec),2,"Ion",kion+20040);
 406
 407   pdgDB->AddParticle("HE3","HE3",3*kGeV+14.931e-3,kFALSE,
 408                      0,2,"Ion",kion+20030);
 409
 410   pdgDB->AddParticle("Cherenkov","Cherenkov",0,kFALSE,
 411                      0,0,"Special",kspe+50);
 412
 413   pdgDB->AddParticle("FeedbackPhoton","FeedbackPhoton",0,kFALSE,
 414                      0,0,"Special",kspe+51);
 415
 416
 417   // To do: define the PDG database extension
 418   // in a common part.
 419   //
 420   // AliMC::ExtendPDGDatabase();
 421   //
 422   // end of "common" implementation
 423   // ======
 424
 425   // map G4 particle names to TDatabasePDG names
 426   // (the map is built only for particles that have not
 427   //  defined standard PDG encoding)
 428
 429   fParticleNameMap.Add("deuteron","Deuteron");
 430   fParticleNameMap.Add("triton",  "Triton");
 431   fParticleNameMap.Add("alpha",   "Alpha");
 432   fParticleNameMap.Add("He3",     "HE3");
 433   fParticleNameMap.Add("opticalphoton","Cherenkov");
 434   // fParticleNameMap.Add("???","FeedbackPhoton");
 435   fParticleNameMap.Add("geantino", "Rootino");
 436
 437   // map G4 particle names to TDatabasePDG encodings
 438   fParticlePDGMap.Add("deuteron", GetPDGEncoding("deuteron"));
 439   fParticlePDGMap.Add("triton", GetPDGEncoding("triton"));
 440   fParticlePDGMap.Add("alpha", GetPDGEncoding("alpha"));
 441   fParticlePDGMap.Add("He3", GetPDGEncoding("He3") );
 442   fParticlePDGMap.Add("opticalphoton", GetPDGEncoding("opticalphoton"));
 443   // fParticlePDGMap.Add("???","FeedbackPhoton");
 444   fParticleNameMap.Add("geantino", GetPDGEncoding("geantino"));
 445
 446   // add verbose
 447   G4cout << "Particle maps have been filled." << G4endl;
 448   //fParticleNameMap.PrintAll();
 449   //fParticlePDGMap.PrintAll();
 450 }
 451
 452 void TG4PhysicsManager::SwitchIsCutVector(TG3Cut cut)
 453 {
 454 // Updates the vector of booleans (fIsCutVector) for the specified cut.
 455 // ---
 456
 457   switch (cut) {
 458     case kCUTGAM:
 459            (*fIsCutVector)[kGamma] = true;
 460            break;
 461     case kBCUTE:
 462            (*fIsCutVector)[kGamma] = true;
 463            break;
 464     case kBCUTM:
 465            (*fIsCutVector)[kGamma] = true;
 466            break;
 467     case kCUTELE:
 468            (*fIsCutVector)[kElectron] = true;
 469            break;
 470     case kDCUTE:
 471            (*fIsCutVector)[kElectron] = true;
 472            break;
 473     case kDCUTM:
 474            (*fIsCutVector)[kElectron] = true;
 475            break;
 476     case kCUTNEU:
 477            (*fIsCutVector)[kNeutralHadron] = true;
 478            break;
 479     case kCUTHAD:
 480            (*fIsCutVector)[kChargedHadron] = true;
 481            break;
 482     case kCUTMUO:
 483            (*fIsCutVector)[kMuon] = true;
 484            break;
 485     default:
 486            break;
 487   }
 488 }
 489
 490 void TG4PhysicsManager::SwitchIsFlagVector(TG3Flag flag)
 491 {
 492 // Updates the vector of booleans (fIsFlagVector) for the specified flag.
 493 // ---
 494
 495   switch (flag) {
 496     case kPAIR:
 497            // gamma
 498            (*fIsFlagVector)[kGamma] = true;
 499            break;
 500     case kCOMP:
 501            // gamma
 502            (*fIsFlagVector)[kGamma] = true;
 503            break;
 504     case kPHOT:
 505            // gamma
 506            (*fIsFlagVector)[kGamma] = true;
 507            break;
 508     case kPFIS:
 509            // gamma
 510            (*fIsFlagVector)[kGamma] = true;
 511            break;
 512     case kDRAY:
 513            // all charged particles
 514            (*fIsFlagVector)[kElectron] = true;
 515            (*fIsFlagVector)[kEplus] = true;
 516            (*fIsFlagVector)[kChargedHadron] = true;
 517            (*fIsFlagVector)[kMuon] = true;
 518            break;
 519     case kANNI:
 520            // e+ only
 521            (*fIsFlagVector)[kEplus] = true;
 522            break;
 523     case kBREM:
 524            // e-/e+, muons
 525            (*fIsFlagVector)[kElectron] = true;
 526            (*fIsFlagVector)[kEplus] = true;
 527            (*fIsFlagVector)[kMuon] = true;
 528            break;
 529     case kHADR:
 530            // hadrons
 531            (*fIsFlagVector)[kNeutralHadron] = true;
 532            (*fIsFlagVector)[kChargedHadron] = true;
 533            break;
 534     case kMUNU:
 535            // muons
 536            (*fIsFlagVector)[kMuon] = true;
 537            break;
 538     case kDCAY:
 539            // any
 540            (*fIsFlagVector)[kAny] = true;
 541            break;
 542     case kLOSS:
 543            // all charged particles
 544            (*fIsFlagVector)[kElectron] = true;
 545            (*fIsFlagVector)[kEplus] = true;
 546            (*fIsFlagVector)[kChargedHadron] = true;
 547            (*fIsFlagVector)[kMuon] = true;
 548            break;
 549     case kMULS:
 550            // all charged particles
 551            (*fIsFlagVector)[kElectron] = true;
 552            (*fIsFlagVector)[kEplus] = true;
 553            (*fIsFlagVector)[kChargedHadron] = true;
 554            (*fIsFlagVector)[kMuon] = true;
 555            break;
 556     default:
 557           break;
 558   }
 559 }
 560
 561 TG3Cut TG4PhysicsManager::GetG3Cut(G4String cutName)
 562 {
 563 // Retrieves corresponding TG3Cut constant from the cutName.
 564 // ---
 565
 566   if      (cutName == fG3CutNameVector[kCUTGAM]) return kCUTGAM;
 567   else if (cutName == fG3CutNameVector[kBCUTE])  return kBCUTE;
 568   else if (cutName == fG3CutNameVector[kBCUTM])  return kBCUTM;
 569   else if (cutName == fG3CutNameVector[kCUTELE]) return kCUTELE;
 570   else if (cutName == fG3CutNameVector[kDCUTE])  return kDCUTE;
 571   else if (cutName == fG3CutNameVector[kDCUTM])  return kDCUTM;
 572   else if (cutName == fG3CutNameVector[kCUTNEU]) return kCUTNEU;
 573   else if (cutName == fG3CutNameVector[kCUTHAD]) return kCUTHAD;
 574   else if (cutName == fG3CutNameVector[kCUTMUO]) return kCUTMUO;
 575   else return kNoG3Cuts;
 576 }
 577
 578 TG3Flag TG4PhysicsManager::GetG3Flag(G4String flagName)
 579 {
 580 // Retrieves corresponding TG3Flag constant from the flagName.
 581 // ---
 582
 583   if      (flagName == fG3FlagNameVector[kPAIR]) return kPAIR;
 584   else if (flagName == fG3FlagNameVector[kCOMP]) return kCOMP;
 585   else if (flagName == fG3FlagNameVector[kPHOT]) return kPHOT;
 586   else if (flagName == fG3FlagNameVector[kPFIS]) return kPFIS;
 587   else if (flagName == fG3FlagNameVector[kDRAY]) return kDRAY;
 588   else if (flagName == fG3FlagNameVector[kANNI]) return kANNI;
 589   else if (flagName == fG3FlagNameVector[kBREM]) return kBREM;
 590   else if (flagName == fG3FlagNameVector[kHADR]) return kHADR;
 591   else if (flagName == fG3FlagNameVector[kMUNU]) return kMUNU;
 592   else if (flagName == fG3FlagNameVector[kDCAY]) return kDCAY;
 593   else if (flagName == fG3FlagNameVector[kLOSS]) return kLOSS;
 594   else if (flagName == fG3FlagNameVector[kMULS]) return kMULS;
 595   else return kNoG3Flags;
 596 }
 597
 598 // public methods
 599
 600 void TG4PhysicsManager::BuildPhysics()
 601 {
 602 // Empty function - not needed in G4.
 603 // (Physics is built within /run/initialize.)
 604
 605   TG4Globals::Warning(
 606     "TG4PhysicsManager::BuildPhysics: is empty function in G4 MC.");
 607 }
 608
 609 void TG4PhysicsManager::SetCut(const char* cutName, Float_t cutValue)
 610 {
 611 // Sets the specified cut.
 612 // ---
 613
 614   if (fLock) LockException();
 615   TG3Cut g3Cut = GetG3Cut(cutName);
 616   if (g3Cut != kNoG3Cuts)
 617     SetCut(g3Cut, cutValue);
 618   else {
 619     G4String text = "TG4PhysicsManager::SetCut:\n";
 620     text = text + "    Parameter " + cutName;
 621     text = text + " is not implemented.";
 622     TG4Globals::Warning(text);
 623   }
 624 }
 625
 626 void TG4PhysicsManager::SetProcess(const char* flagName, Int_t flagValue)
 627 {
 628 // Sets the specified process control.
 629 // ---
 630
 631   if (fLock) LockException();
 632   TG3Flag g3Flag = GetG3Flag(flagName);
 633   if (g3Flag != kNoG3Flags)
 634     SetProcess(g3Flag, flagValue);
 635   else {
 636     G4String text = "TG4PhysicsManager::SetProcess:\n";
 637     text = text + "    Parameter " + flagName;
 638     text = text + " is not implemented.";
 639     TG4Globals::Warning(text);
 640   }
 641 }
 642
 643 void TG4PhysicsManager::SetProcessActivation()
 644 {
 645 // (In)Activates built processes according
 646 // to the setup in fFlagVector.
 647 // ---
 648
 649   if (fPhysicsList) {
 650     // temporarily excluded
 651     // fPhysicsList->SetProcessActivation();
 652   }
 653   else {
 654     G4String text = "TG4PhysicsManager::SetProcessActivation:\n";
 655     text = text +   "   There is no physics list set.";
 656     TG4Globals::Exception(text);
 657   }
 658 }
 659
 660 G4int TG4PhysicsManager::GetPDGEncodingFast(G4ParticleDefinition* particle)
 661 {
 662 // Returns the PDG code of particle;
 663 // if standard PDG code is not defined the preregistred
 664 // fParticlePDGMap is used.
 665 // ---
 666
 667   // get PDG encoding from G4 particle definition
 668   G4int pdgEncoding = particle->GetPDGEncoding();
 669
 670   if (pdgEncoding == 0) {
 671     // use FParticlePDGMap if standard PDG code is not defined
 672     G4String name = particle->GetParticleName();
 673     pdgEncoding = fParticlePDGMap.GetSecond(name);
 674   }
 675
 676   return pdgEncoding;
 677 }
 678
 679 G4bool TG4PhysicsManager::CheckCutWithCutVector(G4String name,
 680                              G4double value, TG3Cut& cut)
 681 {
 682 // Retrieves corresponding TG3Cut from the name and
 683 // in case the value is different from the value in cutVector
 684 // sets true the value of the fIsCutVector element
 685 // corresponding to this cut and returns true;
 686 // returns false otherwise.
 687 // ---
 688
 689   // convert cut name -> TG3Cut
 690   cut = GetG3Cut(name);
 691
 692   // set switch vector element only if the value
 693   // is different from the value in cutVector
 694   if (cut !=kNoG3Cuts) {
 695     // get tolerance from TG4G3Defaults in GeV
 696     G4double tolerance = TG4G3Defaults::CutTolerance()/GeV;
 697     if (!(fCutVector) || (abs(value - (*fCutVector)[cut]) > tolerance)) {
 698       SwitchIsCutVector(cut);
 699       return true;
 700     }
 701     else return false;
 702   }
 703   return false;
 704 }
 705
 706 G4bool TG4PhysicsManager::CheckFlagWithFlagVector(G4String name,
 707                               G4double value, TG3Flag& flag)
 708 {
 709 // Retrieves corresponding TG3Flag from the name and
 710 // in case the value is different from the value in flagVector
 711 // sets true the value of the fIsFlagVector element
 712 // corresponding to this flag and returns true;
 713 // returns false otherwise.
 714 // ---
 715
 716   // convert flag name -> TG3Flag
 717   flag = GetG3Flag(name);
 718
 719   // set switch vector element only if the value
 720   // is different from the value in flagVector
 721   if (flag !=kNoG3Flags) {
 722     if (!(fFlagVector) || (abs(value - (*fFlagVector)[flag]) > 0.01)) {
 723       SwitchIsFlagVector(flag);
 724       return true;
 725     }
 726     else return false;
 727   }
 728   return false;
 729 }
 730
 731 G4bool TG4PhysicsManager::CheckCutWithG3Defaults(G4String name,
 732                               G4double value, TG3Cut& cut)
 733 {
 734 // Retrieves corresponding TG3Cut from the name and
 735 // in case the value is different from the G3 default value
 736 // sets true the value of the SwitchCutVector element
 737 // corresponding to this cut and returns true;
 738 // returns false otherwise.
 739 // ---
 740
 741   // convert cut name -> TG3Cut
 742   cut = GetG3Cut(name);
 743
 744   // set switch vector element only if the value
 745   // is different from G3 default
 746   if (cut !=kNoG3Cuts) {
 747     if (!TG4G3Defaults::IsDefaultCut(cut, value)) {
 748       SwitchIsCutVector(cut);
 749       return true;
 750     }
 751     else return false;
 752   }
 753   return false;
 754 }
 755
 756 G4bool TG4PhysicsManager::CheckFlagWithG3Defaults(G4String name,
 757                               G4double value, TG3Flag& flag)
 758 {
 759 // Retrieves corresponding TG3Flag from the name and
 760 // in case the value is different from the G3 default value
 761 // sets true the value of the SwitchFlagVector element
 762 // corresponding to this flag and returns true;
 763 // returns false otherwise.
 764 // ---
 765
 766   // convert flag name -> TG3Flag
 767   flag = GetG3Flag(name);
 768
 769   // set switch vector element only if the value
 770   // is different from G3 default
 771   if (flag !=kNoG3Flags) {
 772     if (!TG4G3Defaults::IsDefaultFlag(flag, value)) {
 773       SwitchIsFlagVector(flag);
 774       return true;
 775     }
 776     else return false;
 777   }
 778   return false;
 779 }
 780
 781 void TG4PhysicsManager::SetG3DefaultCuts()
 782 {
 783 // Sets G3 default values of kinetic energy cuts.
 784 // ---
 785
 786   if (fLock) LockException();
 787   if (!fCutVector) {
 788     // create vector of kinetic energy cut values
 789     fCutVector = new TG4CutVector();
 790   }
 791   fCutVector->SetG3Defaults();
 792 }
 793
 794 void TG4PhysicsManager::SetG3DefaultProcesses()
 795 {
 796 // Sets G3 default values of control process flags.
 797 // ---
 798
 799   if (fLock) LockException();
 800   if (!fFlagVector) {
 801     // create vector of control process flag values
 802     fFlagVector = new TG4FlagVector;
 803   }
 804   fFlagVector->SetG3Defaults();
 805 }
 806
 807 G4bool TG4PhysicsManager::IsSpecialCuts() const
 808 {
 809 // Returns true if any special cut value is set.
 810 // ---
 811
 812   for (G4int i=0; i<kNofParticlesWSP; i++)
 813   {  if ((*fIsCutVector)[i]) return true; }
 814
 815   return false;
 816 }
 817
 818 G4bool TG4PhysicsManager::IsSpecialFlags() const
 819 {
 820 // Returns true if any special flag value is set.
 821 // ---
 822
 823   for (G4int i=0; i<kNofParticlesWSP; i++)
 824   {  if ((*fIsFlagVector)[i]) return true; }
 825
 826   return false;
 827 }
 828
 829 TG3ParticleWSP TG4PhysicsManager::GetG3ParticleWSP(
 830                                   G4ParticleDefinition* particle) const
 831 {
 832 // Returns TG3ParticleWSP constant for the specified particle.
 833 // (See TG3ParticleWSP.h, too.)
 834 // ---
 835
 836   G4String name = particle->GetParticleName();
 837   G4String pType = particle->GetParticleType();
 838
 839   if (name == "gamma") {
 840     return kGamma;
 841   }
 842   else if (name == "e-") {
 843     return kElectron;
 844   }
 845   else if (name == "e+") {
 846     return kEplus;
 847   }
 848   else if (( pType == "baryon" || pType == "meson" || pType == "nucleus" )) {
 849     if (particle->GetPDGCharge() == 0) {
 850       return kNeutralHadron;
 851     }
 852     else
 853       return kChargedHadron;
 854   }
 855   else if ( name == "mu-" || name == "mu+" ) {
 856     return kMuon;
 857   }
 858   else {
 859     return kNofParticlesWSP;
 860   }
 861 }
 862
 863 void TG4PhysicsManager::GetG3ParticleWSPName(G4int particleWSP,
 864                                              G4String& name) const
 865 {
 866 // Fills the passed name with the name/type of particle specified
 867 // by TG3ParticleWSP constant.
 868 // (See TG3ParticleWSP.h, too.)
 869 // ---
 870
 871   switch (particleWSP) {
 872     case kGamma:
 873       name = "Gamma";
 874       break;
 875     case kElectron:
 876       name = "Electron";
 877       break;
 878     case kEplus:
 879       name = "Eplus";
 880       break;
 881     case kNeutralHadron:
 882       name = "NeutralHadron";
 883       break;
 884     case kChargedHadron:
 885       name = "ChargedHadron";
 886       break;
 887     case kMuon:
 888       name = "Muon";
 889       break;
 890     case kAny:
 891       name = "Any";
 892       break;
 893     case kNofParticlesWSP:
 894       name = "NoSP";
 895       break;
 896     default:
 897       G4String text = "TG4PhysicsList::GetG3ParticleWSPName:\n";
 898       text = text + "   Wrong particleWSP.";
 899       TG4Globals::Exception(text);
 900       break;
 901   }
 902 }
 903