1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
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)
32 // 2009: added polarization (L. Bianchi)
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
39 #include <TObjArray.h>
40 #include <TParticle.h>
42 #include <TVirtualMC.h>
44 #include "AliGenCocktailEventHeader.h"
46 #include "AliGenCocktailEntry.h"
47 #include "AliGenMUONCocktailpp.h"
48 #include "AliGenMUONlib.h"
49 #include "AliGenParam.h"
53 #include "AliDecayer.h"
55 #include "AliGenCorrHF.h"
56 #include "AliDecayerPolarized.h"
58 ClassImp(AliGenMUONCocktailpp)
60 //________________________________________________________________________
61 AliGenMUONCocktailpp::AliGenMUONCocktailpp()
64 fDecayModeResonance(kAll),
65 fDecayModePythia(kAll),
70 fMuonThetaMaxCut(180.),
71 fMuonOriginCut(-999.),
84 // x-sections for pp @ 7 TeV: charmonia from hep-ph/0311048 Tab.9, page 19,
85 // bottomnium as for 10 TeV
87 fSigmaReaction(0.0695),
92 fSigmaUpsilon(0.463e-6),
93 fSigmaUpsilonP(0.154e-6),
94 fSigmaUpsilonPP(0.0886e-6),
96 fSigmaBBbar(0.232e-3),
99 //x-sections for pp @ 10 TeV: charmonia and bottomonia from 14 TeV numbers
100 // scaled down according to ccbar and bbar cross-sections
102 fSigmaReaction(0.0695),
103 fSigmaJPsi(26.06e-6),
104 fSigmaChic1(25.18e-6),
105 fSigmaChic2(41.58e-6),
107 fSigmaUpsilon(0.658e-6),
108 fSigmaUpsilonP(0.218e-6),
109 fSigmaUpsilonPP(0.122e-6),
111 fSigmaBBbar(0.33e-3),
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
118 fSigmaReaction(0.070),
120 fSigmaChic1(31.8e-6),
121 fSigmaChic2(52.5e-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),
134 //_________________________________________________________________________
135 AliGenMUONCocktailpp::~AliGenMUONCocktailpp()
141 //_________________________________________________________________________
142 void AliGenMUONCocktailpp::SetResPolarization(Double_t JpsiPol, Double_t PsiPPol, Double_t UpsPol,
143 Double_t UpsPPol, Double_t UpsPPPol, char *PolFrame){
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;
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;
161 AliInfo(Form("The polarization frame is not valid"));
162 AliInfo(Form("No polarization will be set"));
171 //_________________________________________________________________________
172 void AliGenMUONCocktailpp::CreateCocktail()
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));
185 // Ratios with respect to the reaction cross-section in the
186 // kinematics limit of the MUONCocktail
191 Double_t ratioupsilon;
192 Double_t ratioupsilonP;
193 Double_t ratioupsilonPP;
198 Int_t cmsEnergy = fCMSEnergy;
200 // Cross sections in barns (from PPR Vol. II p: 552) pp - 14 TeV and
201 // corrected from feed down of higher resonances
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;
213 // Cross sections corrected with the BR in mu+mu-
214 // (only in case of use of AliDecayerPolarized)
216 if(fJpsiPol != 0) {sigmajpsi = fSigmaJPsi*0.0593;}
217 if(fChic1Pol != 0) {sigmachic1 = fSigmaChic1*0.;} // tb consistent
218 if(fChic2Pol != 0) {sigmachic2 = fSigmaChic2*0.;} // tb consistent
219 if(fPsiPPol != 0) {sigmapsiP = fSigmaPsiP*0.0075;}
220 if(fUpsPol != 0) {sigmaupsilon = fSigmaUpsilon*0.0248;}
221 if(fUpsPPol != 0) {sigmaupsilonP = fSigmaUpsilonP*0.0193;}
222 if(fUpsPPPol != 0) {sigmaupsilonPP = fSigmaUpsilonPP*0.0218;}
224 // MinBias NN cross section @ pp 14 TeV -PR Vol II p:473
225 Double_t sigmaReaction = fSigmaReaction;
227 Int_t eincStart = 10;
229 AliInfo(Form("the parametrised resonances uses the decay mode %d",fDecayModeResonance));
231 // Generation using CDF scaled distribution for pp@14TeV (from D. Stocco)
232 //----------------------------------------------------------------------
234 AliGenParam * genjpsi;
235 if(cmsEnergy == eincStart){
236 genjpsi = new AliGenParam(1, AliGenMUONlib::kJpsi, "CDF pp 10", "Jpsi");
238 genjpsi = new AliGenParam(1, AliGenMUONlib::kJpsi, "CDF pp ", "Jpsi");
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);
247 genjpsi->Init(); // generation in 4pi
248 ratiojpsi = sigmajpsi / sigmaReaction * genjpsi->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); // get weight
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));
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
259 TVirtualMCDecayer *Jpsidec = 0;
261 AliInfo(Form("******Setting polarized decayer for J/psi"));
263 Jpsidec = new AliDecayerPolarized(fJpsiPol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon);
264 AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fJpsiPol));
267 Jpsidec = new AliDecayerPolarized(fJpsiPol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon);
268 AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fJpsiPol));
270 Jpsidec->SetForceDecay(kAll);
272 genjpsi->SetDecayer(Jpsidec);
275 AddGenerator(genjpsi, "Jpsi", ratiojpsi); // Adding Generator
278 //----------------------------------------------------------------------
280 AliGenParam * genchic1;
281 if(cmsEnergy == eincStart){
282 genchic1 = new AliGenParam(1, AliGenMUONlib::kChic1, "CDF pp 10", "Chic1");
284 genchic1 = new AliGenParam(1, AliGenMUONlib::kChic1, "CDF pp ", "Chic1");
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);
293 genchic1->Init(); // generation in 4pi
294 ratiochic1 = sigmachic1 / sigmaReaction * genchic1->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); // get weight
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));
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
305 TVirtualMCDecayer *Chic1dec = 0;
307 AliInfo(Form("******Setting polarized decayer for Chic1"));
309 Chic1dec = new AliDecayerPolarized(fChic1Pol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon);
310 AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fChic1Pol));
313 Chic1dec = new AliDecayerPolarized(fChic1Pol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon);
314 AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fChic1Pol));
316 Chic1dec->SetForceDecay(kAll);
318 genchic1->SetDecayer(Chic1dec);
321 AddGenerator(genchic1, "Chic1", ratiochic1); // Adding Generator
323 //----------------------------------------------------------------------
325 AliGenParam * genchic2;
326 if(cmsEnergy == eincStart){
327 genchic2 = new AliGenParam(1, AliGenMUONlib::kChic2, "CDF pp 10", "Chic2");
329 genchic2 = new AliGenParam(1, AliGenMUONlib::kChic2, "CDF pp ", "Chic2");
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);
338 genchic2->Init(); // generation in 4pi
339 ratiochic2 = sigmachic2 / sigmaReaction * genchic2->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); // get weight
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));
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
350 TVirtualMCDecayer *Chic2dec = 0;
352 AliInfo(Form("******Setting polarized decayer for Chic2"));
354 Chic2dec = new AliDecayerPolarized(fChic2Pol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon);
355 AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fChic2Pol));
358 Chic2dec = new AliDecayerPolarized(fChic2Pol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon);
359 AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fChic2Pol));
361 Chic2dec->SetForceDecay(kAll);
363 genchic2->SetDecayer(Chic2dec);
366 AddGenerator(genchic2, "Chic2", ratiochic2); // Adding Generator
368 //------------------------------------------------------------------
369 // Psi prime generator
370 AliGenParam * genpsiP;
371 if(cmsEnergy == eincStart){
372 genpsiP = new AliGenParam(1, AliGenMUONlib::kPsiP, "CDF pp 10", "PsiP");
374 genpsiP = new AliGenParam(1, AliGenMUONlib::kPsiP, "CDF pp", "PsiP");
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);
383 genpsiP->Init(); // generation in 4pi
384 ratiopsiP = sigmapsiP / sigmaReaction * genpsiP->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax);
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));
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
395 TVirtualMCDecayer *PsiPdec = 0;
397 AliInfo(Form("******Setting polarized decayer for psi'"));
399 PsiPdec = new AliDecayerPolarized(fPsiPPol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon);
400 AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fPsiPPol));
403 PsiPdec = new AliDecayerPolarized(fPsiPPol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon);
404 AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fPsiPPol));
406 PsiPdec->SetForceDecay(kAll);
408 genpsiP->SetDecayer(PsiPdec);
411 AddGenerator(genpsiP, "PsiP", ratiopsiP); // Adding Generator
413 //------------------------------------------------------------------
414 // Upsilon 1S generator
415 AliGenParam * genupsilon;
416 if(cmsEnergy == eincStart) {
417 genupsilon = new AliGenParam(1, AliGenMUONlib::kUpsilon, "CDF pp 10", "Upsilon");
419 genupsilon = new AliGenParam(1, AliGenMUONlib::kUpsilon, "CDF pp", "Upsilon");
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);
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));
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
439 TVirtualMCDecayer *Upsdec = 0;
441 AliInfo(Form("******Setting polarized decayer for Upsilon"));
443 Upsdec = new AliDecayerPolarized(fUpsPol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon);
444 AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fUpsPol));
447 Upsdec = new AliDecayerPolarized(fUpsPol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon);
448 AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fUpsPol));
450 Upsdec->SetForceDecay(kAll);
452 genupsilon->SetDecayer(Upsdec);
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");
462 genupsilonP = new AliGenParam(1, AliGenMUONlib::kUpsilonP, "CDF pp", "UpsilonP");
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);
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));
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
483 TVirtualMCDecayer *UpsPdec = 0;
485 AliInfo(Form("******Setting polarized decayer for Upsilon'"));
487 UpsPdec = new AliDecayerPolarized(fUpsPPol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon);
488 AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fUpsPPol));
491 UpsPdec = new AliDecayerPolarized(fUpsPPol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon);
492 AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fUpsPPol));
494 UpsPdec->SetForceDecay(kAll);
496 genupsilonP->SetDecayer(UpsPdec);
499 AddGenerator(genupsilonP,"UpsilonP", ratioupsilonP); // Adding Generator
501 //------------------------------------------------------------------
502 // Upsilon 3S generator
503 AliGenParam * genupsilonPP;
504 if(cmsEnergy == eincStart){
505 genupsilonPP = new AliGenParam(1, AliGenMUONlib::kUpsilonPP, "CDF pp 10", "UpsilonPP");
508 genupsilonPP = new AliGenParam(1, AliGenMUONlib::kUpsilonPP, "CDF pp", "UpsilonPP");
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);
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));
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
529 TVirtualMCDecayer *UpsPPdec = 0;
531 AliInfo(Form("******Setting polarized decayer for Upsilon''"));
533 UpsPPdec = new AliDecayerPolarized(fUpsPPPol,AliDecayerPolarized::kColSop,AliDecayerPolarized::kMuon);
534 AliInfo(Form("******Reference frame: %s, alpha: %f","Collins-Soper",fUpsPPPol));
537 UpsPPdec = new AliDecayerPolarized(fUpsPPPol,AliDecayerPolarized::kHelicity,AliDecayerPolarized::kMuon);
538 AliInfo(Form("******Reference frame: %s, alpha: %f","Helicity",fUpsPPPol));
540 UpsPPdec->SetForceDecay(kAll);
542 genupsilonPP->SetDecayer(UpsPPdec);
545 AddGenerator(genupsilonPP,"UpsilonPP", ratioupsilonPP); // Adding Generator
547 //------------------------------------------------------------------
548 // Generator of charm
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);
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));
560 AddGenerator(gencharm,"CorrHFCharm",ratioccbar);
562 //------------------------------------------------------------------
563 // Generator of beauty
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);
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));
575 AddGenerator(genbeauty,"CorrHFBeauty",ratiobbbar);
577 //-------------------------------------------------------------------
580 // This has to go into the Config.C
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();
595 // AddGenerator(pythia,"Pythia",1);
599 void AliGenMUONCocktailpp::Init()
603 TIter next(fEntries);
604 AliGenCocktailEntry *entry;
606 while((entry = (AliGenCocktailEntry*)next())) {
607 entry->Generator()->SetStack(fStack);
612 //_________________________________________________________________________
613 void AliGenMUONCocktailpp::Generate()
616 TIter next(fEntries);
617 AliGenCocktailEntry *entry = 0;
618 AliGenCocktailEntry *preventry = 0;
619 AliGenerator* gen = 0;
621 if (fHeader) delete fHeader;
622 fHeader = new AliGenCocktailEventHeader("MUON Cocktail Header");
624 const TObjArray *partArray = gAlice->GetMCApp()->Particles();
626 // Generate the vertex position used by all generators
627 if(fVertexSmear == kPerEvent) Vertex();
629 // Loop on primordialTrigger:
630 // minimum muon multiplicity above a pt cut in a theta acceptance region
632 Bool_t primordialTrigger = kFALSE;
633 while(!primordialTrigger) {
635 AliRunLoader * runloader = AliRunLoader::Instance();
637 if (runloader->Stack())
638 runloader->Stack()->Clean();
639 // Loop over generators and generate events
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));
648 npart = (strcmp(genName,"Pythia") == 0) ? 1 :
649 gRandom->Poisson(entry->Rate());
653 if (igen == 1) entry->SetFirst(0);
654 else entry->SetFirst((partArray->GetEntriesFast())+1);
656 gen->SetNumberParticles(npart);
658 entry->SetLast(partArray->GetEntriesFast());
664 // Testing primordial trigger: Single muons or dimuons with Pt above a Pt cut
665 // in the muon spectrometer acceptance
668 Int_t numberOfMuons=0;Int_t maxPart = partArray->GetEntriesFast();
669 for(iPart=0; iPart<maxPart; iPart++){
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)) {
682 if (numberOfMuons >= fMuonMultiplicity) {
683 primordialTrigger = kTRUE;
684 fHeader->SetNProduced(maxPart);
692 for (Int_t j=0; j < 3; j++) eventVertex[j] = fVertex[j];
694 fHeader->SetPrimaryVertex(eventVertex);
696 gAlice->SetGenEventHeader(fHeader);
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));