]>
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.), | |
50 | fDebug(0) | |
51 | { | |
52 | // | |
53 | // Default constructor | |
54 | // | |
55 | } | |
56 | ||
57 | //_____________________________________________________________________________ | |
58 | void AliGenPairFlat::Generate() | |
59 | { | |
60 | // | |
61 | // Generate a random pair of particles | |
62 | // | |
63 | ||
64 | Float_t polar[3]= {0,0,0}; | |
65 | Float_t origin[3]; | |
66 | Float_t time; | |
67 | ||
68 | Float_t p[3]; | |
69 | Int_t i, j, nt; | |
70 | Double_t phi, y, mass, mt, e, pt; | |
71 | Float_t weight = 1.; | |
72 | Double_t pt1, pt2, y1, y2; | |
73 | ||
74 | TLorentzVector mother, dau1, dau2; | |
75 | Float_t random[6]; | |
76 | ||
725547b6 | 77 | TF1* pol = new TF1("Pol","1.+[0]*x*x",-1.,1.); |
78 | pol->SetParameter(0, fAlpha); | |
67bfd92b | 79 | |
80 | ||
81 | for (j=0;j<3;j++) origin[j]=fOrigin[j]; | |
82 | time = fTimeOrigin; | |
83 | if(fVertexSmear==kPerEvent) { | |
84 | Vertex(); | |
85 | for (j=0;j<3;j++) origin[j]=fVertex[j]; | |
86 | time = fTime; | |
87 | } | |
88 | ||
89 | printf("\n\n------------------GENERATOR SETTINGS------------------\n\n"); | |
90 | printf("You choosed for the mother the Mass range %f - %f \n",fPairMassMin,fPairMassMax); | |
91 | printf("You choosed for the mother the transverse Momentum range %f - %f \n",fPairPtMin,fPairPtMax); | |
92 | printf("You choosed for the mother the Phi range %f - %f \n",fPairPhiMin,fPairPhiMax); | |
93 | printf("The Particle will decay in (pdg) %i and %i\n \n",fLegPdg1, fLegPdg2); | |
94 | printf("The rest Mass of first daughter (%s) is %f \n", | |
95 | TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg1))->GetTitle(), | |
96 | TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg1))->Mass()); | |
97 | printf("The rest Mass of second daughter (%s) is %f \n", | |
98 | TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg2))->GetTitle(), | |
99 | TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg2))->Mass()); | |
100 | printf("polarization factor is alpha == %lf \n",fAlpha); | |
101 | printf("vertex is at x == %f || y == %f || z == %f \n",origin[0],origin[1],origin[2]); | |
102 | printf("\n----------------------------------------------------------\n"); | |
103 | ||
104 | ||
105 | for(i=0;i<fPairNpart;i++) { | |
106 | ||
107 | // mother properties | |
108 | Rndm(random,4); | |
109 | mass = fPairMassMin+random[0]*(fPairMassMax-fPairMassMin); | |
110 | pt = fPairPtMin+random[1]*(fPairPtMax-fPairPtMin); | |
111 | y = fPairYMin+random[2]*(fPairYMax-fPairYMin); | |
112 | phi = fPairPhiMin+random[3]*(fPairPhiMax-fPairPhiMin); | |
113 | ||
114 | mt = TMath::Sqrt(pt*pt + mass*mass); | |
115 | p[0] = pt*TMath::Cos(phi); | |
116 | p[1] = pt*TMath::Sin(phi); | |
117 | p[2] = mt*TMath::SinH(y); | |
118 | e = mt*TMath::CosH(y); | |
119 | ||
120 | mother.SetPxPyPzE(p[0],p[1],p[2],e); | |
121 | ||
122 | ||
123 | 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); | |
124 | ||
125 | //decay procedure | |
725547b6 | 126 | if(!Decay(mother,dau1,dau2,pol))continue; |
67bfd92b | 127 | |
128 | pt1 = dau1.Pt(); | |
129 | pt2 = dau2.Pt(); | |
130 | y1 = dau1.Rapidity(); | |
131 | y2 = dau2.Rapidity(); | |
132 | ||
133 | //first daughter | |
134 | p[0] = dau1.Px(); | |
135 | p[1] = dau1.Py(); | |
136 | p[2] = dau1.Pz(); | |
137 | e = dau1.E(); | |
138 | ||
139 | if(fVertexSmear==kPerTrack) { | |
140 | Rndm(random,6); | |
141 | for (j=0;j<3;j++) { | |
142 | origin[j]=fOrigin[j]+fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())* | |
143 | TMath::Sqrt(-2*TMath::Log(random[2*j+1])); | |
144 | } | |
145 | ||
146 | Rndm(random,2); | |
147 | time = fTimeOrigin + fOsigma[2]/TMath::Ccgs()* | |
148 | TMath::Cos(2*random[0]*TMath::Pi())* | |
149 | TMath::Sqrt(-2*TMath::Log(random[1])); | |
150 | } | |
151 | ||
152 | PushTrack(fTrackIt,-1,fLegPdg1,p,origin,polar,time,kPPrimary,nt, weight); | |
153 | 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); | |
154 | ||
155 | //second daughter | |
156 | p[0] = dau2.Px(); | |
157 | p[1] = dau2.Py(); | |
158 | p[2] = dau2.Pz(); | |
159 | e = dau2.E(); | |
160 | ||
161 | if(fVertexSmear==kPerTrack) { | |
162 | Rndm(random,6); | |
163 | for (j=0;j<3;j++) { | |
164 | origin[j]=fOrigin[j]+fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())* | |
165 | TMath::Sqrt(-2*TMath::Log(random[2*j+1])); | |
166 | } | |
167 | ||
168 | Rndm(random,2); | |
169 | time = fTimeOrigin + fOsigma[2]/TMath::Ccgs()* | |
170 | TMath::Cos(2*random[0]*TMath::Pi())* | |
171 | TMath::Sqrt(-2*TMath::Log(random[1])); | |
172 | } | |
173 | ||
174 | PushTrack(fTrackIt,-1,fLegPdg2,p,origin,polar,time,kPPrimary,nt,weight); | |
175 | 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); | |
176 | ||
177 | }//particle loop | |
178 | ||
179 | } | |
180 | ||
181 | //_____________________________________________________________________________ | |
182 | ||
183 | void AliGenPairFlat::Init() | |
184 | { | |
185 | // | |
186 | // Init | |
187 | // | |
188 | ||
189 | printf("AliGenPairFlat::Init() not implemented!!!\n"); | |
190 | ||
191 | } | |
192 | ||
193 | ||
194 | //_____________________________________________________________________________ | |
195 | ||
196 | ||
df8949d4 | 197 | Bool_t AliGenPairFlat::Decay(TLorentzVector& mother, TLorentzVector &dau1, TLorentzVector &dau2, TF1* polfactor) |
67bfd92b | 198 | { |
199 | // | |
200 | // decay procedure | |
201 | // | |
202 | ||
203 | Double_t mp, md1, md2, ed1, ed2, pd1, pd2; | |
204 | Double_t costheta, sintheta, phi; | |
205 | ||
206 | mp = mother.M(); | |
207 | md1 = TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg1))->Mass(); | |
208 | md2 = TDatabasePDG::Instance()->GetParticle(TMath::Abs(fLegPdg2))->Mass(); | |
209 | ||
210 | if(mp < md1+md2){ | |
211 | printf("Daughters are heavier than Mother!! Check Kinematics!! \n"); | |
212 | return 0; | |
213 | } | |
214 | ||
215 | ed1 = (mp*mp + md1*md1 - md2*md2)/(2.*mp); | |
216 | ed2 = mp-ed1; | |
217 | pd1 = TMath::Sqrt((ed1+md1)*(ed1-md1)); | |
218 | pd2 = TMath::Sqrt((mp*mp-(md1+md2)*(md1+md2))*(mp*mp -(md1-md2)*(md1-md2)))/(2.*mp); | |
219 | ||
725547b6 | 220 | costheta = polfactor->GetRandom(); //polarization |
67bfd92b | 221 | sintheta = TMath::Sqrt((1.+costheta)*(1.-costheta)); |
222 | phi = 2.0*TMath::Pi()*gRandom->Rndm(); | |
223 | ||
224 | dau1.SetPxPyPzE(pd1*sintheta*TMath::Cos(phi), | |
225 | pd1*sintheta*TMath::Sin(phi), | |
226 | pd1*costheta, | |
227 | ed1); | |
228 | ||
229 | dau2.SetPxPyPzE((-1.)*pd1*sintheta*TMath::Cos(phi), | |
230 | (-1.)*pd1*sintheta*TMath::Sin(phi), | |
231 | (-1.)*pd1*costheta, | |
232 | ed2); | |
233 | ||
234 | TVector3 boost = mother.BoostVector(); | |
235 | ||
236 | dau1.Boost(boost); | |
237 | dau2.Boost(boost); | |
238 | ||
239 | return 1; | |
240 | ||
241 | } | |
242 | ||
243 | ||
244 |