Geant4 Cross Reference |
1 // 1 2 // ******************************************* 3 // * License and Disclaimer 4 // * 5 // * The Geant4 software is copyright of th 6 // * the Geant4 Collaboration. It is provided 7 // * conditions of the Geant4 Software License 8 // * LICENSE and available at http://cern.ch/ 9 // * include a list of copyright holders. 10 // * 11 // * Neither the authors of this software syst 12 // * institutes,nor the agencies providing fin 13 // * work make any representation or warran 14 // * regarding this software system or assum 15 // * use. Please see the license in the file 16 // * for the full disclaimer and the limitatio 17 // * 18 // * This code implementation is the result 19 // * technical work of the GEANT4 collaboratio 20 // * By using, copying, modifying or distri 21 // * any work based on the software) you ag 22 // * use in resulting scientific publicati 23 // * acceptance of all terms of the Geant4 Sof 24 // ******************************************* 25 // 26 // ------------------------------------------- 27 // 28 // GEANT4 Class header file 29 // 30 // 31 // File name: G4VEmModel 32 // 33 // Author: Vladimir Ivanchenko 34 // 35 // Creation date: 03.01.2002 36 // 37 // Modifications: 38 // 39 // 23-12-02 V.Ivanchenko change interface befo 40 // 24-01-03 Cut per region (V.Ivanchenko) 41 // 13-02-03 Add name (V.Ivanchenko) 42 // 25-02-03 Add sample theta and displacement 43 // 23-07-03 Replace G4Material by G4MaterialCu 44 // calculation (V.Ivanchenko) 45 // 01-03-04 L.Urban signature changed in Sampl 46 // 23-04-04 L.urban signature of SampleCosineT 47 // 17-11-04 Add method CrossSectionPerAtom (V. 48 // 14-03-05 Reduce number of pure virtual meth 49 // separate (V.Ivanchenko) 50 // 24-03-05 Remove IsInCharge and add G4VParti 51 // 08-04-05 Major optimisation of internal int 52 // 15-04-05 optimize internal interface for ms 53 // 08-05-05 A -> N (V.Ivanchenko) 54 // 25-07-05 Move constructor and destructor to 55 // 02-02-06 ComputeCrossSectionPerAtom: defaul 56 // 06-02-06 add method ComputeMeanFreePath() ( 57 // 07-03-06 Optimize msc methods (V.Ivanchenko 58 // 29-06-06 Add member currentElement and Get/ 59 // 29-10-07 Added SampleScattering (V.Ivanchen 60 // 15-07-08 Reorder class members and improve 61 // 21-07-08 Added vector of G4ElementSelector 62 // 12-09-08 Added methods GetParticleCharge, G 63 // CorrectionsAlongStep, ActivateNucl 64 // 16-02-09 Moved implementations of virtual m 65 // 07-04-09 Moved msc methods from G4VEmModel 66 // 13-10-10 Added G4VEmAngularDistribution (VI 67 // 68 // Class Description: 69 // 70 // Abstract interface to energy loss models 71 72 // ------------------------------------------- 73 // 74 75 #ifndef G4VEmModel_h 76 #define G4VEmModel_h 1 77 78 #include "globals.hh" 79 #include "G4DynamicParticle.hh" 80 #include "G4ParticleDefinition.hh" 81 #include "G4MaterialCutsCouple.hh" 82 #include "G4Material.hh" 83 #include "G4Element.hh" 84 #include "G4ElementVector.hh" 85 #include "G4Isotope.hh" 86 #include "G4DataVector.hh" 87 #include "G4VEmFluctuationModel.hh" 88 #include "G4VEmAngularDistribution.hh" 89 #include "G4EmElementSelector.hh" 90 #include <CLHEP/Random/RandomEngine.h> 91 #include <vector> 92 93 class G4ElementData; 94 class G4PhysicsTable; 95 class G4Region; 96 class G4VParticleChange; 97 class G4ParticleChangeForLoss; 98 class G4ParticleChangeForGamma; 99 class G4Track; 100 class G4LossTableManager; 101 102 class G4VEmModel 103 { 104 105 public: 106 107 explicit G4VEmModel(const G4String& nam); 108 109 virtual ~G4VEmModel(); 110 111 //------------------------------------------ 112 // Virtual methods to be implemented for any 113 //------------------------------------------ 114 115 virtual void Initialise(const G4ParticleDefi 116 117 virtual void SampleSecondaries(std::vector<G 118 const G4Mater 119 const G4Dynam 120 G4double tmin 121 G4double tmax 122 123 //------------------------------------------ 124 // Methods for initialisation of MT; may be 125 //------------------------------------------ 126 127 // initialisation in local thread 128 virtual void InitialiseLocal(const G4Particl 129 G4VEmModel* mas 130 131 // initialisation of a new material at run t 132 virtual void InitialiseForMaterial(const G4P 133 const G4M 134 135 // initialisation of a new element at run ti 136 virtual void InitialiseForElement(const G4Pa 137 G4int Z); 138 139 //------------------------------------------ 140 // Methods with standard implementation; may 141 //------------------------------------------ 142 143 // main method to compute dEdx 144 virtual G4double ComputeDEDXPerVolume(const 145 const 146 G4doub 147 G4doub 148 149 // main method to compute cross section per 150 virtual G4double CrossSectionPerVolume(const 151 const 152 G4dou 153 G4dou 154 G4dou 155 156 // method to get partial cross section 157 virtual G4double GetPartialCrossSection(cons 158 G4in 159 cons 160 G4do 161 162 // main method to compute cross section per 163 virtual G4double ComputeCrossSectionPerAtom( 164 165 166 167 168 169 170 // main method to compute cross section per 171 virtual G4double ComputeCrossSectionPerShell 172 173 174 175 176 177 // Compute effective ion charge square 178 virtual G4double ChargeSquareRatio(const G4T 179 180 // Compute effective ion charge square 181 virtual G4double GetChargeSquareRatio(const 182 const 183 G4doub 184 185 // Compute ion charge 186 virtual G4double GetParticleCharge(const G4P 187 const G4M 188 G4double 189 190 // Initialisation for a new track 191 virtual void StartTracking(G4Track*); 192 193 // add correction to energy loss and compute 194 virtual void CorrectionsAlongStep(const G4Ma 195 const G4Dy 196 const G4do 197 G4double& 198 199 // value which may be tabulated (by default 200 virtual G4double Value(const G4MaterialCutsC 201 const G4ParticleDefin 202 G4double kineticEnerg 203 204 // threshold for zero value 205 virtual G4double MinPrimaryEnergy(const G4Ma 206 const G4Pa 207 G4double c 208 209 // model can define low-energy limit for the 210 virtual G4double MinEnergyCut(const G4Partic 211 const G4Materi 212 213 // initialisation at run time for a given ma 214 virtual void SetupForMaterial(const G4Partic 215 const G4Materi 216 G4double kinet 217 218 // add a region for the model 219 virtual void DefineForRegion(const G4Region* 220 221 // fill number of different type of secondar 222 virtual void FillNumberOfSecondaries(G4int& 223 G4int& 224 225 // for automatic documentation 226 virtual void ModelDescription(std::ostream& 227 228 protected: 229 230 // initialisation of the ParticleChange for 231 G4ParticleChangeForLoss* GetParticleChangeFo 232 233 // initialisation of the ParticleChange for 234 G4ParticleChangeForGamma* GetParticleChangeF 235 236 // kinematically allowed max kinetic energy 237 virtual G4double MaxSecondaryEnergy(const G4 238 G4double 239 240 public: 241 242 //------------------------------------------ 243 // Generic methods common to all models 244 //------------------------------------------ 245 246 // should be called at initialisation to bui 247 void InitialiseElementSelectors(const G4Part 248 const G4Data 249 250 // should be called at initialisation to acc 251 inline std::vector<G4EmElementSelector*>* Ge 252 253 // should be called at initialisation to set 254 inline void SetElementSelectors(std::vector< 255 256 // dEdx per unit length, base material appro 257 inline G4double ComputeDEDX( const G4Materia 258 const G4Particl 259 G4double kineti 260 G4double cutEne 261 262 // cross section per volume, base material a 263 inline G4double CrossSection(const G4Materia 264 const G4Particl 265 G4double kineti 266 G4double cutEne 267 G4double maxEne 268 269 // compute mean free path via cross section 270 inline G4double ComputeMeanFreePath(const G4 271 G4double 272 const G4 273 G4double 274 G4double 275 276 // generic cross section per element 277 inline G4double ComputeCrossSectionPerAtom(c 278 c 279 G 280 G 281 G 282 283 // atom can be selected effitiantly if eleme 284 inline const G4Element* SelectRandomAtom(con 285 con 286 G4d 287 G4d 288 G4d 289 // same as SelectRandomAtom above but more e 290 inline const G4Element* SelectTargetAtom(con 291 con 292 G4d 293 G4d 294 G4d 295 G4d 296 297 // to select atom cross section per volume i 298 const G4Element* SelectRandomAtom(const G4Ma 299 const G4Pa 300 G4double k 301 G4double c 302 G4double m 303 304 // to select atom if cross section is propor 305 const G4Element* GetCurrentElement(const G4M 306 G4int SelectRandomAtomNumber(const G4Materia 307 308 // select isotope in order to have precise m 309 const G4Isotope* GetCurrentIsotope(const G4E 310 G4int SelectIsotopeNumber(const G4Element*) 311 312 //------------------------------------------ 313 // Get/Set methods 314 //------------------------------------------ 315 316 void SetParticleChange(G4VParticleChange*, G 317 318 void SetCrossSectionTable(G4PhysicsTable*, G 319 320 inline G4ElementData* GetElementData(); 321 322 inline G4PhysicsTable* GetCrossSectionTable( 323 324 inline G4VEmFluctuationModel* GetModelOfFluc 325 326 inline G4VEmAngularDistribution* GetAngularD 327 328 inline G4VEmModel* GetTripletModel(); 329 330 inline void SetTripletModel(G4VEmModel*); 331 332 inline void SetAngularDistribution(G4VEmAngu 333 334 inline G4double HighEnergyLimit() const; 335 336 inline G4double LowEnergyLimit() const; 337 338 inline G4double HighEnergyActivationLimit() 339 340 inline G4double LowEnergyActivationLimit() c 341 342 inline G4double PolarAngleLimit() const; 343 344 inline G4double SecondaryThreshold() const; 345 346 inline G4bool DeexcitationFlag() const; 347 348 inline G4bool ForceBuildTableFlag() const; 349 350 inline G4bool UseAngularGeneratorFlag() cons 351 352 inline void SetAngularGeneratorFlag(G4bool); 353 354 inline void SetHighEnergyLimit(G4double); 355 356 inline void SetLowEnergyLimit(G4double); 357 358 inline void SetActivationHighEnergyLimit(G4d 359 360 inline void SetActivationLowEnergyLimit(G4do 361 362 inline G4bool IsActive(G4double kinEnergy) c 363 364 inline void SetPolarAngleLimit(G4double); 365 366 inline void SetSecondaryThreshold(G4double); 367 368 inline void SetDeexcitationFlag(G4bool val); 369 370 inline void SetForceBuildTable(G4bool val); 371 372 inline void SetFluctuationFlag(G4bool val); 373 374 inline G4bool IsMaster() const; 375 376 inline void SetUseBaseMaterials(G4bool val); 377 378 inline G4bool UseBaseMaterials() const; 379 380 inline G4double MaxSecondaryKinEnergy(const 381 382 inline const G4String& GetName() const; 383 384 inline void SetCurrentCouple(const G4Materia 385 386 inline G4bool IsLocked() const; 387 388 inline void SetLocked(G4bool); 389 390 // obsolete methods 391 [[deprecated("Use G4EmParameters::Instance() 392 void SetLPMFlag(G4bool); 393 394 void SetMasterThread(G4bool); 395 396 // hide assignment operator 397 G4VEmModel & operator=(const G4VEmModel &ri 398 G4VEmModel(const G4VEmModel&) = delete; 399 400 protected: 401 402 inline const G4MaterialCutsCouple* CurrentCo 403 404 inline void SetCurrentElement(const G4Elemen 405 406 private: 407 408 // ======== Parameters of the class fixed at 409 410 G4VEmFluctuationModel* flucModel = null 411 G4VEmAngularDistribution* anglModel = null 412 G4VEmModel* fTripletModel = 413 const G4MaterialCutsCouple* fCurrentCouple = 414 const G4Element* fCurrentElement 415 std::vector<G4EmElementSelector*>* elmSelect 416 G4LossTableManager* fEmManager; 417 418 protected: 419 420 G4ElementData* fElementData = 421 G4VParticleChange* pParticleChange 422 G4PhysicsTable* xSectionTable = 423 const G4Material* pBaseMaterial = 424 const std::vector<G4double>* theDensityFacto 425 const std::vector<G4int>* theDensityIdx = 426 427 G4double inveplus; 428 G4double pFactor = 1.0; 429 430 private: 431 432 G4double lowLimit; 433 G4double highLimit; 434 G4double eMinActive = 0.0; 435 G4double eMaxActive = DBL_MAX; 436 G4double secondaryThreshold = DBL_MAX; 437 G4double polarAngleLimit; 438 439 G4int nSelectors = 0; 440 G4int nsec = 5; 441 442 protected: 443 444 std::size_t currentCoupleIndex = 0; 445 std::size_t basedCoupleIndex = 0; 446 G4bool lossFlucFlag = true; 447 448 private: 449 450 G4bool flagDeexcitation = false; 451 G4bool flagForceBuildTable = false; 452 G4bool isMaster = true; 453 454 G4bool localTable = true; 455 G4bool localElmSelectors = true; 456 G4bool useAngularGenerator = false; 457 G4bool useBaseMaterials = false; 458 G4bool isLocked = false; 459 460 const G4String name; 461 std::vector<G4double> xsec; 462 463 }; 464 465 // ======== Run time inline methods ========== 466 467 inline void G4VEmModel::SetCurrentCouple(const 468 { 469 if(fCurrentCouple != ptr) { 470 fCurrentCouple = ptr; 471 basedCoupleIndex = currentCoupleIndex = pt 472 pBaseMaterial = ptr->GetMaterial(); 473 pFactor = 1.0; 474 if(useBaseMaterials) { 475 basedCoupleIndex = (*theDensityIdx)[curr 476 if(nullptr != pBaseMaterial->GetBaseMate 477 pBaseMaterial = pBaseMaterial->GetBaseMateri 478 pFactor = (*theDensityFactor)[currentCou 479 } 480 } 481 } 482 483 //....oooOO0OOooo........oooOO0OOooo........oo 484 485 inline const G4MaterialCutsCouple* G4VEmModel: 486 { 487 return fCurrentCouple; 488 } 489 490 //....oooOO0OOooo........oooOO0OOooo........oo 491 492 inline void G4VEmModel::SetCurrentElement(cons 493 { 494 fCurrentElement = elm; 495 } 496 497 //....oooOO0OOooo........oooOO0OOooo........oo 498 499 inline 500 G4double G4VEmModel::MaxSecondaryKinEnergy(con 501 { 502 return MaxSecondaryEnergy(dynPart->GetPartic 503 dynPart->GetKineti 504 } 505 506 //....oooOO0OOooo........oooOO0OOooo........oo 507 508 inline G4double G4VEmModel::ComputeDEDX(const 509 const 510 G4doub 511 G4doub 512 { 513 SetCurrentCouple(couple); 514 return pFactor*ComputeDEDXPerVolume(pBaseMat 515 } 516 517 //....oooOO0OOooo........oooOO0OOooo........oo 518 519 inline G4double G4VEmModel::CrossSection(const 520 const 521 G4dou 522 G4dou 523 G4dou 524 { 525 SetCurrentCouple(couple); 526 return pFactor*CrossSectionPerVolume(pBaseMa 527 cutEner 528 } 529 530 //....oooOO0OOooo........oooOO0OOooo........oo 531 532 inline 533 G4double G4VEmModel::ComputeMeanFreePath(const 534 G4dou 535 const 536 G4dou 537 G4dou 538 { 539 G4double cross = CrossSectionPerVolume(mater 540 return (cross > 0.0) ? 1./cross : DBL_MAX; 541 } 542 543 //....oooOO0OOooo........oooOO0OOooo........oo 544 545 inline G4double 546 G4VEmModel::ComputeCrossSectionPerAtom(const G 547 const G 548 G4doubl 549 G4doubl 550 G4doubl 551 { 552 fCurrentElement = elm; 553 return ComputeCrossSectionPerAtom(part,kinEn 554 cutEnergy, 555 } 556 557 //....oooOO0OOooo........oooOO0OOooo........oo 558 559 inline const G4Element* 560 G4VEmModel::SelectRandomAtom(const G4MaterialC 561 const G4ParticleD 562 G4double kinEnerg 563 G4double cutEnerg 564 G4double maxEnerg 565 { 566 SetCurrentCouple(couple); 567 fCurrentElement = (nSelectors > 0) ? 568 ((*elmSelectors)[couple->GetIndex()])->Sel 569 SelectRandomAtom(pBaseMaterial,part,kinEne 570 return fCurrentElement; 571 } 572 573 //....oooOO0OOooo........oooOO0OOooo........oo 574 575 inline const G4Element* 576 G4VEmModel::SelectTargetAtom(const G4MaterialC 577 const G4ParticleD 578 G4double kinEnerg 579 G4double logKinE, 580 G4double cutEnerg 581 G4double maxEnerg 582 { 583 SetCurrentCouple(couple); 584 fCurrentElement = (nSelectors > 0) 585 ? ((*elmSelectors)[couple->GetIndex()])->Se 586 : SelectRandomAtom(pBaseMaterial,part,kinEn 587 return fCurrentElement; 588 } 589 590 // ======== Get/Set inline methods used at ini 591 592 inline G4VEmFluctuationModel* G4VEmModel::GetM 593 { 594 return flucModel; 595 } 596 597 //....oooOO0OOooo........oooOO0OOooo........oo 598 599 inline G4VEmAngularDistribution* G4VEmModel::G 600 { 601 return anglModel; 602 } 603 604 //....oooOO0OOooo........oooOO0OOooo........oo 605 606 inline void G4VEmModel::SetAngularDistribution 607 { 608 if(p != anglModel) { 609 delete anglModel; 610 anglModel = p; 611 } 612 } 613 614 //....oooOO0OOooo........oooOO0OOooo........oo 615 616 inline G4VEmModel* G4VEmModel::GetTripletModel 617 { 618 return fTripletModel; 619 } 620 621 //....oooOO0OOooo........oooOO0OOooo........oo 622 623 inline void G4VEmModel::SetTripletModel(G4VEmM 624 { 625 if(p != fTripletModel) { 626 delete fTripletModel; 627 fTripletModel = p; 628 } 629 } 630 631 //....oooOO0OOooo........oooOO0OOooo........oo 632 633 inline G4double G4VEmModel::HighEnergyLimit() 634 { 635 return highLimit; 636 } 637 638 //....oooOO0OOooo........oooOO0OOooo........oo 639 640 inline G4double G4VEmModel::LowEnergyLimit() c 641 { 642 return lowLimit; 643 } 644 645 //....oooOO0OOooo........oooOO0OOooo........oo 646 647 inline G4double G4VEmModel::HighEnergyActivati 648 { 649 return eMaxActive; 650 } 651 652 //....oooOO0OOooo........oooOO0OOooo........oo 653 654 inline G4double G4VEmModel::LowEnergyActivatio 655 { 656 return eMinActive; 657 } 658 659 //....oooOO0OOooo........oooOO0OOooo........oo 660 661 inline G4double G4VEmModel::PolarAngleLimit() 662 { 663 return polarAngleLimit; 664 } 665 666 //....oooOO0OOooo........oooOO0OOooo........oo 667 668 inline G4double G4VEmModel::SecondaryThreshold 669 { 670 return secondaryThreshold; 671 } 672 673 //....oooOO0OOooo........oooOO0OOooo........oo 674 675 inline G4bool G4VEmModel::DeexcitationFlag() c 676 { 677 return flagDeexcitation; 678 } 679 680 //....oooOO0OOooo........oooOO0OOooo........oo 681 682 inline G4bool G4VEmModel::ForceBuildTableFlag( 683 { 684 return flagForceBuildTable; 685 } 686 687 //....oooOO0OOooo........oooOO0OOooo........oo 688 689 inline G4bool G4VEmModel::UseAngularGeneratorF 690 { 691 return useAngularGenerator; 692 } 693 694 //....oooOO0OOooo........oooOO0OOooo........oo 695 696 inline void G4VEmModel::SetAngularGeneratorFla 697 { 698 useAngularGenerator = val; 699 } 700 701 //....oooOO0OOooo........oooOO0OOooo........oo 702 703 inline void G4VEmModel::SetFluctuationFlag(G4b 704 { 705 lossFlucFlag = val; 706 } 707 708 //....oooOO0OOooo........oooOO0OOooo........oo 709 710 inline G4bool G4VEmModel::IsMaster() const 711 { 712 return isMaster; 713 } 714 715 //....oooOO0OOooo........oooOO0OOooo........oo 716 717 inline void G4VEmModel::SetUseBaseMaterials(G4 718 { 719 useBaseMaterials = val; 720 } 721 722 //....oooOO0OOooo........oooOO0OOooo........oo 723 724 inline G4bool G4VEmModel::UseBaseMaterials() c 725 { 726 return useBaseMaterials; 727 } 728 729 //....oooOO0OOooo........oooOO0OOooo........oo 730 731 inline void G4VEmModel::SetHighEnergyLimit(G4d 732 { 733 highLimit = val; 734 } 735 736 //....oooOO0OOooo........oooOO0OOooo........oo 737 738 inline void G4VEmModel::SetLowEnergyLimit(G4do 739 { 740 lowLimit = val; 741 } 742 743 //....oooOO0OOooo........oooOO0OOooo........oo 744 745 inline void G4VEmModel::SetActivationHighEnerg 746 { 747 eMaxActive = val; 748 } 749 750 //....oooOO0OOooo........oooOO0OOooo........oo 751 752 inline void G4VEmModel::SetActivationLowEnergy 753 { 754 eMinActive = val; 755 } 756 757 //....oooOO0OOooo........oooOO0OOooo........oo 758 759 inline G4bool G4VEmModel::IsActive(G4double ki 760 { 761 return (kinEnergy >= eMinActive && kinEnergy 762 } 763 764 //....oooOO0OOooo........oooOO0OOooo........oo 765 766 inline void G4VEmModel::SetPolarAngleLimit(G4d 767 { 768 if(!isLocked) { polarAngleLimit = val; } 769 } 770 771 //....oooOO0OOooo........oooOO0OOooo........oo 772 773 inline void G4VEmModel::SetSecondaryThreshold( 774 { 775 secondaryThreshold = val; 776 } 777 778 //....oooOO0OOooo........oooOO0OOooo........oo 779 780 inline void G4VEmModel::SetDeexcitationFlag(G4 781 { 782 flagDeexcitation = val; 783 } 784 785 //....oooOO0OOooo........oooOO0OOooo........oo 786 787 inline void G4VEmModel::SetForceBuildTable(G4b 788 { 789 flagForceBuildTable = val; 790 } 791 792 //....oooOO0OOooo........oooOO0OOooo........oo 793 794 inline const G4String& G4VEmModel::GetName() c 795 { 796 return name; 797 } 798 799 //....oooOO0OOooo........oooOO0OOooo........oo 800 801 inline std::vector<G4EmElementSelector*>* G4VE 802 { 803 return elmSelectors; 804 } 805 806 //....oooOO0OOooo........oooOO0OOooo........oo 807 808 inline void 809 G4VEmModel::SetElementSelectors(std::vector<G4 810 { 811 if(p != elmSelectors) { 812 elmSelectors = p; 813 nSelectors = (nullptr != elmSelectors) ? G 814 localElmSelectors = false; 815 } 816 } 817 818 //....oooOO0OOooo........oooOO0OOooo........oo 819 820 inline G4ElementData* G4VEmModel::GetElementDa 821 { 822 return fElementData; 823 } 824 825 //....oooOO0OOooo........oooOO0OOooo........oo 826 827 inline G4PhysicsTable* G4VEmModel::GetCrossSec 828 { 829 return xSectionTable; 830 } 831 832 //....oooOO0OOooo........oooOO0OOooo........oo 833 834 inline G4bool G4VEmModel::IsLocked() const 835 { 836 return isLocked; 837 } 838 839 //....oooOO0OOooo........oooOO0OOooo........oo 840 841 inline void G4VEmModel::SetLocked(G4bool val) 842 { 843 isLocked = val; 844 } 845 846 //....oooOO0OOooo........oooOO0OOooo........oo 847 848 #endif 849