]>
Commit | Line | Data |
---|---|---|
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 | /* $Id$ */ | |
17 | ||
18 | // | |
19 | // Classe to create the coktail for physics with muons for pp at 14 TeV | |
20 | // The followoing sources: | |
21 | // jpsi,psiP, upsilon, upsilonP, upsilonPP are added to Pythia | |
22 | // The free parameeters are : | |
23 | // pp reaction cross-section | |
24 | // production cross-sections in pp collisions | |
25 | // July 07:added heavy quark production from AliGenCorrHF and heavy quark | |
26 | // production switched off in Pythia | |
27 | // Aug. 07: added trigger cut on total momentum | |
28 | // 2009: added possibility to hide x-sections (B. Vulpescu) | |
29 | // 2009: added possibility to have the cocktail (fast generator and param.) | |
30 | // for pp @ 10 TeV or pp @ 14 TeV (N. Bastid) | |
31 | //----------------- | |
32 | // 2009: added polarization (L. Bianchi) | |
33 | //------------------ | |
34 | // 11/2009: added chi_c1 & chi_c2 (P.Crochet & N.Bastid). | |
35 | // Cross-sections for charmonia are now directly taken from the Yellow Report | |
36 | // (hep-ph/0311048) Tab.9, page 19. See below for details w.r.t. beam energy. | |
37 | // usage: see example of Config in $ALICE_ROOT/prod/LHC09a10/Config.C | |
38 | ||
39 | #include <TObjArray.h> | |
40 | #include <TParticle.h> | |
41 | #include <TF1.h> | |
42 | #include <TVirtualMC.h> | |
43 | ||
44 | #include "AliGenCocktailEventHeader.h" | |
45 | ||
46 | #include "AliGenCocktailEntry.h" | |
47 | #include "AliGenMUONCocktailpp.h" | |
48 | #include "AliGenMUONlib.h" | |
49 | #include "AliGenParam.h" | |
50 | #include "AliMC.h" | |
51 | #include "AliRun.h" | |
52 | #include "AliStack.h" | |
53 | #include "AliDecayer.h" | |
54 | #include "AliLog.h" | |
55 | #include "AliGenCorrHF.h" | |
56 | #include "AliDecayerPolarized.h" | |
57 | ||
58 | ClassImp(AliGenMUONCocktailpp) | |
59 | ||
60 | //________________________________________________________________________ | |
61 | AliGenMUONCocktailpp::AliGenMUONCocktailpp() | |
62 | :AliGenCocktail(), | |
63 | fDecayer(0), | |
64 | fDecayModeResonance(kAll), | |
65 | fDecayModePythia(kAll), | |
66 | fMuonMultiplicity(0), | |
67 | fMuonPtCut(0.), | |
68 | fMuonPCut(0.), | |
69 | fMuonThetaMinCut(0.), | |
70 | fMuonThetaMaxCut(180.), | |
71 | fMuonOriginCut(-999.), | |
72 | fNSucceded(0), | |
73 | fNGenerated(0), | |
74 | ||
75 | fJpsiPol(0), | |
76 | fChic1Pol(0), | |
77 | fChic2Pol(0), | |
78 | fPsiPPol(0), | |
79 | fUpsPol(0), | |
80 | fUpsPPol(0), | |
81 | fUpsPPPol(0), | |
82 | fPolFrame(0), | |
83 | ||
84 | // x-sections for pp @ 7 TeV: charmonia from hep-ph/0311048 Tab.9, page 19, | |
85 | // bottomnium as for 10 TeV | |
86 | /* fCMSEnergy(7), | |
87 | fSigmaReaction(0.0695), | |
88 | fSigmaJPsi(21.8e-6), | |
89 | fSigmaChic1(21.1e-6), | |
90 | fSigmaChic2(34.9e-6), | |
91 | fSigmaPsiP(4.93e-6), | |
92 | fSigmaUpsilon(0.463e-6), | |
93 | fSigmaUpsilonP(0.154e-6), | |
94 | fSigmaUpsilonPP(0.0886e-6), | |
95 | fSigmaCCbar(6.91e-3), | |
96 | fSigmaBBbar(0.232e-3), | |
97 | */ | |
98 | ||
99 | //x-sections for pp @ 10 TeV: charmonia and bottomonia from 14 TeV numbers | |
100 | // scaled down according to ccbar and bbar cross-sections | |
101 | fCMSEnergy(10), | |
102 | fSigmaReaction(0.0695), | |
103 | fSigmaJPsi(26.06e-6), | |
104 | fSigmaChic1(25.18e-6), | |
105 | fSigmaChic2(41.58e-6), | |
106 | fSigmaPsiP(5.88e-6), | |
107 | fSigmaUpsilon(0.658e-6), | |
108 | fSigmaUpsilonP(0.218e-6), | |
109 | fSigmaUpsilonPP(0.122e-6), | |
110 | fSigmaCCbar(8.9e-3), | |
111 | fSigmaBBbar(0.33e-3), | |
112 | ||
113 | ||
114 | //x-sections for pp @ 14 TeV: charmonia from hep-ph/0311048 Tab.9, page 19, | |
115 | // bottomonium from hep-ph/0311048 Tab.9, page 19 taken inton account that | |
116 | // feed-down from chib is included | |
117 | /* fCMSEnergy(14), | |
118 | fSigmaReaction(0.070), | |
119 | fSigmaJPsi(32.9e-6), | |
120 | fSigmaChic1(31.8e-6), | |
121 | fSigmaChic2(52.5e-6), | |
122 | fSigmaPsiP(7.43e-6), | |
123 | fSigmaUpsilon(0.989e-6), | |
124 | fSigmaUpsilonP(0.502e-6), | |
125 | fSigmaUpsilonPP(0.228e-6), | |
126 | fSigmaCCbar(11.2e-3), | |
127 | fSigmaBBbar(0.51e-3), | |
128 | */ | |
129 | fSigmaSilent(kFALSE) | |
130 | { | |
131 | // Constructor | |
132 | } | |
133 | ||
134 | //_________________________________________________________________________ | |
135 | AliGenMUONCocktailpp::~AliGenMUONCocktailpp() | |
136 | // Destructor | |
137 | { | |
138 | ||
139 | } | |
140 | ||
141 | //_________________________________________________________________________ | |
142 | void AliGenMUONCocktailpp::SetResPolarization(Double_t JpsiPol, Double_t PsiPPol, Double_t UpsPol, | |
143 | Double_t UpsPPol, Double_t UpsPPPol, char *PolFrame){ | |
144 | ||
145 | if (strcmp(PolFrame,"kColSop")==0){ | |
146 | fJpsiPol = (JpsiPol>=-1 && JpsiPol<=1) ? JpsiPol : 0; | |
147 | fPsiPPol = (PsiPPol>=-1 && PsiPPol<=1) ? PsiPPol : 0; | |
148 | fUpsPol = (UpsPol>=-1 && UpsPol<=1) ? UpsPol : 0; | |
149 | fUpsPPol = (UpsPPol>=-1 && UpsPPol<=1) ? UpsPPol : 0; | |
150 | fUpsPPPol = (UpsPPPol>=-1 && UpsPPPol<=1) ? UpsPPPol : 0; | |
151 | fPolFrame = 0; | |
152 | } else if (strcmp(PolFrame,"kHelicity")==0){ | |
153 | fJpsiPol = (JpsiPol>=-1 && JpsiPol<=1) ? JpsiPol : 0; | |
154 | fPsiPPol = (PsiPPol>=-1 && PsiPPol<=1) ? PsiPPol : 0; | |
155 | fUpsPol = (UpsPol>=-1 && UpsPol<=1) ? UpsPol : 0; | |
156 | fUpsPPol = (UpsPPol>=-1 && UpsPPol<=1) ? UpsPPol : 0; | |
157 | fUpsPPPol = (UpsPPPol>=-1 && UpsPPPol<=1) ? UpsPPPol : 0; | |
158 | fPolFrame = 1; | |
159 | ||
160 | } else { | |
161 | AliInfo(Form("The polarization frame is not valid")); | |
162 | AliInfo(Form("No polarization will be set")); | |
163 | fJpsiPol=0.; | |
164 | fPsiPPol=0.; | |
165 | fUpsPol=0.; | |
166 | fUpsPPol=0.; | |
167 | fUpsPPPol=0.; | |
168 | } | |
169 | } | |
170 | ||
171 | //_________________________________________________________________________ | |
172 | void AliGenMUONCocktailpp::CreateCocktail() | |
173 | { | |
174 | ||
175 | // These limits are only used for renormalization of quarkonia cross section | |
176 | // Pythia events are generated in 4pi | |
177 | Double_t ptMin = fPtMin; | |
178 | Double_t ptMax = fPtMax; | |
179 | Double_t yMin = fYMin;; | |
180 | Double_t yMax = fYMax;; | |
181 | Double_t phiMin = fPhiMin*180./TMath::Pi(); | |
182 | Double_t phiMax = fPhiMax*180./TMath::Pi(); | |
183 | AliInfo(Form("Ranges pT:%4.1f : %4.1f GeV/c, y:%4.2f : %4.2f, Phi:%5.1f : %5.1f degres",ptMin,ptMax,yMin,yMax,phiMin,phiMax)); | |
184 | ||
185 | // Ratios with respect to the reaction cross-section in the | |
186 | // kinematics limit of the MUONCocktail | |
187 | Double_t ratiojpsi; | |
188 | Double_t ratiochic1; | |
189 | Double_t ratiochic2; | |
190 | Double_t ratiopsiP; | |
191 | Double_t ratioupsilon; | |
192 | Double_t ratioupsilonP; | |
193 | Double_t ratioupsilonPP; | |
194 | Double_t ratioccbar; | |
195 | Double_t ratiobbbar; | |
196 | ||
197 | // Beam energy | |
198 | Int_t cmsEnergy = fCMSEnergy; | |
199 | ||
200 | // Cross sections in barns (from PPR Vol. II p: 552) pp - 14 TeV and | |
201 | // corrected from feed down of higher resonances | |
202 | ||
203 | Double_t sigmajpsi = fSigmaJPsi; | |
204 | Double_t sigmachic1 = fSigmaChic1; | |
205 | Double_t sigmachic2 = fSigmaChic2; | |
206 | Double_t sigmapsiP = fSigmaPsiP; | |
207 | Double_t sigmaupsilon = fSigmaUpsilon; | |
208 | Double_t sigmaupsilonP = fSigmaUpsilonP; | |
209 | Double_t sigmaupsilonPP = fSigmaUpsilonPP; | |
210 | Double_t sigmaccbar = fSigmaCCbar; | |
211 | Double_t sigmabbbar = fSigmaBBbar; | |
212 | ||
213 | // Cross sections corrected with the BR in mu+mu- | |
214 | // (only in case of use of AliDecayerPolarized) | |
215 | ||
216 | if(TMath::Abs(fJpsiPol) > 1.e-30) {sigmajpsi = fSigmaJPsi*0.0593;} | |
217 | if(TMath::Abs(fChic1Pol) > 1.e-30) {sigmachic1 = fSigmaChic1*0.;} // tb consistent | |
218 | if(TMath::Abs(fChic2Pol) > 1.e-30) {sigmachic2 = fSigmaChic2*0.;} // tb consistent | |
219 | if(TMath::Abs(fPsiPPol) > 1.e-30) {sigmapsiP = fSigmaPsiP*0.0075;} | |
220 | if(TMath::Abs(fUpsPol) > 1.e-30) {sigmaupsilon = fSigmaUpsilon*0.0248;} | |
221 | if(TMath::Abs(fUpsPPol) > 1.e-30) {sigmaupsilonP = fSigmaUpsilonP*0.0193;} | |
222 | if(TMath::Abs(fUpsPPPol) > 1.e-30) {sigmaupsilonPP = fSigmaUpsilonPP*0.0218;} | |
223 | ||
224 | // MinBias NN cross section @ pp 14 TeV -PR Vol II p:473 | |
225 | Double_t sigmaReaction = fSigmaReaction; | |
226 | ||
227 | Int_t eincStart = 10; | |
228 | ||
229 | AliInfo(Form("the parametrised resonances uses the decay mode %d",fDecayModeResonance)); | |
230 | ||
231 | // Generation using CDF scaled distribution for pp@14TeV (from D. Stocco) | |
232 | //---------------------------------------------------------------------- | |
233 | // Jpsi generator | |
234 | AliGenParam * genjpsi; | |
235 | if(cmsEnergy == eincStart){ | |
236 | genjpsi = new AliGenParam(1, AliGenMUONlib::kJpsi, "CDF pp 10", "Jpsi"); | |
237 | } else { | |
238 | genjpsi = new AliGenParam(1, AliGenMUONlib::kJpsi, "CDF pp ", "Jpsi"); | |
239 | } | |
240 | // first step: generation in 4pi | |
241 | genjpsi->SetPtRange(0.,100.); | |
242 | genjpsi->SetYRange(-8.,8.); | |
243 | genjpsi->SetPhiRange(0.,360.); | |
244 | genjpsi->SetForceDecay(fDecayModeResonance); | |
245 | //if (!gMC) genjpsi->SetDecayer(fDecayer); | |
246 | ||
247 | genjpsi->Init(); // generation in 4pi | |
248 | ratiojpsi = sigmajpsi / sigmaReaction * genjpsi->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); // get weight | |
249 | if (!fSigmaSilent) { | |
250 | AliInfo(Form("jpsi prod. cross-section in pp %5.3g b",sigmajpsi)); | |
251 | AliInfo(Form("jpsi prod. probability per collision in acceptance %5.3g",ratiojpsi)); | |
252 | } | |
253 | // second step: generation in selected kinematical range | |
254 | genjpsi->SetPtRange(ptMin, ptMax); | |
255 | genjpsi->SetYRange(yMin, yMax); | |
256 | genjpsi->SetPhiRange(phiMin, phiMax); | |
257 | genjpsi->Init(); // generation in selected kinematic range | |
258 | ||
259 | TVirtualMCDecayer *jpsiDec = 0; | |
260 | if(TMath::Abs(fJpsiPol) > 1.e-30){ | |
261 | AliInfo(Form("******Setting polarized decayer for J/psi")); | |
262 | if(fPolFrame==0) { | |
263 | jpsiDec = new AliDecayerPolarized(fJpsiPol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon); | |
264 | AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fJpsiPol)); | |
265 | } | |
266 | if(fPolFrame==1) { | |
267 | jpsiDec = new AliDecayerPolarized(fJpsiPol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon); | |
268 | AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fJpsiPol)); | |
269 | } | |
270 | jpsiDec->SetForceDecay(kAll); | |
271 | jpsiDec->Init(); | |
272 | genjpsi->SetDecayer(jpsiDec); | |
273 | } | |
274 | ||
275 | AddGenerator(genjpsi, "Jpsi", ratiojpsi); // Adding Generator | |
276 | ||
277 | ||
278 | //---------------------------------------------------------------------- | |
279 | // Chic1 generator | |
280 | AliGenParam * genchic1; | |
281 | if(cmsEnergy == eincStart){ | |
282 | genchic1 = new AliGenParam(1, AliGenMUONlib::kChic1, "CDF pp 10", "Chic1"); | |
283 | } else { | |
284 | genchic1 = new AliGenParam(1, AliGenMUONlib::kChic1, "CDF pp ", "Chic1"); | |
285 | } | |
286 | // first step: generation in 4pi | |
287 | genchic1->SetPtRange(0.,100.); | |
288 | genchic1->SetYRange(-8.,8.); | |
289 | genchic1->SetPhiRange(0.,360.); | |
290 | genchic1->SetForceDecay(fDecayModeResonance); | |
291 | //if (!gMC) genchic1->SetDecayer(fDecayer); | |
292 | ||
293 | genchic1->Init(); // generation in 4pi | |
294 | ratiochic1 = sigmachic1 / sigmaReaction * genchic1->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); // get weight | |
295 | if (!fSigmaSilent) { | |
296 | AliInfo(Form("chic1 prod. cross-section in pp %5.3g b",sigmachic1)); | |
297 | AliInfo(Form("chic1 prod. probability per collision in acceptance %5.3g",ratiochic1)); | |
298 | } | |
299 | // second step: generation in selected kinematical range | |
300 | genchic1->SetPtRange(ptMin, ptMax); | |
301 | genchic1->SetYRange(yMin, yMax); | |
302 | genchic1->SetPhiRange(phiMin, phiMax); | |
303 | genchic1->Init(); // generation in selected kinematic range | |
304 | ||
305 | TVirtualMCDecayer *chic1Dec = 0; | |
306 | if(fChic1Pol != 0){ | |
307 | AliInfo(Form("******Setting polarized decayer for Chic1")); | |
308 | if(fPolFrame==0) { | |
309 | chic1Dec = new AliDecayerPolarized(fChic1Pol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon); | |
310 | AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fChic1Pol)); | |
311 | } | |
312 | if(fPolFrame==1) { | |
313 | chic1Dec = new AliDecayerPolarized(fChic1Pol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon); | |
314 | AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fChic1Pol)); | |
315 | } | |
316 | chic1Dec->SetForceDecay(kAll); | |
317 | chic1Dec->Init(); | |
318 | genchic1->SetDecayer(chic1Dec); | |
319 | } | |
320 | ||
321 | AddGenerator(genchic1, "Chic1", ratiochic1); // Adding Generator | |
322 | ||
323 | //---------------------------------------------------------------------- | |
324 | // Chic2 generator | |
325 | AliGenParam * genchic2; | |
326 | if(cmsEnergy == eincStart){ | |
327 | genchic2 = new AliGenParam(1, AliGenMUONlib::kChic2, "CDF pp 10", "Chic2"); | |
328 | } else { | |
329 | genchic2 = new AliGenParam(1, AliGenMUONlib::kChic2, "CDF pp ", "Chic2"); | |
330 | } | |
331 | // first step: generation in 4pi | |
332 | genchic2->SetPtRange(0.,100.); | |
333 | genchic2->SetYRange(-8.,8.); | |
334 | genchic2->SetPhiRange(0.,360.); | |
335 | genchic2->SetForceDecay(fDecayModeResonance); | |
336 | //if (!gMC) genchic1->SetDecayer(fDecayer); | |
337 | ||
338 | genchic2->Init(); // generation in 4pi | |
339 | ratiochic2 = sigmachic2 / sigmaReaction * genchic2->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); // get weight | |
340 | if (!fSigmaSilent) { | |
341 | AliInfo(Form("chic2 prod. cross-section in pp %5.3g b",sigmachic2)); | |
342 | AliInfo(Form("chic2 prod. probability per collision in acceptance %5.3g",ratiochic2)); | |
343 | } | |
344 | // second step: generation in selected kinematical range | |
345 | genchic2->SetPtRange(ptMin, ptMax); | |
346 | genchic2->SetYRange(yMin, yMax); | |
347 | genchic2->SetPhiRange(phiMin, phiMax); | |
348 | genchic2->Init(); // generation in selected kinematic range | |
349 | ||
350 | TVirtualMCDecayer *chic2Dec = 0; | |
351 | if(fChic2Pol != 0){ | |
352 | AliInfo(Form("******Setting polarized decayer for Chic2")); | |
353 | if(fPolFrame==0) { | |
354 | chic2Dec = new AliDecayerPolarized(fChic2Pol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon); | |
355 | AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fChic2Pol)); | |
356 | } | |
357 | if(fPolFrame==1) { | |
358 | chic2Dec = new AliDecayerPolarized(fChic2Pol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon); | |
359 | AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fChic2Pol)); | |
360 | } | |
361 | chic2Dec->SetForceDecay(kAll); | |
362 | chic2Dec->Init(); | |
363 | genchic2->SetDecayer(chic2Dec); | |
364 | } | |
365 | ||
366 | AddGenerator(genchic2, "Chic2", ratiochic2); // Adding Generator | |
367 | ||
368 | //------------------------------------------------------------------ | |
369 | // Psi prime generator | |
370 | AliGenParam * genpsiP; | |
371 | if(cmsEnergy == eincStart){ | |
372 | genpsiP = new AliGenParam(1, AliGenMUONlib::kPsiP, "CDF pp 10", "PsiP"); | |
373 | } else { | |
374 | genpsiP = new AliGenParam(1, AliGenMUONlib::kPsiP, "CDF pp", "PsiP"); | |
375 | } | |
376 | // first step: generation in 4pi | |
377 | genpsiP->SetPtRange(0.,100.); | |
378 | genpsiP->SetYRange(-8.,8.); | |
379 | genpsiP->SetPhiRange(0.,360.); | |
380 | genpsiP->SetForceDecay(fDecayModeResonance); | |
381 | //if (!gMC) genpsiP->SetDecayer(fDecayer); | |
382 | ||
383 | genpsiP->Init(); // generation in 4pi | |
384 | ratiopsiP = sigmapsiP / sigmaReaction * genpsiP->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); | |
385 | if (!fSigmaSilent) { | |
386 | AliInfo(Form("psiP prod. cross-section in pp %5.3g b",sigmapsiP)); | |
387 | AliInfo(Form("psiP prod. probability per collision in acceptance %5.3g",ratiopsiP)); | |
388 | } | |
389 | // second step: generation in selected kinematical range | |
390 | genpsiP->SetPtRange(ptMin, ptMax); | |
391 | genpsiP->SetYRange(yMin, yMax); | |
392 | genpsiP->SetPhiRange(phiMin, phiMax); | |
393 | genpsiP->Init(); // generation in selected kinematic range | |
394 | ||
395 | TVirtualMCDecayer *psipDec = 0; | |
396 | if(TMath::Abs(fPsiPPol) > 1.e-30){ | |
397 | AliInfo(Form("******Setting polarized decayer for psi'")); | |
398 | if(fPolFrame==0) { | |
399 | psipDec = new AliDecayerPolarized(fPsiPPol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon); | |
400 | AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fPsiPPol)); | |
401 | } | |
402 | if(fPolFrame==1) { | |
403 | psipDec = new AliDecayerPolarized(fPsiPPol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon); | |
404 | AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fPsiPPol)); | |
405 | } | |
406 | psipDec->SetForceDecay(kAll); | |
407 | psipDec->Init(); | |
408 | genpsiP->SetDecayer(psipDec); | |
409 | } | |
410 | ||
411 | AddGenerator(genpsiP, "PsiP", ratiopsiP); // Adding Generator | |
412 | ||
413 | //------------------------------------------------------------------ | |
414 | // Upsilon 1S generator | |
415 | AliGenParam * genupsilon; | |
416 | if(cmsEnergy == eincStart) { | |
417 | genupsilon = new AliGenParam(1, AliGenMUONlib::kUpsilon, "CDF pp 10", "Upsilon"); | |
418 | } else { | |
419 | genupsilon = new AliGenParam(1, AliGenMUONlib::kUpsilon, "CDF pp", "Upsilon"); | |
420 | } | |
421 | // first step: generation in 4pi | |
422 | genupsilon->SetPtRange(0.,100.); | |
423 | genupsilon->SetYRange(-8.,8.); | |
424 | genupsilon->SetPhiRange(0.,360.); | |
425 | genupsilon->SetForceDecay(fDecayModeResonance); | |
426 | //if (!gMC) genupsilon->SetDecayer(fDecayer); | |
427 | genupsilon->Init(); // generation in 4pi | |
428 | ratioupsilon = sigmaupsilon / sigmaReaction * genupsilon->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); | |
429 | if (!fSigmaSilent) { | |
430 | AliInfo(Form("Upsilon 1S prod. cross-section in pp %5.3g b",sigmaupsilon)); | |
431 | AliInfo(Form("Upsilon 1S prod. probability per collision in acceptance %5.3g",ratioupsilon)); | |
432 | } | |
433 | // second step: generation in selected kinematical range | |
434 | genupsilon->SetPtRange(ptMin, ptMax); | |
435 | genupsilon->SetYRange(yMin, yMax); | |
436 | genupsilon->SetPhiRange(phiMin, phiMax); | |
437 | genupsilon->Init(); // generation in selected kinematical range | |
438 | ||
439 | TVirtualMCDecayer *upsDec = 0; | |
440 | if(TMath::Abs(fUpsPol) > 1.e-30){ | |
441 | AliInfo(Form("******Setting polarized decayer for Upsilon")); | |
442 | if(fPolFrame==0) { | |
443 | upsDec = new AliDecayerPolarized(fUpsPol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon); | |
444 | AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fUpsPol)); | |
445 | } | |
446 | if(fPolFrame==1) { | |
447 | upsDec = new AliDecayerPolarized(fUpsPol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon); | |
448 | AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fUpsPol)); | |
449 | } | |
450 | upsDec->SetForceDecay(kAll); | |
451 | upsDec->Init(); | |
452 | genupsilon->SetDecayer(upsDec); | |
453 | } | |
454 | ||
455 | AddGenerator(genupsilon,"Upsilon", ratioupsilon); // Adding Generator | |
456 | //------------------------------------------------------------------ | |
457 | // Upsilon 2S generator | |
458 | AliGenParam * genupsilonP; | |
459 | if(cmsEnergy == eincStart){ | |
460 | genupsilonP = new AliGenParam(1, AliGenMUONlib::kUpsilonP, "CDF pp 10", "UpsilonP"); | |
461 | } else { | |
462 | genupsilonP = new AliGenParam(1, AliGenMUONlib::kUpsilonP, "CDF pp", "UpsilonP"); | |
463 | } | |
464 | ||
465 | // first step: generation in 4pi | |
466 | genupsilonP->SetPtRange(0.,100.); | |
467 | genupsilonP->SetYRange(-8.,8.); | |
468 | genupsilonP->SetPhiRange(0.,360.); | |
469 | genupsilonP->SetForceDecay(fDecayModeResonance); | |
470 | //if (!gMC) genupsilonP->SetDecayer(fDecayer); | |
471 | genupsilonP->Init(); // generation in 4pi | |
472 | ratioupsilonP = sigmaupsilonP / sigmaReaction * genupsilonP->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); | |
473 | if (!fSigmaSilent) { | |
474 | AliInfo(Form("Upsilon 2S prod. cross-section in pp %5.3g b",sigmaupsilonP)); | |
475 | AliInfo(Form("Upsilon 2S prod. probability per collision in acceptance %5.3g",ratioupsilonP)); | |
476 | } | |
477 | // second step: generation in the kinematical range | |
478 | genupsilonP->SetPtRange(ptMin, ptMax); | |
479 | genupsilonP->SetYRange(yMin, yMax); | |
480 | genupsilonP->SetPhiRange(phiMin, phiMax); | |
481 | genupsilonP->Init(); // generation in selected kinematical range | |
482 | ||
483 | TVirtualMCDecayer *upspDec = 0; | |
484 | if(TMath::Abs(fUpsPPol) > 1.e-30){ | |
485 | AliInfo(Form("******Setting polarized decayer for Upsilon'")); | |
486 | if(fPolFrame==0) { | |
487 | upspDec = new AliDecayerPolarized(fUpsPPol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon); | |
488 | AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fUpsPPol)); | |
489 | } | |
490 | if(fPolFrame==1) { | |
491 | upspDec = new AliDecayerPolarized(fUpsPPol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon); | |
492 | AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fUpsPPol)); | |
493 | } | |
494 | upspDec->SetForceDecay(kAll); | |
495 | upspDec->Init(); | |
496 | genupsilonP->SetDecayer(upspDec); | |
497 | } | |
498 | ||
499 | AddGenerator(genupsilonP,"UpsilonP", ratioupsilonP); // Adding Generator | |
500 | ||
501 | //------------------------------------------------------------------ | |
502 | // Upsilon 3S generator | |
503 | AliGenParam * genupsilonPP; | |
504 | if(cmsEnergy == eincStart){ | |
505 | genupsilonPP = new AliGenParam(1, AliGenMUONlib::kUpsilonPP, "CDF pp 10", "UpsilonPP"); | |
506 | } | |
507 | else { | |
508 | genupsilonPP = new AliGenParam(1, AliGenMUONlib::kUpsilonPP, "CDF pp", "UpsilonPP"); | |
509 | } | |
510 | ||
511 | // first step: generation in 4pi | |
512 | genupsilonPP->SetPtRange(0.,100.); | |
513 | genupsilonPP->SetYRange(-8.,8.); | |
514 | genupsilonPP->SetPhiRange(0.,360.); | |
515 | genupsilonPP->SetForceDecay(fDecayModeResonance); | |
516 | //if (!gMC) genupsilonPP->SetDecayer(fDecayer); | |
517 | genupsilonPP->Init(); // generation in 4pi | |
518 | ratioupsilonPP = sigmaupsilonPP / sigmaReaction * genupsilonPP->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); | |
519 | if (!fSigmaSilent) { | |
520 | AliInfo(Form("Upsilon 3S prod. cross-section in pp %5.3g b",sigmaupsilonPP)); | |
521 | AliInfo(Form("Upsilon 3S prod. probability per collision in acceptance %5.3g",ratioupsilonPP)); | |
522 | } | |
523 | // second step: generation in selected kinematical range | |
524 | genupsilonPP->SetPtRange(ptMin, ptMax); | |
525 | genupsilonPP->SetYRange(yMin, yMax); | |
526 | genupsilonPP->SetPhiRange(phiMin, phiMax); | |
527 | genupsilonPP->Init(); // generation in selected kinematical range | |
528 | ||
529 | TVirtualMCDecayer *upsppDec = 0; | |
530 | if(fUpsPPPol != 0){ | |
531 | AliInfo(Form("******Setting polarized decayer for Upsilon''")); | |
532 | if(fPolFrame==0) { | |
533 | upsppDec = new AliDecayerPolarized(fUpsPPPol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon); | |
534 | AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fUpsPPPol)); | |
535 | } | |
536 | if(fPolFrame==1) { | |
537 | upsppDec = new AliDecayerPolarized(fUpsPPPol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon); | |
538 | AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fUpsPPPol)); | |
539 | } | |
540 | upsppDec->SetForceDecay(kAll); | |
541 | upsppDec->Init(); | |
542 | genupsilonPP->SetDecayer(upsppDec); | |
543 | } | |
544 | ||
545 | AddGenerator(genupsilonPP,"UpsilonPP", ratioupsilonPP); // Adding Generator | |
546 | ||
547 | //------------------------------------------------------------------ | |
548 | // Generator of charm | |
549 | ||
550 | AliGenCorrHF *gencharm = new AliGenCorrHF(1, 4, cmsEnergy); | |
551 | gencharm->SetMomentumRange(0,9999); | |
552 | gencharm->SetForceDecay(kAll); | |
553 | ratioccbar = sigmaccbar/sigmaReaction; | |
554 | //if (!gMC) gencharm->SetDecayer(fDecayer); | |
555 | gencharm->Init(); | |
556 | if (!fSigmaSilent) { | |
557 | AliInfo(Form("c-cbar prod. cross-section in pp %5.3g b",sigmaccbar)); | |
558 | AliInfo(Form("c-cbar prod. probability per collision in acceptance %5.3g",ratioccbar)); | |
559 | } | |
560 | AddGenerator(gencharm,"CorrHFCharm",ratioccbar); | |
561 | ||
562 | //------------------------------------------------------------------ | |
563 | // Generator of beauty | |
564 | ||
565 | AliGenCorrHF *genbeauty = new AliGenCorrHF(1, 5, cmsEnergy); | |
566 | genbeauty->SetMomentumRange(0,9999); | |
567 | genbeauty->SetForceDecay(kAll); | |
568 | ratiobbbar = sigmabbbar/sigmaReaction; | |
569 | //if (!gMC) genbeauty->SetDecayer(fDecayer); | |
570 | genbeauty->Init(); | |
571 | if (!fSigmaSilent) { | |
572 | AliInfo(Form("b-bbar prod. cross-section in pp %5.3g b",sigmabbbar)); | |
573 | AliInfo(Form("b-bbar prod. probability per collision in acceptance %5.3g",ratiobbbar)); | |
574 | } | |
575 | AddGenerator(genbeauty,"CorrHFBeauty",ratiobbbar); | |
576 | ||
577 | //------------------------------------------------------------------- | |
578 | // Pythia generator | |
579 | // | |
580 | // This has to go into the Config.C | |
581 | // | |
582 | // AliGenPythia *pythia = new AliGenPythia(1); | |
583 | // pythia->SetProcess(kPyMbMSEL1); | |
584 | // pythia->SetStrucFunc(kCTEQ5L); | |
585 | // pythia->SetEnergyCMS(14000.); | |
586 | // AliInfo(Form("\n\npythia uses the decay mode %d", GetDecayModePythia())); | |
587 | // Decay_t dt = gener->GetDecayModePythia(); | |
588 | // pythia->SetForceDecay(dt); | |
589 | // pythia->SetPtRange(0.,100.); | |
590 | // pythia->SetYRange(-8.,8.); | |
591 | // pythia->SetPhiRange(0.,360.); | |
592 | // pythia->SetPtHard(2.76,-1.0); | |
593 | // pythia->SwitchHFOff(); | |
594 | // pythia->Init(); | |
595 | // AddGenerator(pythia,"Pythia",1); | |
596 | ||
597 | } | |
598 | ||
599 | void AliGenMUONCocktailpp::Init() | |
600 | { | |
601 | // | |
602 | // Initialisation | |
603 | TIter next(fEntries); | |
604 | AliGenCocktailEntry *entry; | |
605 | if (fStack) { | |
606 | while((entry = (AliGenCocktailEntry*)next())) { | |
607 | entry->Generator()->SetStack(fStack); | |
608 | } | |
609 | } | |
610 | } | |
611 | ||
612 | //_________________________________________________________________________ | |
613 | void AliGenMUONCocktailpp::Generate() | |
614 | { | |
615 | // Generate event | |
616 | TIter next(fEntries); | |
617 | AliGenCocktailEntry *entry = 0; | |
618 | AliGenCocktailEntry *preventry = 0; | |
619 | AliGenerator* gen = 0; | |
620 | ||
621 | if (fHeader) delete fHeader; | |
622 | fHeader = new AliGenCocktailEventHeader("MUON Cocktail Header"); | |
623 | ||
624 | const TObjArray *partArray = gAlice->GetMCApp()->Particles(); | |
625 | ||
626 | // Generate the vertex position used by all generators | |
627 | if(fVertexSmear == kPerEvent) Vertex(); | |
628 | ||
629 | // Loop on primordialTrigger: | |
630 | // minimum muon multiplicity above a pt cut in a theta acceptance region | |
631 | ||
632 | Bool_t primordialTrigger = kFALSE; | |
633 | while(!primordialTrigger) { | |
634 | //Reseting stack | |
635 | AliRunLoader * runloader = AliRunLoader::Instance(); | |
636 | if (runloader) | |
637 | if (runloader->Stack()) | |
638 | runloader->Stack()->Clean(); | |
639 | // Loop over generators and generate events | |
640 | Int_t igen = 0; | |
641 | Int_t npart = 0; | |
642 | const char* genName = 0; | |
643 | while((entry = (AliGenCocktailEntry*)next())) { | |
644 | gen = entry->Generator(); | |
645 | genName = entry->GetName(); | |
646 | gen->SetVertex(fVertex.At(0), fVertex.At(1), fVertex.At(2)); | |
647 | ||
648 | npart = (strcmp(genName,"Pythia") == 0) ? 1 : | |
649 | gRandom->Poisson(entry->Rate()); | |
650 | ||
651 | if (npart > 0) { | |
652 | igen++; | |
653 | if (igen == 1) entry->SetFirst(0); | |
654 | else entry->SetFirst((partArray->GetEntriesFast())+1); | |
655 | ||
656 | gen->SetNumberParticles(npart); | |
657 | gen->Generate(); | |
658 | entry->SetLast(partArray->GetEntriesFast()); | |
659 | preventry = entry; | |
660 | } | |
661 | } | |
662 | next.Reset(); | |
663 | ||
664 | // Testing primordial trigger: Single muons or dimuons with Pt above a Pt cut | |
665 | // in the muon spectrometer acceptance | |
666 | Int_t iPart; | |
667 | fNGenerated++; | |
668 | Int_t numberOfMuons=0;Int_t maxPart = partArray->GetEntriesFast(); | |
669 | for(iPart=0; iPart<maxPart; iPart++){ | |
670 | ||
671 | TParticle *part = gAlice->GetMCApp()->Particle(iPart); | |
672 | if ( TMath::Abs(part->GetPdgCode()) == 13 ){ | |
673 | if((part->Vz() > fMuonOriginCut) && //take only the muons that decayed before the abs + 1 int. length in C abs | |
674 | (part->Theta()*180./TMath::Pi()>fMuonThetaMinCut) && | |
675 | (part->Theta()*180./TMath::Pi()<fMuonThetaMaxCut) && | |
676 | (part->Pt()>fMuonPtCut) && | |
677 | (part->P()>fMuonPCut)) { | |
678 | numberOfMuons++; | |
679 | } | |
680 | } | |
681 | } | |
682 | if (numberOfMuons >= fMuonMultiplicity) { | |
683 | primordialTrigger = kTRUE; | |
684 | fHeader->SetNProduced(maxPart); | |
685 | } | |
686 | ||
687 | } | |
688 | fNSucceded++; | |
689 | ||
690 | TArrayF eventVertex; | |
691 | eventVertex.Set(3); | |
692 | for (Int_t j=0; j < 3; j++) eventVertex[j] = fVertex[j]; | |
693 | ||
694 | fHeader->SetPrimaryVertex(eventVertex); | |
695 | ||
696 | gAlice->SetGenEventHeader(fHeader); | |
697 | ||
698 | // AliInfo(Form("Generated Events are %d and Succeeded Events are %d",fNGenerated,fNSucceded)); | |
699 | AliDebug(5,Form("Generated Events are %d and Succeeded Events are %d",fNGenerated,fNSucceded)); | |
700 | } | |
701 | ||
702 |