Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/examples/advanced/composite_calorimeter/src/CCalPrimaryGeneratorAction.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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 // File: CCalPrimaryGeneratorAction.cc
 28 // Description: CCalPrimaryGeneratorAction Sets up particle beam
 29 ///////////////////////////////////////////////////////////////////////////////
 30 
 31 #include <CLHEP/Random/RandFlat.h>
 32 
 33 #include "CCalPrimaryGeneratorAction.hh"
 34 #include "CCalPrimaryGeneratorMessenger.hh"
 35 #include "G4HEPEvtInterface.hh"
 36 
 37 #include "G4Exp.hh"
 38 #include "G4PhysicalConstants.hh"
 39 #include "G4SystemOfUnits.hh"
 40 #include "G4Event.hh"
 41 #include "G4ParticleGun.hh"
 42 #include "G4ParticleTable.hh"
 43 #include "G4ParticleDefinition.hh"
 44 #include "G4HEPEvtInterface.hh"
 45 #include "G4RunManager.hh"
 46 
 47 //#define debug
 48 
 49 CCalPrimaryGeneratorAction::CCalPrimaryGeneratorAction(): particleGun(0),
 50   generatorInput(singleFixed),  verboseLevel(0), n_particle(1), 
 51   particleName("pi-"), particleEnergy(100*GeV), particlePosition(0.,0.,0.),
 52   particleDir(1.,1.,0.1), isInitialized(0), scanSteps(0) {
 53   
 54   //Initialise the messenger
 55   gunMessenger = new CCalPrimaryGeneratorMessenger(this);
 56     
 57   //Default settings:
 58   SetMinimumEnergy(1.*GeV);
 59   SetMaximumEnergy(1.*TeV);
 60   SetMinimumEta(0.);
 61   SetMaximumEta(3.5);
 62   SetMinimumPhi(0.*degree);
 63   SetMaximumPhi(360.*degree);
 64   SetStepsPhi(1);
 65   SetStepsEta(1);
 66     
 67   // Number of particles per gun
 68   particleGun = new G4ParticleGun(n_particle);
 69     
 70   // Getting particle definition
 71   G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
 72   G4ParticleDefinition* particle = particleTable->FindParticle(particleName);
 73   particleGun->SetParticleDefinition(particle);
 74     
 75   // Setting particle properties
 76   particleGun->SetParticleEnergy(particleEnergy);
 77   particleGun->SetParticleMomentumDirection(particleDir);
 78   particleGun->SetParticlePosition(particlePosition);
 79   print(0);
 80 }
 81 
 82 CCalPrimaryGeneratorAction::~CCalPrimaryGeneratorAction() {
 83   if (gunMessenger)
 84     delete gunMessenger;
 85   if (particleGun)
 86     delete particleGun;
 87 }
 88 
 89 void CCalPrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) {
 90 
 91   if (isInitialized == 0) initialize();
 92 
 93   if (generatorInput == singleRandom) {
 94     particleEnergy = CLHEP::RandFlat::shoot(energyMin,energyMax);
 95     particleGun->SetParticleEnergy(particleEnergy);
 96 
 97     G4double eta = CLHEP::RandFlat::shoot(etaMin,etaMax);
 98     G4double phi = CLHEP::RandFlat::shoot(phiMin,phiMax);
 99     G4double theta = std::atan(G4Exp(-eta))*2.;
100     G4double randomX = std::sin(theta)*std::cos(phi);
101     G4double randomY = std::sin(theta)*std::sin(phi);
102     G4double randomZ = std::cos(theta);
103   
104     particleDir = G4ThreeVector(randomX,randomY,randomZ);
105     particleGun->SetParticleMomentumDirection(particleDir);
106     if (verboseLevel >= 2 ) {
107       G4cout << "Energy " << particleEnergy/GeV << " GeV; Theta " 
108              << theta/deg << " degree; Phi " << phi/deg << " degree" << G4endl;
109       G4cout << "Shooting in " << particleDir << " direction "<< G4endl;
110     }
111   } else if (generatorInput == singleScan) {
112     G4double scanEtaStep, scanPhiStep;
113     if (scanSteps == 0) { 
114       scanPhiStep = (phiMax - phiMin) / phiSteps;
115       phiValue = phiMin - scanPhiStep; //first step is going to change scanCurrentPhiValue
116       etaValue = etaMin;
117     }
118 
119     scanEtaStep = (etaMax - etaMin) / etaSteps;
120     scanPhiStep = (phiMax - phiMin) / phiSteps;
121 #ifdef debug
122     if (verboseLevel > 2 ) {
123       G4cout << " scanEtaStep " << scanEtaStep << " # of Steps " << etaSteps 
124              << G4endl;
125       G4cout << " scanPhiStep " << scanPhiStep << " # of Steps " << phiSteps 
126              << G4endl;
127     }
128 #endif
129 
130     //----- First scan in phi, then in eta
131     if (phiMax - phiValue < 1.E-6 * scanPhiStep) { // !only <= 1.E6 steps allowed
132       if (etaMax - etaValue < 1.E-6 * scanEtaStep) { // !only <= 1.E6 steps allowed
133         G4cout << " Scan completed!" << G4endl;
134         return;
135       } else {
136         etaValue += scanEtaStep; 
137         phiValue  = phiMin;
138       }
139     } else {
140       phiValue += scanPhiStep;
141     }    
142     G4double theta = std::atan(G4Exp(-etaValue))*2.;
143 
144     G4double scanX = std::sin(theta)*std::cos(phiValue);
145     G4double scanY = std::sin(theta)*std::sin(phiValue);
146     G4double scanZ = std::cos(theta);
147     if (verboseLevel >= 2 ) {
148       G4cout << "Scan eta " << etaValue << " Phi " << phiValue/deg
149              << " theta " << theta/deg << G4endl;
150     }
151     particleDir = G4ThreeVector(scanX,scanY,scanZ);
152     particleGun->SetParticleMomentumDirection(particleDir);
153 #ifdef debug
154     if (verboseLevel > 2 ) {
155       G4cout  << "Shooting in " << particleDir << " direction "<< G4endl;
156     }
157 #endif
158     scanSteps++;
159   }
160   
161   // Generate GEANT4 primary vertex
162   particleGun->GeneratePrimaryVertex(anEvent);
163 } 
164 
165 
166 void CCalPrimaryGeneratorAction::SetVerboseLevel(G4int val){
167   verboseLevel = val;
168 }
169 
170 
171 void CCalPrimaryGeneratorAction::SetRandom(G4String val) { 
172 
173   if (val=="on") {
174     generatorInput = singleRandom;
175     print (1);
176   } else {
177     generatorInput = singleFixed;
178     print (1);
179   }
180 }
181 
182 
183 void CCalPrimaryGeneratorAction::SetScan(G4String val) { 
184 
185   if (val=="on") {
186     generatorInput = singleScan;
187     scanSteps = 0;
188     print (1);
189   } else {
190     generatorInput = singleFixed;
191     print (1);
192   }  
193 }
194 
195 
196 void CCalPrimaryGeneratorAction::SetMinimumEnergy(G4double p){
197 
198   if (p <= 0.) {
199     G4cerr << "CCalPrimaryGeneratorAction::SetMinimumEnergy: value " << p/GeV 
200            << "GeV is out of bounds, it will not be used" << G4endl;
201     G4cerr << " Should be  >0.  Please check" << G4endl; 
202   } else {
203     energyMin = p;
204 #ifdef debug
205     if (verboseLevel >= 1 ) {
206       G4cout << " CCalPrimaryGeneratorAction: setting min. value of energy to "
207              << energyMin/GeV << " GeV " << G4endl;
208     }
209 #endif
210   }
211 }
212 
213 
214 void CCalPrimaryGeneratorAction::SetMaximumEnergy(G4double p){
215 
216   if (p <= 0.) {
217     G4cerr << "CCalPrimaryGeneratorAction::SetMaximumEnergy: value " << p/GeV 
218            << "GeV is out of bounds, it will not be used" << G4endl;
219     G4cerr << " Should be  >0.  Please check" << G4endl; 
220   } else {
221     energyMax = p;
222 #ifdef debug
223     if (verboseLevel >= 1 ) {
224       G4cout << " CCalPrimaryGeneratorAction: setting max. value of energy to "
225              << energyMax/GeV << " GeV " << G4endl;
226     }
227 #endif
228   }
229 }
230 
231 
232 void CCalPrimaryGeneratorAction::SetMinimumPhi(G4double p){
233 
234   if (std::fabs(p)>2.*pi) {
235     G4cerr << "CCalPrimaryGeneratorAction::SetMinimumPhi: setting value quite "
236            << "large" << G4endl;
237     G4cerr << " Should be given in radians - Please check" << G4endl;
238   } else {
239     phiMin = std::fabs(p);
240 #ifdef debug
241     if (verboseLevel >= 1 ) {
242       G4cout << " CCalPrimaryGeneratorAction: setting min. value of phi to "
243              << phiMin << G4endl;
244     }
245 #endif
246   }
247 }
248 
249 
250 void CCalPrimaryGeneratorAction::SetMaximumPhi(G4double p){
251 
252   if (std::fabs(p)>2.*pi) {
253     G4cerr << "CCalPrimaryGeneratorAction::SetMaximumPhi: setting value quite "
254            << "large" << G4endl;
255     G4cerr << " Should be given in radians - Please check" << G4endl;
256   } else {
257     phiMax = std::fabs(p);
258 #ifdef debug
259     if (verboseLevel >= 1 ) {
260       G4cout << " CCalPrimaryGeneratorAction: setting max. value of phi to "
261              << phiMax << G4endl;
262     }
263 #endif
264   }
265 }
266 
267 
268 void CCalPrimaryGeneratorAction::SetStepsPhi(G4int val){
269 
270   if (val <= 0) {
271     G4cerr << "CCalPrimaryGeneratorAction::SetStepsPhi: value " << val 
272            << " is out of bounds, it will not be used" << G4endl;
273     G4cerr << " Should be  > 0  Please check" << G4endl; 
274   } else {
275     phiSteps = val;
276 #ifdef debug
277     if (verboseLevel >= 1 ) {
278       G4cout << " CCalPrimaryGeneratorAction: setting no. of steps in phi to "
279              << phiSteps << G4endl;
280     }
281 #endif
282   }
283 }
284 
285 
286 void CCalPrimaryGeneratorAction::SetMinimumEta(G4double p){
287 
288   etaMin = p;
289 #ifdef debug
290   if (verboseLevel >= 1 ) {
291     G4cout << " CCalPrimaryGeneratorAction: setting min. value of eta to "
292            << etaMin << G4endl;
293   }
294 #endif
295 }
296 
297 
298 void CCalPrimaryGeneratorAction::SetMaximumEta(G4double p){
299 
300   etaMax = p;
301 #ifdef debug
302   if (verboseLevel >= 1 ) {
303     G4cout << " CCalPrimaryGeneratorAction: setting max. value of eta to "
304            << etaMax << G4endl;
305   }
306 #endif
307 }
308 
309 
310 void CCalPrimaryGeneratorAction::SetStepsEta(G4int val){
311 
312   if (val <= 0) {
313     G4cerr<<"CCalPrimaryGeneratorAction::SetStepsEta: value " << val << " is out of bounds, it will not be used"<<G4endl;
314     G4cerr<<" Should be  > 0  Please check"<<G4endl; 
315   } else {
316     etaSteps = val;
317 #ifdef debug
318     if (verboseLevel >= 1 ) {
319       G4cout << " CCalPrimaryGeneratorAction: setting no. of steps in eta to "
320              << etaSteps << G4endl;
321     }
322 #endif
323   }
324 }
325 
326 void CCalPrimaryGeneratorAction::SetGunPosition(const G4ThreeVector & pos) const {
327 
328   particleGun->SetParticlePosition(pos);
329 }
330 
331 void CCalPrimaryGeneratorAction::SetRunNo(G4int val){
332   G4RunManager::GetRunManager()->SetRunIDCounter( val );
333 }
334 
335 void CCalPrimaryGeneratorAction::initialize(){
336 
337   isInitialized = 1;
338 
339   print (1);
340 }
341 
342 
343 void CCalPrimaryGeneratorAction::print(G4int val){
344 
345   if (verboseLevel >= val) {
346 
347     if (generatorInput == singleRandom) {
348       G4cout << G4endl
349              << "**********************************************************************" << G4endl
350              << "*                                                                    *" << G4endl  
351              << "* CCalPrimaryGeneratorAction DEFAULT Random Energy/Direction setting:*" << G4endl
352              << "*                                                                    *" << G4endl  
353              << "*                                                                    *" << G4endl  
354              << "*   Energy in    [ "<< energyMin/GeV   << " - " << energyMax/GeV   << "] (GeV) "<< G4endl
355              << "*   Phi angle in [ "<< phiMin          << " - " << phiMax          << "] (rad) "<< G4endl
356              << "*                [ "<< phiMin/degree   << " - " << phiMax/degree   << "] (deg) "<< G4endl 
357              << "*   Eta in       [ "<< etaMin          << " - " << etaMax          << "]       "<< G4endl
358              << "*                                                                    *" << G4endl  
359              << "*                                                                    *" << G4endl  
360                    << "**********************************************************************" << G4endl;
361     } else if (generatorInput == singleScan) {
362       G4cout << G4endl
363              << "**********************************************************************" << G4endl
364              << "*                                                                    *" << G4endl  
365              << "* CCalPrimaryGeneratorAction DEFAULT Scan Direction settings :       *" << G4endl
366              << "*                                                                    *" << G4endl  
367              << "*                                                                    *" << G4endl  
368              << "*   Phi angle in [ " << phiMin/degree   << " - " << phiMax/degree << "] (deg) " << G4endl
369              << "*   Eta in       [ " << etaMin          << " - " << etaMax        << "]       " << G4endl
370              << "*   Steps along eta " << etaSteps << " and along phi " << phiSteps << G4endl
371              << "*                                                                    *" << G4endl  
372              << "*                                                                    *" << G4endl  
373                    << "**********************************************************************" << G4endl;
374     } else  if (generatorInput == singleFixed) {
375       G4cout << G4endl
376              << "*******************************************************************" << G4endl
377              << "*                                                                 *" << G4endl  
378              << "* CCalPrimaryGeneratorAction: Current settings :                  *" << G4endl
379              << "*                                                                 *" << G4endl  
380              << "* " << particleGun->GetNumberOfParticles() 
381              << "  " << particleGun->GetParticleDefinition()->GetParticleName() 
382              << " of " << particleGun->GetParticleEnergy()/GeV << " GeV" << G4endl
383              << "* will be shot from " << particleGun->GetParticlePosition() << G4endl;
384       G4cout << "* in direction " << particleGun->GetParticleMomentumDirection() << G4endl;
385       G4cout << "*                                                                 *" << G4endl  
386              << "*                                                                 *" << G4endl  
387              << "*******************************************************************" << G4endl;
388     }
389   } 
390 }
391