Geant4 Cross Reference |
1 // 2 // ******************************************************************** 3 // * License and Disclaimer * 4 // * * 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. * 10 // * * 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitation of liability. * 17 // * * 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************************************** 25 // 26 // 27 // ------------------------------------------------------------------- 28 // 29 // GEANT4 Class file 30 // 31 // 32 // File name: G4RDGenerator2BS 33 // 34 // Author: Andreia Trindade (andreia@lip.pt) 35 // Pedro Rodrigues (psilva@lip.pt) 36 // Luis Peralta (luis@lip.pt) 37 // 38 // Creation date: 2 June 2003 39 // 40 // Modifications: 41 // 02 Jun 2003 First implementation acording with new design 42 // 05 Nov 2003 MGP Fixed std namespace 43 // 17 Nov 2003 MGP Fixed compilation problem on Windows 44 // 45 // Class Description: 46 // 47 // Concrete base class for Bremsstrahlung Angular Distribution Generation - 2BS Distribution 48 // 49 // Class Description: End 50 // 51 // ------------------------------------------------------------------- 52 // 53 // 54 55 #include "G4RDGenerator2BS.hh" 56 #include "Randomize.hh" 57 #include "G4PhysicalConstants.hh" 58 // 59 60 G4RDGenerator2BS::G4RDGenerator2BS(const G4String& name):G4RDVBremAngularDistribution(name) 61 {;} 62 63 // 64 65 G4RDGenerator2BS::~G4RDGenerator2BS() 66 {;} 67 68 // 69 70 G4double G4RDGenerator2BS::PolarAngle(const G4double initial_energy, 71 const G4double final_energy, 72 const G4int Z) 73 { 74 75 // Adapted from "Improved bremsstrahlung photon angular sampling in the EGS4 code system" 76 // by Alex F. Bielajew, Rahde Mohan anc Chen-Shou Chui, PIRS-0203 77 // Ionizing Radiation Standards 78 // Institute for National Measurement Standards 79 // National Research Council of Canada 80 // Departement of Medical Physics, Memorial Sloan-Kettering Cancer Center, New York 81 82 83 G4double theta = 0; 84 85 G4double initialTotalEnergy = (initial_energy+electron_mass_c2)/electron_mass_c2; 86 G4double finalTotalEnergy = (final_energy+electron_mass_c2)/electron_mass_c2; 87 EnergyRatio = finalTotalEnergy/initialTotalEnergy; 88 G4double gMaxEnergy = (pi*initialTotalEnergy)*(pi*initialTotalEnergy); 89 90 G4double Zeff = std::sqrt(static_cast<G4double>(Z) * (static_cast<G4double>(Z) + 1.0)); 91 z = (0.00008116224*(std::pow(Zeff,0.3333333))); 92 93 // Rejection arguments 94 rejection_argument1 = (1.0+EnergyRatio*EnergyRatio); 95 rejection_argument2 = - 2.0*EnergyRatio + 3.0*rejection_argument1; 96 rejection_argument3 = ((1-EnergyRatio)/(2.0*initialTotalEnergy*EnergyRatio))* 97 ((1-EnergyRatio)/(2.0*initialTotalEnergy*EnergyRatio)); 98 99 // Calculate rejection function at 0, 1 and Emax 100 G4double gfunction0 = RejectionFunction(0); 101 G4double gfunction1 = RejectionFunction(1); 102 G4double gfunctionEmax = RejectionFunction(gMaxEnergy); 103 104 105 // Calculate Maximum value 106 G4double gMaximum = std::max(gfunction0,gfunction1); 107 gMaximum = std::max(gMaximum,gfunctionEmax); 108 109 G4double rand, gfunctionTest, randTest; 110 111 do{ 112 rand = G4UniformRand(); 113 rand = rand/(1-rand+1.0/gMaxEnergy); 114 gfunctionTest = RejectionFunction(rand); 115 randTest = G4UniformRand(); 116 117 }while(randTest > (gfunctionTest/gMaximum)); 118 119 theta = std::sqrt(rand)/initialTotalEnergy; 120 121 122 return theta; 123 } 124 // 125 126 G4double G4RDGenerator2BS::RejectionFunction(G4double value) const 127 { 128 129 G4double argument = (1+value)*(1+value); 130 131 G4double gfunction = (4+std::log(rejection_argument3+(z/argument)))* 132 ((4*EnergyRatio*value/argument)-rejection_argument1)+rejection_argument2; 133 134 return gfunction; 135 136 } 137 138 void G4RDGenerator2BS::PrintGeneratorInformation() const 139 { 140 G4cout << "\n" << G4endl; 141 G4cout << "Bremsstrahlung Angular Generator is 2BS Generator from 2BS Koch & Motz distribution (Rev Mod Phys 31(4), 920 (1959))" << G4endl; 142 G4cout << "Sampling algorithm adapted from PIRS-0203" << G4endl; 143 G4cout << "\n" << G4endl; 144 } 145 146