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 // Author: Alfonso Mantero (Alfonso.Mantero@ge 29 // 30 // History: 31 // ----------- 32 // 02 Sep 2003 Alfonso Mantero created 33 // 34 // ------------------------------------------- 35 36 #include "XrayFluoMercuryPrimaryGeneratorActio 37 #include "XrayFluoMercuryDetectorConstruction. 38 #include "XrayFluoMercuryPrimaryGeneratorMesse 39 #include "XrayFluoRunAction.hh" 40 #include "XrayFluoAnalysisManager.hh" 41 #include "XrayFluoDataSet.hh" 42 #include "G4PhysicalConstants.hh" 43 #include "G4SystemOfUnits.hh" 44 #include "G4DataVector.hh" 45 #include "G4Event.hh" 46 #include "G4ParticleGun.hh" 47 #include "G4ParticleTable.hh" 48 #include "G4ParticleDefinition.hh" 49 #include "Randomize.hh" 50 51 //....oooOO0OOooo........oooOO0OOooo........oo 52 53 XrayFluoMercuryPrimaryGeneratorAction::XrayFlu 54 :globalFlag(false),spectrum("off") 55 { 56 57 XrayFluoDetector = XrayFluoDC; 58 59 G4int n_particle = 1; 60 particleGun = new G4ParticleGun(n_particle) 61 62 //create a messenger for this class 63 gunMessenger = new XrayFluoMercuryPrimaryGen 64 runManager = new XrayFluoRunAction(); 65 66 // default particle kinematic 67 68 G4ParticleTable* particleTable = G4ParticleT 69 G4String particleName; 70 G4ParticleDefinition* particle 71 = particleTable->FindParticle(particleName 72 particleGun->SetParticleDefinition(particle) 73 particleGun->SetParticleMomentumDirection(G4 74 75 76 particleGun->SetParticleEnergy(10.*keV); 77 G4double position = -0.5*(XrayFluoDetector-> 78 particleGun->SetParticlePosition(G4ThreeVect 79 80 G4cout << "XrayFluoMercuryPrimaryGeneratorAc 81 82 } 83 84 85 //....oooOO0OOooo........oooOO0OOooo........oo 86 87 XrayFluoMercuryPrimaryGeneratorAction::~XrayFl 88 { 89 delete particleGun; 90 delete gunMessenger; 91 delete runManager; 92 93 G4cout << "XrayFluoMercuryPrimaryGeneratorAc 94 95 } 96 97 //....oooOO0OOooo........oooOO0OOooo........oo 98 99 void XrayFluoMercuryPrimaryGeneratorAction::Ge 100 { 101 //this function is called at the begining of 102 // 103 104 // Conidering the sunas a Poin-like source. 105 106 G4double z0 = -0.5*(XrayFluoDetector->GetWor 107 G4double y0 = 0.*m, x0 = 0.*m; 108 109 110 // Let's try to illuminate only the prtion o 111 112 G4double spacecraftLatitude = XrayFluoDetect 113 G4double mercuryDia = XrayFluoDetector->GetM 114 G4double sunDia = XrayFluoDetector->GetSunDi 115 G4double opticField = XrayFluoDetector->GetO 116 117 118 G4double a = 2*std::tan(opticField/2); 119 120 // if (!pointLikeFlag) { 121 122 // let's decide from wich point of the sun 123 124 G4double theta = std::acos(2.*G4UniformRand( 125 G4double phi = 2. * pi * G4UniformRand(); 126 G4double rho = sunDia/2; 127 128 G4double sunPosX = x0 + rho * std::sin(the 129 G4double sunPosY = y0 + rho * std::sin(the 130 G4double sunPosZ = z0 + rho * std::cos(the 131 132 particleGun->SetParticlePosition(G4ThreeVe 133 134 // the angle at the center of Mercury subt 135 G4double alpha = 2 * a/mercuryDia; 136 137 if(!globalFlag){ 138 theta = alpha * G4UniformRand() + (180.* 139 phi = alpha * G4UniformRand() + 90. * de 140 } 141 142 else if(globalFlag){ 143 theta = pi/2. * rad * G4UniformRand() + 144 phi = 2*pi*rad * G4UniformRand() ; 145 } 146 147 rho = mercuryDia/2.; 148 149 G4double mercuryPosX = rho * std::sin(thet 150 G4double mercuryPosY = rho * std::sin(thet 151 G4double mercuryPosZ = rho * std::cos(thet 152 153 particleGun->SetParticleMomentumDirection( 154 G4ThreeVector(mercuryPosX-sunPosX ,m 155 156 // } 157 // if (pointLikeFlag) { 158 159 // // theta is the angle that the mean direc 160 // // point of the surface of Mercury) makes 161 // G4double theta = std::asin( mercuryDia/2. 162 // std::sqrt(std::pow(z0,2)+std::pow(mer 163 164 // // on the y axis, the light emitted from 165 // G4double phi = std::asin( mercuryDia/2.*s 166 // std::sqrt( std::pow(mercuryDia/2. 167 // std::pow(z0 - mercuryDia/2.*std:: 168 // - theta; 169 170 // // on the x axis, the light emitted from 171 // G4double zeta = std::atan( a/std::sqrt(st 172 173 174 175 // //alpha in [-zeta;zeta] 176 // G4double alpha = (2*zeta)*G4UniformRand() 177 // //beta in [theta-phi;theta+phi] 178 // G4double beta = (G4UniformRand()*2*phi) - 179 180 // G4double dirY = std::sin(beta); 181 // G4double dirX = std::sin(alpha); 182 183 // particleGun->SetParticleMomentumDirection 184 185 // particleGun->SetParticlePosition(G4ThreeV 186 187 // } 188 189 190 191 //shoot particles according to a certain spe 192 if (spectrum =="on") 193 { 194 G4String particle = particleGun->GetPar 195 ->GetParticleName(); 196 if(particle == "proton"|| particle == "a 197 { 198 G4DataVector* energies = runManager->GetE 199 G4DataVector* data = runManager->GetData( 200 201 G4double sum = runManager->GetDataSum(); 202 G4double partSum = 0; 203 G4int j = 0; 204 G4double random= sum*G4UniformRand(); 205 while (partSum<random) 206 { 207 partSum += (*data)[j]; 208 j++; 209 } 210 211 particleGun->SetParticleEnergy((*energies) 212 213 } 214 else if (particle == "gamma") 215 { 216 const XrayFluoDataSet* dataSet = runManage 217 218 G4int i = 0; 219 G4int id = 0; 220 G4double minEnergy = 0. * keV; 221 G4double particleEnergy= 0.; 222 G4double maxEnergy = 10. * keV; 223 G4double energyRange = maxEnergy - minEner 224 225 while ( i == 0) 226 { 227 G4double random = G4UniformRand(); 228 229 G4double randomNum = G4UniformRand(); 230 231 particleEnergy = (random*energyRange) 232 233 if ((dataSet->FindValue(particleEnergy 234 { 235 i = 1; 236 237 } 238 } 239 particleGun->SetParticleEnergy(particleEn 240 } 241 } 242 243 244 #ifdef G4ANALYSIS_USE 245 246 G4double partEnergy = particleGun->GetPartic 247 XrayFluoAnalysisManager* analysis = XrayFlu 248 analysis->analysePrimaryGenerator(partEnergy 249 250 #endif 251 252 particleGun->GeneratePrimaryVertex(anEvent); 253 } 254 255 //....oooOO0OOooo........oooOO0OOooo........oo 256 257 258 259 260 261 262 263 264 265