]>
Commit | Line | Data |
---|---|---|
67bfd92b | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | ||
17 | // Generator for particle pairs in a preset | |
18 | // kinematic range | |
19 | // ranges can be set for invariant mass, pair pT, pair rapidity | |
20 | // and pair azimuth | |
21 | ||
22 | // | |
23 | // Comments and suggestions: markus.konrad.kohler@cern.ch | |
24 | // | |
25 | ||
26 | #include "TPDGCode.h" | |
27 | #include <TDatabasePDG.h> | |
28 | #include "AliGenEventHeader.h" | |
29 | #include "TF1.h" | |
30 | ||
31 | #include "AliGenPairFlat.h" | |
32 | ||
33 | ClassImp(AliGenPairFlat) | |
34 | ||
35 | //_____________________________________________________________________________ | |
36 | AliGenPairFlat::AliGenPairFlat() | |
37 | :AliGenerator(), | |
38 | fPairNpart(10), | |
39 | fPairYMin(-10.), | |
40 | fPairYMax(10.), | |
41 | fPairPhiMin(0.), | |
42 | fPairPhiMax(TMath::TwoPi()), | |
43 | fPairPtMin(0.001), | |
44 | fPairPtMax(50.), | |
45 | fPairMassMin(0.), | |
46 | fPairMassMax(10.), | |
47 | fLegPdg1(11), | |
48 | fLegPdg2(-11), | |
49 | fAlpha(0.), | |
8c7c80ad | 50 | fDebug(0), |
51 | fPol(0) | |
67bfd92b | 52 | { |
53 | // | |
54 | // Default constructor | |
55 | // | |
56 | } | |
57 | ||
8c7c80ad | 58 | |
59 | //_____________________________________________________________________________ | |
60 | AliGenPairFlat::~AliGenPairFlat() | |
61 | { | |
62 | // | |
63 | // Destructor | |
64 | // | |
65 | delete fPol; | |
66 | ||
67 | } | |
68 | ||
67bfd92b | 69 | //_____________________________________________________________________________ |
70 | void AliGenPairFlat::Generate() | |
71 | { | |
72 | // | |
73 | // Generate a random pair of particles | |
74 | // | |
75 | ||
76 | Float_t polar[3]= {0,0,0}; | |
77 | Float_t origin[3]; | |
78 | Float_t time; | |
79 | ||
80 | Float_t p[3]; | |
81 | Int_t i, j, nt; | |
82 | Double_t phi, y, mass, mt, e, pt; | |
83 | Float_t weight = 1.; | |
84 | Double_t pt1, pt2, y1, y2; | |
85 | ||
86 | TLorentzVector mother, dau1, dau2; | |
87 | Float_t random[6]; | |
88 | ||
8c7c80ad | 89 | fPol = new TF1("fPol","1.+[0]*x*x",-1.,1.); |
90 | fPol->SetParameter(0, fAlpha); | |
67bfd92b | 91 | |
92 | ||
93 | for (j=0;j<3;j++) origin[j]=fOrigin[j]; | |
94 | time = fTimeOrigin; | |
95 | if(fVertexSmear==kPerEvent) { | |
96 | Vertex(); | |
97 | for (j=0;j<3;j++) origin[j]=fVertex[j]; | |
98 | time = fTime; | |
99 | } | |
100 | ||
8c7c80ad | 101 | if(fDebug == 2){ |
67bfd92b | 102 | printf("\n\n------------------GENERATOR SETTINGS------------------\n\n"); |
103 | printf("You choosed for the mother the Mass range %f - %f \n",fPairMassMin,fPairMassMax); | |
104 | printf("You choosed for the mother the transverse Momentum range %f - %f \n",fPairPtMin,fPairPtMax); | |
105 | printf("You choosed for the mother the Phi range %f - %f \n",fPairPhiMin,fPairPhiMax); | |
106 | printf("The Particle will decay in (pdg) %i and %i\n \n",fLegPdg1, fLegPdg2); | |
107 | printf("The rest Mass of first daughter (%s) is %f \n", | |
108 | TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg1))->GetTitle(), | |
109 | TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg1))->Mass()); | |
110 | printf("The rest Mass of second daughter (%s) is %f \n", | |
111 | TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg2))->GetTitle(), | |
112 | TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg2))->Mass()); | |
113 | printf("polarization factor is alpha == %lf \n",fAlpha); | |
114 | printf("vertex is at x == %f || y == %f || z == %f \n",origin[0],origin[1],origin[2]); | |
115 | printf("\n----------------------------------------------------------\n"); | |
8c7c80ad | 116 | } |
67bfd92b | 117 | |
118 | for(i=0;i<fPairNpart;i++) { | |
119 | ||
120 | // mother properties | |
121 | Rndm(random,4); | |
122 | mass = fPairMassMin+random[0]*(fPairMassMax-fPairMassMin); | |
123 | pt = fPairPtMin+random[1]*(fPairPtMax-fPairPtMin); | |
124 | y = fPairYMin+random[2]*(fPairYMax-fPairYMin); | |
125 | phi = fPairPhiMin+random[3]*(fPairPhiMax-fPairPhiMin); | |
126 | ||
127 | mt = TMath::Sqrt(pt*pt + mass*mass); | |
128 | p[0] = pt*TMath::Cos(phi); | |
129 | p[1] = pt*TMath::Sin(phi); | |
130 | p[2] = mt*TMath::SinH(y); | |
131 | e = mt*TMath::CosH(y); | |
132 | ||
133 | mother.SetPxPyPzE(p[0],p[1],p[2],e); | |
134 | ||
135 | ||
136 | if (fDebug == 2) printf("p = (%+11.4e,%+11.4e,%+11.4e) GeV || pt = %+11.4e || y = %+11.4e || weight=%+11.4e || E= %+11.4e\n",p[0],p[1],p[2],pt, y,weight,e); | |
137 | ||
138 | //decay procedure | |
8c7c80ad | 139 | if(!Decay(mother,dau1,dau2,fPol))continue; |
67bfd92b | 140 | |
141 | pt1 = dau1.Pt(); | |
142 | pt2 = dau2.Pt(); | |
143 | y1 = dau1.Rapidity(); | |
144 | y2 = dau2.Rapidity(); | |
145 | ||
146 | //first daughter | |
147 | p[0] = dau1.Px(); | |
148 | p[1] = dau1.Py(); | |
149 | p[2] = dau1.Pz(); | |
150 | e = dau1.E(); | |
151 | ||
152 | if(fVertexSmear==kPerTrack) { | |
153 | Rndm(random,6); | |
154 | for (j=0;j<3;j++) { | |
155 | origin[j]=fOrigin[j]+fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())* | |
156 | TMath::Sqrt(-2*TMath::Log(random[2*j+1])); | |
157 | } | |
158 | ||
159 | Rndm(random,2); | |
160 | time = fTimeOrigin + fOsigma[2]/TMath::Ccgs()* | |
161 | TMath::Cos(2*random[0]*TMath::Pi())* | |
162 | TMath::Sqrt(-2*TMath::Log(random[1])); | |
163 | } | |
164 | ||
165 | PushTrack(fTrackIt,-1,fLegPdg1,p,origin,polar,time,kPPrimary,nt, weight); | |
166 | if (fDebug == 2) printf("dau1-->id=%+3d, p = (%+11.4e,%+11.4e,%+11.4e) GeV || pt = %+11.4e || weight=%+11.4e || E= %+11.4e\n",fLegPdg1,p[0],p[1],p[2],pt1,weight,e); | |
167 | ||
168 | //second daughter | |
169 | p[0] = dau2.Px(); | |
170 | p[1] = dau2.Py(); | |
171 | p[2] = dau2.Pz(); | |
172 | e = dau2.E(); | |
173 | ||
174 | if(fVertexSmear==kPerTrack) { | |
175 | Rndm(random,6); | |
176 | for (j=0;j<3;j++) { | |
177 | origin[j]=fOrigin[j]+fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())* | |
178 | TMath::Sqrt(-2*TMath::Log(random[2*j+1])); | |
179 | } | |
180 | ||
181 | Rndm(random,2); | |
182 | time = fTimeOrigin + fOsigma[2]/TMath::Ccgs()* | |
183 | TMath::Cos(2*random[0]*TMath::Pi())* | |
184 | TMath::Sqrt(-2*TMath::Log(random[1])); | |
185 | } | |
186 | ||
187 | PushTrack(fTrackIt,-1,fLegPdg2,p,origin,polar,time,kPPrimary,nt,weight); | |
188 | if (fDebug == 2) printf("dau2-->id=%+3d, p = (%+11.4e,%+11.4e,%+11.4e) GeV || pt = %+11.4e || weight=%+11.4e || E= %+11.4e\n",fLegPdg2,p[0],p[1],p[2],pt2, weight,e); | |
189 | ||
190 | }//particle loop | |
191 | ||
192 | } | |
193 | ||
194 | //_____________________________________________________________________________ | |
67bfd92b | 195 | void AliGenPairFlat::Init() |
196 | { | |
197 | // | |
198 | // Init | |
199 | // | |
200 | ||
201 | printf("AliGenPairFlat::Init() not implemented!!!\n"); | |
202 | ||
203 | } | |
204 | ||
205 | ||
206 | //_____________________________________________________________________________ | |
8c7c80ad | 207 | Bool_t AliGenPairFlat::Decay(const TLorentzVector& mother, TLorentzVector &dau1, TLorentzVector &dau2, TF1* polfactor) |
67bfd92b | 208 | { |
209 | // | |
210 | // decay procedure | |
211 | // | |
212 | ||
213 | Double_t mp, md1, md2, ed1, ed2, pd1, pd2; | |
214 | Double_t costheta, sintheta, phi; | |
215 | ||
216 | mp = mother.M(); | |
217 | md1 = TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg1))->Mass(); | |
218 | md2 = TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg2))->Mass(); | |
219 | ||
220 | if(mp < md1+md2){ | |
221 | printf("Daughters are heavier than Mother!! Check Kinematics!! \n"); | |
222 | return 0; | |
223 | } | |
224 | ||
225 | ed1 = (mp*mp + md1*md1 - md2*md2)/(2.*mp); | |
226 | ed2 = mp-ed1; | |
227 | pd1 = TMath::Sqrt((ed1+md1)*(ed1-md1)); | |
228 | pd2 = TMath::Sqrt((mp*mp-(md1+md2)*(md1+md2))*(mp*mp -(md1-md2)*(md1-md2)))/(2.*mp); | |
229 | ||
725547b6 | 230 | costheta = polfactor->GetRandom(); //polarization |
67bfd92b | 231 | sintheta = TMath::Sqrt((1.+costheta)*(1.-costheta)); |
232 | phi = 2.0*TMath::Pi()*gRandom->Rndm(); | |
233 | ||
234 | dau1.SetPxPyPzE(pd1*sintheta*TMath::Cos(phi), | |
235 | pd1*sintheta*TMath::Sin(phi), | |
236 | pd1*costheta, | |
237 | ed1); | |
238 | ||
239 | dau2.SetPxPyPzE((-1.)*pd1*sintheta*TMath::Cos(phi), | |
240 | (-1.)*pd1*sintheta*TMath::Sin(phi), | |
241 | (-1.)*pd1*costheta, | |
242 | ed2); | |
243 | ||
244 | TVector3 boost = mother.BoostVector(); | |
245 | ||
246 | dau1.Boost(boost); | |
247 | dau2.Boost(boost); | |
248 | ||
249 | return 1; | |
250 | ||
251 | } | |
252 | ||
253 | ||
254 |