/[PAMELA software]/PamVMC/src/PamVMCPrimaryGenerator.cxx
ViewVC logotype

Contents of /PamVMC/src/PamVMCPrimaryGenerator.cxx

Parent Directory Parent Directory | Revision Log Revision Log

Revision 1.5 - (show annotations) (download)
Fri Jun 12 18:39:32 2009 UTC (16 years, 4 months ago) by pam-rm2
Branch: MAIN
CVS Tags: v1r0
Changes since 1.2: +17 -19 lines 
- Introduced user-defined names of output files and random seeds number.
Users can do it use options of PamVMCApplication constructor:
PamVMCApplication(const char* name,  const char *title, const char*
filename="pamtest", Int_t seed=0).
The Random object that I use is TRandom3 object which has astronomical
large period (in case of default initialization 0). All random generators
in the code use this object by calling of gRandom singleton which keeps
it.

- Corrected TOF digitization routine. No problems with TDC hits due to
hadronic interactions anymore.

- Some small changes was done to compile code under Root 5.23. +
geant4_vmc v. 2.6 without any warnings

- Some classes of PamG4RunConfiguartion was changed for geant4_vmc v.
2.6.Some obsolete classes was deleted as soon as developers implemented
regions.

- Navigation was changed from "geomRootToGeant4" to "geomRoot", because on
VMC web page written that as soon as Geant4 has no option ONLY/MANY
translation of overlapped geometry to Geant4 through VGM could be wrong.
I'd like to stay with Root navigation:
http://root.cern.ch/root/vmc/Geant4VMC.html. This should be default
option.

- New Tracker digitization routine written by Sergio was implemented

- PamVMC again became compatible with geant4_vmc v.2.5 and ROOT 5.20.
 The problem was that ROOT developers introduced in TVirtualMC class a new
