Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer << 3 // * DISCLAIMER * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th << 5 // * The following disclaimer summarizes all the specific disclaimers * 6 // * the Geant4 Collaboration. It is provided << 6 // * of contributors to this software. The specific disclaimers,which * 7 // * conditions of the Geant4 Software License << 7 // * govern, are listed with their locations in: * 8 // * LICENSE and available at http://cern.ch/ << 8 // * http://cern.ch/geant4/license * 9 // * include a list of copyright holders. << 10 // * 9 // * * 11 // * Neither the authors of this software syst 10 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 11 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 12 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 13 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file << 14 // * use. * 16 // * for the full disclaimer and the limitatio << 17 // * 15 // * * 18 // * This code implementation is the result << 16 // * This code implementation is the intellectual property of the * 19 // * technical work of the GEANT4 collaboratio << 17 // * GEANT4 collaboration. * 20 // * By using, copying, modifying or distri << 18 // * By copying, distributing or modifying the Program (or any work * 21 // * any work based on the software) you ag << 19 // * based on the Program) you indicate your acceptance of this * 22 // * use in resulting scientific publicati << 20 // * statement, and all its terms. * 23 // * acceptance of all terms of the Geant4 Sof << 24 // ******************************************* 21 // ******************************************************************** 25 // 22 // 26 // G4SingleParticleSource source implementatio << 27 // << 28 // Author: Fan Lei, QinetiQ ltd. - 05/02/2004 << 29 // Customer: ESA/ESTEC << 30 // Revision: Andrea Dotti, SLAC << 31 // ------------------------------------------- << 32 << 33 #include <cmath> << 34 23 35 #include "G4SingleParticleSource.hh" << 24 /////////////////////////////////////////////////////////////////////////////// 36 << 25 // 37 #include "G4SystemOfUnits.hh" << 26 // MODULE: G4SingleParticleSource.hh >> 27 // >> 28 // Version: 1.0 >> 29 // Date: 5/02/04 >> 30 // Author: Fan Lei >> 31 // Organisation: QinetiQ ltd. >> 32 // Customer: ESA/ESTEC >> 33 // >> 34 /////////////////////////////////////////////////////////////////////////////// >> 35 // >> 36 // CHANGE HISTORY >> 37 // -------------- >> 38 // >> 39 // Version 1.0, 05/02/2004, Fan Lei, Created. >> 40 // Based on the G4GeneralParticleSource class in Geant4 v6.0 >> 41 // >> 42 /////////////////////////////////////////////////////////////////////////////// >> 43 // 38 #include "G4PrimaryParticle.hh" 44 #include "G4PrimaryParticle.hh" 39 #include "G4Event.hh" 45 #include "G4Event.hh" 40 #include "Randomize.hh" 46 #include "Randomize.hh" >> 47 #include <cmath> 41 #include "G4ParticleTable.hh" 48 #include "G4ParticleTable.hh" 42 #include "G4Geantino.hh" 49 #include "G4Geantino.hh" 43 #include "G4ParticleDefinition.hh" 50 #include "G4ParticleDefinition.hh" 44 #include "G4IonTable.hh" 51 #include "G4IonTable.hh" 45 #include "G4Ions.hh" 52 #include "G4Ions.hh" 46 #include "G4TrackingManager.hh" 53 #include "G4TrackingManager.hh" 47 #include "G4Track.hh" 54 #include "G4Track.hh" 48 #include "G4AutoLock.hh" << 55 #include "G4SingleParticleSource.hh" 49 << 50 G4SingleParticleSource::part_prop_t::part_prop << 51 { << 52 momentum_direction = G4ParticleMomentum(1,0, << 53 energy = 1.*MeV; << 54 position = G4ThreeVector(); << 55 } << 56 56 57 G4SingleParticleSource::G4SingleParticleSource 57 G4SingleParticleSource::G4SingleParticleSource() 58 { 58 { 59 // Initialise all variables 59 // Initialise all variables 60 // Position distribution Variables 60 // Position distribution Variables 61 61 62 NumberOfParticlesToBeGenerated = 1; 62 NumberOfParticlesToBeGenerated = 1; 63 definition = G4Geantino::GeantinoDefinition( << 63 particle_definition = G4Geantino::GeantinoDefinition(); 64 << 64 G4ThreeVector zero; 65 charge = 0.0; << 65 particle_momentum_direction = G4ParticleMomentum(1,0,0); 66 time = 0; << 66 particle_energy = 1.0*MeV; 67 polarization = G4ThreeVector(); << 67 particle_position = zero; >> 68 particle_time = 0.0; >> 69 particle_polarization = zero; >> 70 particle_charge = 0.0; >> 71 particle_weight = 1.0; 68 72 69 biasRndm = new G4SPSRandomGenerator(); 73 biasRndm = new G4SPSRandomGenerator(); 70 posGenerator = new G4SPSPosDistribution(); 74 posGenerator = new G4SPSPosDistribution(); 71 posGenerator->SetBiasRndm(biasRndm); 75 posGenerator->SetBiasRndm(biasRndm); 72 angGenerator = new G4SPSAngDistribution(); 76 angGenerator = new G4SPSAngDistribution(); 73 angGenerator->SetPosDistribution(posGenerato 77 angGenerator->SetPosDistribution(posGenerator); 74 angGenerator->SetBiasRndm(biasRndm); 78 angGenerator->SetBiasRndm(biasRndm); 75 eneGenerator = new G4SPSEneDistribution(); 79 eneGenerator = new G4SPSEneDistribution(); 76 eneGenerator->SetBiasRndm(biasRndm); 80 eneGenerator->SetBiasRndm(biasRndm); 77 81 >> 82 // verbosity 78 verbosityLevel = 0; 83 verbosityLevel = 0; 79 << 84 80 G4MUTEXINIT(mutex); << 81 } 85 } 82 86 83 G4SingleParticleSource::~G4SingleParticleSourc 87 G4SingleParticleSource::~G4SingleParticleSource() 84 { << 88 {} 85 delete biasRndm; << 86 delete posGenerator; << 87 delete angGenerator; << 88 delete eneGenerator; << 89 << 90 G4MUTEXDESTROY(mutex); << 91 } << 92 89 93 void G4SingleParticleSource::SetVerbosity(G4in << 90 void G4SingleParticleSource::SetVerbosity(int vL) 94 { 91 { 95 G4AutoLock l(&mutex); << 96 verbosityLevel = vL; 92 verbosityLevel = vL; 97 posGenerator->SetVerbosity(vL); 93 posGenerator->SetVerbosity(vL); 98 angGenerator->SetVerbosity(vL); 94 angGenerator->SetVerbosity(vL); 99 eneGenerator->SetVerbosity(vL); 95 eneGenerator->SetVerbosity(vL); >> 96 G4cout << "Verbosity Set to: " << verbosityLevel << G4endl; 100 } 97 } 101 98 102 void G4SingleParticleSource:: << 99 void G4SingleParticleSource::SetParticleDefinition 103 SetParticleDefinition(G4ParticleDefinition* aP << 100 (G4ParticleDefinition* aParticleDefinition) 104 { 101 { 105 definition = aParticleDefinition; << 102 particle_definition = aParticleDefinition; 106 charge = aParticleDefinition->GetPDGCharge() << 103 particle_charge = particle_definition->GetPDGCharge(); 107 } 104 } 108 105 109 void G4SingleParticleSource::GeneratePrimaryVe << 106 void G4SingleParticleSource::GeneratePrimaryVertex(G4Event *evt) 110 { 107 { 111 if (definition == nullptr) << 108 if(particle_definition==NULL) return; 112 { << 113 // TODO: Should this rise an exception??? << 114 return; << 115 } << 116 << 117 if (verbosityLevel > 1) << 118 { << 119 G4cout << " NumberOfParticlesToBeGenerated << 120 << NumberOfParticlesToBeGenerated < << 121 } << 122 109 123 part_prop_t& pp = ParticleProperties.Get(); << 110 if(verbosityLevel > 1) >> 111 G4cout << " NumberOfParticlesToBeGenerated: "<<NumberOfParticlesToBeGenerated << G4endl; 124 112 125 // Position stuff 113 // Position stuff 126 pp.position = posGenerator->GenerateOne(); << 114 particle_position = posGenerator->GenerateOne(); 127 115 128 // Create a new vertex << 116 // create a new vertex 129 auto* vertex = new G4PrimaryVertex(pp.positi << 117 G4PrimaryVertex* vertex = new G4PrimaryVertex(particle_position,particle_time); 130 118 131 for (G4int i=0; i<NumberOfParticlesToBeGener << 119 for( G4int i=0; i<NumberOfParticlesToBeGenerated; i++ ) { 132 { << 133 // Angular stuff 120 // Angular stuff 134 pp.momentum_direction = angGenerator->Gene << 121 particle_momentum_direction = angGenerator->GenerateOne(); 135 << 136 // Energy stuff 122 // Energy stuff 137 pp.energy = eneGenerator->GenerateOne(defi << 123 particle_energy = eneGenerator->GenerateOne(particle_definition); 138 << 124 139 if (verbosityLevel >= 2) << 125 if(verbosityLevel >= 2) 140 { << 141 G4cout << "Creating primaries and assign 126 G4cout << "Creating primaries and assigning to vertex" << G4endl; >> 127 // create new primaries and set them to the vertex >> 128 G4double mass = particle_definition->GetPDGMass(); >> 129 G4double energy = particle_energy + mass; >> 130 G4double pmom = std::sqrt(energy*energy-mass*mass); >> 131 G4double px = pmom*particle_momentum_direction.x(); >> 132 G4double py = pmom*particle_momentum_direction.y(); >> 133 G4double pz = pmom*particle_momentum_direction.z(); >> 134 >> 135 if(verbosityLevel > 1){ >> 136 G4cout << "Particle name: "<<particle_definition->GetParticleName() << G4endl; >> 137 G4cout << " Energy: "<<particle_energy << G4endl; >> 138 G4cout << " Position: "<<particle_position<< G4endl; >> 139 G4cout << " Direction: "<<particle_momentum_direction << G4endl; 142 } 140 } 143 << 141 G4PrimaryParticle* particle = 144 // Create new primaries and set them to th << 142 new G4PrimaryParticle(particle_definition,px,py,pz); 145 // << 146 G4double mass = definition->GetPDGMass(); << 147 auto* particle = new G4PrimaryParticle(def << 148 particle->SetKineticEnergy(pp.energy ); << 149 particle->SetMass( mass ); 143 particle->SetMass( mass ); 150 particle->SetMomentumDirection( pp.momentu << 144 particle->SetCharge( particle_charge ); 151 particle->SetCharge( charge ); << 145 particle->SetPolarization(particle_polarization.x(), 152 particle->SetPolarization(polarization.x() << 146 particle_polarization.y(), 153 polarization.y() << 147 particle_polarization.z()); 154 polarization.z() << 155 if (verbosityLevel > 1) << 156 { << 157 G4cout << "Particle name: " << definitio << 158 G4cout << " Energy: " << pp.energy << 159 G4cout << " Position: " << pp.positi << 160 G4cout << " Direction: " << pp.moment << 161 } << 162 << 163 // Set bweight equal to the multiple of al 148 // Set bweight equal to the multiple of all non-zero weights 164 // << 149 particle_weight = biasRndm->GetBiasWeight(); 165 G4double weight = eneGenerator->GetWeight( << 150 // pass it to primary particle >> 151 particle->SetWeight(particle_weight); 166 152 167 if(eneGenerator->IfApplyEnergyWeight()) << 153 vertex->SetPrimary( particle ); 168 { << 154 169 weight *= eneGenerator->GetArbEneWeight( << 170 } << 171 << 172 // Pass it to primary particle << 173 // << 174 particle->SetWeight(weight); << 175 vertex->SetPrimary(particle); << 176 } 155 } 177 << 156 // now pass the weight to the primary vertex. CANNOT be used here! 178 // Now pass the weight to the primary vertex << 179 // vertex->SetWeight(particle_weight); 157 // vertex->SetWeight(particle_weight); 180 evt->AddPrimaryVertex(vertex); << 158 evt->AddPrimaryVertex( vertex ); 181 << 159 if(verbosityLevel > 1) 182 if (verbosityLevel > 1) << 160 G4cout << " Primary Vetex generated !"<< G4endl; 183 { << 184 G4cout << " Primary Vetex generated !" << << 185 } << 186 } 161 } >> 162 >> 163 >> 164 >> 165 >> 166 >> 167 >> 168 >> 169 >> 170 187 171