method SetMagField and new base class:TVirtualMagField from which
user-defined classes shoukd be derived
1 // $Id: PamVMCPrimaryGenerator.cxx,v 1.0 2006/06/03
2
3
4 #include <TVirtualMC.h>
5 #include <TVirtualMCStack.h>
6 #include <TPDGCode.h>
7 #include <TDatabasePDG.h>
8 #include <TParticlePDG.h>
9 #include <TVector3.h>
10 #include <TMath.h>
11 #include <Riostream.h>
12
13 #include "PamVMCPrimaryGenerator.h"
14
15 using namespace TMath;
16
17 ClassImp(PamVMCPrimary)
18
19 PamVMCPrimary & operator+=(PamVMCPrimary &a, const PamVMCPrimary &b)
20 {
21 a.fPDG=b.fPDG;
22 a.fX0=b.fX0;
23 a.fY0=b.fY0;
24 a.fZ0=b.fZ0;
25 a.fTHETA=b.fTHETA;
26 a.fPHI=b.fPHI;
27 a.fP0=b.fP0;
28 a.fGOOD=b.fGOOD;
29
30 return a;
31 }
32
33
34 ClassImp(PamVMCPrimaryGenerator)
35
36 PamVMCPrimaryGenerator::PamVMCPrimaryGenerator(TVirtualMCStack* stack)
37 : TObject(),
38 fStack(stack),
39 fevno(0),
40 fmass(0.),
41 fcharge(0.),
42 frandom(0)
43 {
44 // Standard constructor
45
46 ftheta = new TF1("ftheta","sin(x)*cos(x)",0.,acos(-1.)/4.);
47 ftheta->SetNpx(1000);
48
49 fprimColl = new TClonesArray("PamVMCPrimary");
50 fprim.fPDG=kProton;
51 fprim.fX0=1.;
52 fprim.fY0=1.;
53 fprim.fZ0=130.;
54 fprim.fTHETA=0.;
55 fprim.fPHI=0.;
56 fprim.fP0=1.; //1GV
57
58 }
59
60 PamVMCPrimaryGenerator::PamVMCPrimaryGenerator()
61 : TObject(),
62 fStack(0),
63 fevno(0),
64 fmass(0.),
65 fcharge(0.),
66 fprimColl(0),
67 frandom(0)
68 {
69 // Default constructor
70 //Default primary proton
71 ftheta = new TF1("ftheta","sin(x)*cos(x)",0.,acos(-1.)/4.);
72 ftheta->SetNpx(1000);
73
74 fprim.fPDG=kProton;
75 fprim.fX0=1.;
76 fprim.fY0=1.;
77 fprim.fZ0=130.;
78 fprim.fTHETA=0.;
79 fprim.fPHI=0.;
80 fprim.fP0=1.; //1GV
81 }
82
83 PamVMCPrimaryGenerator::~PamVMCPrimaryGenerator()
84 {
85 // Destructor
86 delete ftheta;
87 delete fprimColl;
88 }
89
90 // private methods
91
92
93 void PamVMCPrimaryGenerator::GeneratePrimary()
94 {
95 // Add one primary particle to the user stack (derived from TVirtualMCStack).
96
97 // Track ID (filled by stack)
98 Int_t ntr;
99
100 // Option: to be tracked
101 Int_t toBeDone = 1;
102
103 // Particle type
104 Int_t pdg = fprim.fPDG;
105
106 Double_t fvx, fvy, fvz;
107 fvx=fprim.fX0;
108 fvy=fprim.fY0;
109 fvz=fprim.fZ0;
110
111 // Position
112
113 Double_t tof = 0.;
114
115 // Energy (in GeV)
116 Double_t kinEnergy = MomentumToKinE(fprim.fP0);
117 Double_t e = fmass + kinEnergy;
118
119 // Particle momentum
120 Double_t px, py, pz;
121
122 px = fprim.fP0*Sin(fprim.fTHETA)*Cos(fprim.fPHI);
123 py = fprim.fP0*Sin(fprim.fTHETA)*Sin(fprim.fPHI);
124 pz = -fprim.fP0*Cos(fprim.fTHETA);
125
126 // Polarization
127 TVector3 polar;
128
129 // Add particle to stack
130 fStack->PushTrack(toBeDone, -1, pdg, px, py, pz, e, fvx, fvy, fvz, tof,
131 polar.X(), polar.Y(), polar.Z(),
132 kPPrimary, ntr, 1., 0);
133
134 PamVMCPrimary * pc = (PamVMCPrimary *)fprimColl->New(fevno++);
135
136 *pc = fprim;
137 }
138
139
140 void PamVMCPrimaryGenerator::SetParticle(Int_t pdg){
141 fprim.fPDG=pdg;
142 //TParticlePDG* particlePDG = TDatabasePDG::Instance()->GetParticle(fprim.fPDG);
143 fmass = (TDatabasePDG::Instance()->GetParticle(fprim.fPDG))->Mass();
144 fcharge = ((TDatabasePDG::Instance()->GetParticle(fprim.fPDG))->Charge())/3.;
145 }
146
147 void PamVMCPrimaryGenerator::SetMomentum(
148 Double_t px, Double_t py, Double_t pz)
149 {
150 fprim.fP0= Sqrt(px*px+py*py+pz*pz);
151 fprim.fTHETA=ATan(Sqrt(px*px+py*py)/pz);
152 fprim.fPHI=ATan(py/px);
153 }
154
155 void PamVMCPrimaryGenerator::GenSpe(Double_t PEmin, Double_t PEmax, Bool_t isEnergy)
156 {
157 if(isEnergy) {
158 fprim.fP0=frandom->Uniform(KinEToMomentum(PEmin),KinEToMomentum(PEmax));
159 } else{
160 fprim.fP0=frandom->Uniform(PEmin,PEmax);
161 }
162
163 }
164
165 void PamVMCPrimaryGenerator::GenSpe(Double_t PEmin, Double_t PEmax, Double_t gamma, Bool_t isEnergy)
166 {
167 Double_t alpha = 1.+gamma; //integral spectral index
168 if(alpha==0.){
169 fprim.fP0=Exp(Log(PEmin)+frandom->Uniform(0.,1.)*(Log(PEmax)-Log(PEmin)));
170 } else {
171 if(PEmin==0.) PEmin=1.E-10;
172 fprim.fP0=Power((frandom->Uniform(0.,1.)*(Power(PEmax,alpha)-Power(PEmin,alpha))+Power(PEmin,alpha)),1./alpha);
173 }
174
175 if(isEnergy) fprim.fP0=KinEToMomentum(fprim.fP0);
176
177 }
  ViewVC Help
Powered by ViewVC 1.1.23