]>
Commit | Line | Data |
---|---|---|
103ac317 | 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 | |
b44c3901 | 25 | // July 07:added heavy quark production from AliGenCorrHF and heavy quark |
26 | // production switched off in Pythia | |
f81a4aca | 27 | // Aug. 07: added trigger cut on total momentum |
103ac317 | 28 | |
29 | #include <TObjArray.h> | |
30 | #include <TParticle.h> | |
31 | #include <TF1.h> | |
aaa95f22 | 32 | #include <TVirtualMC.h> |
103ac317 | 33 | |
34 | #include "AliGenCocktailEventHeader.h" | |
35 | ||
36 | #include "AliGenCocktailEntry.h" | |
37 | #include "AliGenMUONCocktailpp.h" | |
38 | #include "AliGenMUONlib.h" | |
39 | #include "AliGenParam.h" | |
40 | #include "AliMC.h" | |
41 | #include "AliRun.h" | |
42 | #include "AliStack.h" | |
aaa95f22 | 43 | #include "AliDecayer.h" |
103ac317 | 44 | #include "AliLog.h" |
b44c3901 | 45 | #include "AliGenCorrHF.h" |
103ac317 | 46 | |
47 | ClassImp(AliGenMUONCocktailpp) | |
48 | ||
49 | //________________________________________________________________________ | |
50 | AliGenMUONCocktailpp::AliGenMUONCocktailpp() | |
1c56e311 | 51 | :AliGenCocktail(), |
aaa95f22 | 52 | fDecayer(0), |
53 | fDecayModeResonance(kAll), | |
54 | fDecayModePythia(kAll), | |
1c56e311 | 55 | fMuonMultiplicity(0), |
56 | fMuonPtCut(0.), | |
f81a4aca | 57 | fMuonPCut(0.), |
18bc0c8e | 58 | fMuonThetaMinCut(0.), |
f81a4aca | 59 | fMuonThetaMaxCut(180.), |
aaa95f22 | 60 | fMuonOriginCut(-999.), |
18bc0c8e | 61 | fNSucceded(0), |
9ff768ee | 62 | fNGenerated(0), |
63 | fSigmaJPsi(49.44e-6), | |
64 | fSigmaPsiP(7.67e-6), | |
65 | fSigmaUpsilon(0.989e-6), | |
66 | fSigmaUpsilonP(0.502e-6), | |
67 | fSigmaUpsilonPP(0.228e-6), | |
68 | fSigmaCCbar(11.2e-3), | |
69 | fSigmaBBbar(0.51e-3), | |
70 | fSigmaSilent(kFALSE) | |
103ac317 | 71 | { |
72 | // Constructor | |
103ac317 | 73 | } |
93a2041b | 74 | |
103ac317 | 75 | //_________________________________________________________________________ |
76 | AliGenMUONCocktailpp::~AliGenMUONCocktailpp() | |
77 | // Destructor | |
78 | { | |
79 | ||
80 | } | |
81 | ||
82 | //_________________________________________________________________________ | |
8058ead1 | 83 | void AliGenMUONCocktailpp::CreateCocktail() |
103ac317 | 84 | { |
85 | // MinBias NN cross section @ pp 14 TeV -PR Vol II p:473 | |
86 | Double_t sigmaReaction = 0.070; | |
87 | ||
88 | // These limits are only used for renormalization of quarkonia cross section | |
89 | // Pythia events are generated in 4pi | |
90 | Double_t ptMin = fPtMin; | |
91 | Double_t ptMax = fPtMax; | |
92 | Double_t yMin = fYMin;; | |
93 | Double_t yMax = fYMax;; | |
94 | Double_t phiMin = fPhiMin*180./TMath::Pi(); | |
95 | Double_t phiMax = fPhiMax*180./TMath::Pi(); | |
96 | 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)); | |
97 | ||
98 | // Ratios with respect to the reaction cross-section in the | |
99 | // kinematics limit of the MUONCocktail | |
100 | Double_t ratiojpsi; | |
101 | Double_t ratiopsiP; | |
102 | Double_t ratioupsilon; | |
103 | Double_t ratioupsilonP; | |
104 | Double_t ratioupsilonPP; | |
b44c3901 | 105 | Double_t ratioccbar; |
106 | Double_t ratiobbbar; | |
107 | ||
18bc0c8e | 108 | // Cross sections in barns (from PPR Vol. II p: 552) pp - 14 TeV and |
109 | // corrected from feed down of higher resonances | |
103ac317 | 110 | |
9ff768ee | 111 | Double_t sigmajpsi = fSigmaJPsi; |
112 | Double_t sigmapsiP = fSigmaPsiP; | |
113 | Double_t sigmaupsilon = fSigmaUpsilon; | |
114 | Double_t sigmaupsilonP = fSigmaUpsilonP; | |
115 | Double_t sigmaupsilonPP = fSigmaUpsilonPP; | |
116 | Double_t sigmaccbar = fSigmaCCbar; | |
117 | Double_t sigmabbbar = fSigmaBBbar; | |
3285b419 | 118 | |
aaa95f22 | 119 | AliInfo(Form("the parametrised resonances uses the decay mode %d",fDecayModeResonance)); |
120 | ||
18bc0c8e | 121 | // Generation using CDF scaled distribution for pp@14TeV (from D. Stocco) |
122 | //---------------------------------------------------------------------- | |
103ac317 | 123 | // Jpsi generator |
18bc0c8e | 124 | AliGenParam * genjpsi = new AliGenParam(1, AliGenMUONlib::kJpsi, "CDF pp", "Jpsi"); |
103ac317 | 125 | // first step: generation in 4pi |
126 | genjpsi->SetPtRange(0.,100.); | |
127 | genjpsi->SetYRange(-8.,8.); | |
128 | genjpsi->SetPhiRange(0.,360.); | |
aaa95f22 | 129 | genjpsi->SetForceDecay(fDecayModeResonance); |
130 | if (!gMC) genjpsi->SetDecayer(fDecayer); | |
131 | ||
103ac317 | 132 | genjpsi->Init(); // generation in 4pi |
133 | ratiojpsi = sigmajpsi / sigmaReaction * genjpsi->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); // get weight | |
9ff768ee | 134 | if (!fSigmaSilent) { |
3285b419 | 135 | AliInfo(Form("jpsi prod. cross-section in pp(14 TeV) %5.3g b",sigmajpsi)); |
136 | AliInfo(Form("jpsi prod. probability per collision in acceptance %5.3g",ratiojpsi)); | |
137 | } | |
103ac317 | 138 | // second step: generation in selected kinematical range |
139 | genjpsi->SetPtRange(ptMin, ptMax); | |
140 | genjpsi->SetYRange(yMin, yMax); | |
141 | genjpsi->SetPhiRange(phiMin, phiMax); | |
142 | genjpsi->Init(); // generation in selected kinematic range | |
143 | AddGenerator(genjpsi, "Jpsi", ratiojpsi); // Adding Generator | |
103ac317 | 144 | //------------------------------------------------------------------ |
145 | // Psi prime generator | |
18bc0c8e | 146 | AliGenParam * genpsiP = new AliGenParam(1, AliGenMUONlib::kPsiP, "CDF pp", "PsiP"); |
103ac317 | 147 | // first step: generation in 4pi |
148 | genpsiP->SetPtRange(0.,100.); | |
149 | genpsiP->SetYRange(-8.,8.); | |
150 | genpsiP->SetPhiRange(0.,360.); | |
aaa95f22 | 151 | genpsiP->SetForceDecay(fDecayModeResonance); |
152 | if (!gMC) genpsiP->SetDecayer(fDecayer); | |
153 | ||
103ac317 | 154 | genpsiP->Init(); // generation in 4pi |
155 | ratiopsiP = sigmapsiP / sigmaReaction * genpsiP->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); | |
9ff768ee | 156 | if (!fSigmaSilent) { |
3285b419 | 157 | AliInfo(Form("psiP prod. cross-section in pp(14 TeV) %5.3g b",sigmapsiP)); |
158 | AliInfo(Form("psiP prod. probability per collision in acceptance %5.3g",ratiopsiP)); | |
159 | } | |
103ac317 | 160 | // second step: generation in selected kinematical range |
161 | genpsiP->SetPtRange(ptMin, ptMax); | |
162 | genpsiP->SetYRange(yMin, yMax); | |
163 | genpsiP->SetPhiRange(phiMin, phiMax); | |
164 | genpsiP->Init(); // generation in selected kinematic range | |
165 | AddGenerator(genpsiP, "PsiP", ratiopsiP); // Adding Generator | |
103ac317 | 166 | //------------------------------------------------------------------ |
167 | // Upsilon 1S generator | |
18bc0c8e | 168 | AliGenParam * genupsilon = new AliGenParam(1, AliGenMUONlib::kUpsilon, "CDF pp", "Upsilon"); |
103ac317 | 169 | // first step: generation in 4pi |
170 | genupsilon->SetPtRange(0.,100.); | |
171 | genupsilon->SetYRange(-8.,8.); | |
172 | genupsilon->SetPhiRange(0.,360.); | |
aaa95f22 | 173 | genupsilon->SetForceDecay(fDecayModeResonance); |
174 | if (!gMC) genupsilon->SetDecayer(fDecayer); | |
103ac317 | 175 | genupsilon->Init(); // generation in 4pi |
176 | ratioupsilon = sigmaupsilon / sigmaReaction * genupsilon->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); | |
9ff768ee | 177 | if (!fSigmaSilent) { |
3285b419 | 178 | AliInfo(Form("Upsilon 1S prod. cross-section in pp(14 TeV) %5.3g b",sigmaupsilon)); |
179 | AliInfo(Form("Upsilon 1S prod. probability per collision in acceptance %5.3g",ratioupsilon)); | |
180 | } | |
103ac317 | 181 | // second step: generation in selected kinematical range |
182 | genupsilon->SetPtRange(ptMin, ptMax); | |
183 | genupsilon->SetYRange(yMin, yMax); | |
184 | genupsilon->SetPhiRange(phiMin, phiMax); | |
185 | genupsilon->Init(); // generation in selected kinematical range | |
186 | AddGenerator(genupsilon,"Upsilon", ratioupsilon); // Adding Generator | |
103ac317 | 187 | //------------------------------------------------------------------ |
188 | // Upsilon 2S generator | |
18bc0c8e | 189 | AliGenParam * genupsilonP = new AliGenParam(1, AliGenMUONlib::kUpsilonP, "CDF pp", "UpsilonP"); |
103ac317 | 190 | // first step: generation in 4pi |
191 | genupsilonP->SetPtRange(0.,100.); | |
192 | genupsilonP->SetYRange(-8.,8.); | |
193 | genupsilonP->SetPhiRange(0.,360.); | |
aaa95f22 | 194 | genupsilonP->SetForceDecay(fDecayModeResonance); |
195 | if (!gMC) genupsilonP->SetDecayer(fDecayer); | |
103ac317 | 196 | genupsilonP->Init(); // generation in 4pi |
197 | ratioupsilonP = sigmaupsilonP / sigmaReaction * genupsilonP->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); | |
9ff768ee | 198 | if (!fSigmaSilent) { |
3285b419 | 199 | AliInfo(Form("Upsilon 2S prod. cross-section in pp(14 TeV) %5.3g b",sigmaupsilonP)); |
200 | AliInfo(Form("Upsilon 2S prod. probability per collision in acceptance %5.3g",ratioupsilonP)); | |
201 | } | |
103ac317 | 202 | // second step: generation in the kinematical range |
203 | genupsilonP->SetPtRange(ptMin, ptMax); | |
204 | genupsilonP->SetYRange(yMin, yMax); | |
205 | genupsilonP->SetPhiRange(phiMin, phiMax); | |
206 | genupsilonP->Init(); // generation in selected kinematical range | |
207 | AddGenerator(genupsilonP,"UpsilonP", ratioupsilonP); // Adding Generator | |
103ac317 | 208 | |
209 | //------------------------------------------------------------------ | |
210 | // Upsilon 3S generator | |
18bc0c8e | 211 | AliGenParam * genupsilonPP = new AliGenParam(1, AliGenMUONlib::kUpsilonPP, "CDF pp", "UpsilonPP"); |
103ac317 | 212 | // first step: generation in 4pi |
213 | genupsilonPP->SetPtRange(0.,100.); | |
214 | genupsilonPP->SetYRange(-8.,8.); | |
215 | genupsilonPP->SetPhiRange(0.,360.); | |
aaa95f22 | 216 | genupsilonPP->SetForceDecay(fDecayModeResonance); |
217 | if (!gMC) genupsilonPP->SetDecayer(fDecayer); | |
103ac317 | 218 | genupsilonPP->Init(); // generation in 4pi |
219 | ratioupsilonPP = sigmaupsilonPP / sigmaReaction * genupsilonPP->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); | |
9ff768ee | 220 | if (!fSigmaSilent) { |
3285b419 | 221 | AliInfo(Form("Upsilon 3S prod. cross-section in pp(14 TeV) %5.3g b",sigmaupsilonPP)); |
222 | AliInfo(Form("Upsilon 3S prod. probability per collision in acceptance %5.3g",ratioupsilonPP)); | |
223 | } | |
103ac317 | 224 | // second step: generation in selected kinematical range |
225 | genupsilonPP->SetPtRange(ptMin, ptMax); | |
226 | genupsilonPP->SetYRange(yMin, yMax); | |
227 | genupsilonPP->SetPhiRange(phiMin, phiMax); | |
228 | genupsilonPP->Init(); // generation in selected kinematical range | |
229 | AddGenerator(genupsilonPP,"UpsilonPP", ratioupsilonPP); // Adding Generator | |
b44c3901 | 230 | |
231 | //------------------------------------------------------------------ | |
232 | // Generator of charm | |
233 | ||
234 | AliGenCorrHF *gencharm = new AliGenCorrHF(1, 4); | |
8058ead1 | 235 | gencharm->SetMomentumRange(0,9999); |
236 | gencharm->SetForceDecay(kAll); | |
237 | ratioccbar = sigmaccbar/sigmaReaction; | |
238 | if (!gMC) gencharm->SetDecayer(fDecayer); | |
239 | gencharm->Init(); | |
9ff768ee | 240 | if (!fSigmaSilent) { |
3285b419 | 241 | AliInfo(Form("c-cbar prod. cross-section in pp(14 TeV) %5.3g b",sigmaccbar)); |
242 | AliInfo(Form("c-cbar prod. probability per collision in acceptance %5.3g",ratioccbar)); | |
243 | } | |
8058ead1 | 244 | AddGenerator(gencharm,"CorrHFCharm",ratioccbar); |
b44c3901 | 245 | |
103ac317 | 246 | //------------------------------------------------------------------ |
b44c3901 | 247 | // Generator of beauty |
248 | ||
8058ead1 | 249 | AliGenCorrHF *genbeauty = new AliGenCorrHF(1, 5); |
250 | genbeauty->SetMomentumRange(0,9999); | |
251 | genbeauty->SetForceDecay(kAll); | |
252 | ratiobbbar = sigmabbbar/sigmaReaction; | |
253 | if (!gMC) genbeauty->SetDecayer(fDecayer); | |
254 | genbeauty->Init(); | |
9ff768ee | 255 | if (!fSigmaSilent) { |
3285b419 | 256 | AliInfo(Form("b-bbar prod. cross-section in pp(14 TeV) %5.3g b",sigmabbbar)); |
257 | AliInfo(Form("b-bbar prod. probability per collision in acceptance %5.3g",ratiobbbar)); | |
258 | } | |
8058ead1 | 259 | AddGenerator(genbeauty,"CorrHFBeauty",ratiobbbar); |
b44c3901 | 260 | |
8058ead1 | 261 | //------------------------------------------------------------------- |
103ac317 | 262 | // Pythia generator |
8058ead1 | 263 | // |
264 | // This has to go into the Config.C | |
265 | // | |
266 | // AliGenPythia *pythia = new AliGenPythia(1); | |
267 | // pythia->SetProcess(kPyMbMSEL1); | |
268 | // pythia->SetStrucFunc(kCTEQ5L); | |
269 | // pythia->SetEnergyCMS(14000.); | |
270 | // AliInfo(Form("\n\npythia uses the decay mode %d", GetDecayModePythia())); | |
271 | // Decay_t dt = gener->GetDecayModePythia(); | |
272 | // pythia->SetForceDecay(dt); | |
273 | // pythia->SetPtRange(0.,100.); | |
274 | // pythia->SetYRange(-8.,8.); | |
275 | // pythia->SetPhiRange(0.,360.); | |
276 | // pythia->SetPtHard(2.76,-1.0); | |
277 | // pythia->SwitchHFOff(); | |
278 | // pythia->Init(); | |
279 | // AddGenerator(pythia,"Pythia",1); | |
280 | } | |
103ac317 | 281 | |
8058ead1 | 282 | void AliGenMUONCocktailpp::Init() |
283 | { | |
284 | // | |
285 | // Initialisation | |
aaa95f22 | 286 | TIter next(fEntries); |
287 | AliGenCocktailEntry *entry; | |
288 | if (fStack) { | |
289 | while((entry = (AliGenCocktailEntry*)next())) { | |
290 | entry->Generator()->SetStack(fStack); | |
291 | } | |
292 | } | |
103ac317 | 293 | } |
294 | ||
295 | //_________________________________________________________________________ | |
296 | void AliGenMUONCocktailpp::Generate() | |
297 | { | |
298 | // Generate event | |
299 | TIter next(fEntries); | |
300 | AliGenCocktailEntry *entry = 0; | |
301 | AliGenCocktailEntry *preventry = 0; | |
302 | AliGenerator* gen = 0; | |
303 | ||
d94af0c1 | 304 | const TObjArray *partArray = gAlice->GetMCApp()->Particles(); |
103ac317 | 305 | |
306 | // Generate the vertex position used by all generators | |
307 | if(fVertexSmear == kPerEvent) Vertex(); | |
308 | ||
309 | // Loop on primordialTrigger: | |
310 | // minimum muon multiplicity above a pt cut in a theta acceptance region | |
311 | ||
312 | Bool_t primordialTrigger = kFALSE; | |
313 | while(!primordialTrigger) { | |
314 | //Reseting stack | |
33c3c91a | 315 | AliRunLoader * runloader = AliRunLoader::Instance(); |
103ac317 | 316 | if (runloader) |
317 | if (runloader->Stack()) | |
34da8494 | 318 | runloader->Stack()->Clean(); |
103ac317 | 319 | // Loop over generators and generate events |
320 | Int_t igen = 0; | |
321 | Int_t npart = 0; | |
322 | const char* genName = 0; | |
323 | while((entry = (AliGenCocktailEntry*)next())) { | |
324 | gen = entry->Generator(); | |
325 | genName = entry->GetName(); | |
326 | gen->SetVertex(fVertex.At(0), fVertex.At(1), fVertex.At(2)); | |
327 | ||
328 | npart = (strcmp(genName,"Pythia") == 0) ? 1 : | |
329 | gRandom->Poisson(entry->Rate()); | |
330 | ||
331 | if (npart > 0) { | |
332 | igen++; | |
333 | if (igen == 1) entry->SetFirst(0); | |
334 | else entry->SetFirst((partArray->GetEntriesFast())+1); | |
335 | ||
336 | gen->SetNumberParticles(npart); | |
337 | gen->Generate(); | |
338 | entry->SetLast(partArray->GetEntriesFast()); | |
339 | preventry = entry; | |
340 | } | |
341 | } | |
342 | next.Reset(); | |
343 | ||
344 | // Testing primordial trigger: Single muons or dimuons with Pt above a Pt cut | |
345 | // in the muon spectrometer acceptance | |
346 | Int_t iPart; | |
347 | fNGenerated++; | |
b44c3901 | 348 | Int_t numberOfMuons=0;Int_t maxPart = partArray->GetEntriesFast(); |
aaa95f22 | 349 | for(iPart=0; iPart<maxPart; iPart++){ |
103ac317 | 350 | |
aaa95f22 | 351 | TParticle *part = gAlice->GetMCApp()->Particle(iPart); |
352 | if ( TMath::Abs(part->GetPdgCode()) == 13 ){ | |
353 | if((part->Vz() > fMuonOriginCut) && //take only the muons that decayed before the abs + 1 int. length in C abs | |
354 | (part->Theta()*180./TMath::Pi()>fMuonThetaMinCut) && | |
355 | (part->Theta()*180./TMath::Pi()<fMuonThetaMaxCut) && | |
f81a4aca | 356 | (part->Pt()>fMuonPtCut) && |
357 | (part->P()>fMuonPCut)) { | |
aaa95f22 | 358 | numberOfMuons++; |
103ac317 | 359 | } |
aaa95f22 | 360 | } |
361 | } | |
103ac317 | 362 | if (numberOfMuons >= fMuonMultiplicity) primordialTrigger = kTRUE; |
363 | } | |
364 | fNSucceded++; | |
365 | ||
aaa95f22 | 366 | // AliInfo(Form("Generated Events are %d and Succeeded Events are %d",fNGenerated,fNSucceded)); |
103ac317 | 367 | AliDebug(5,Form("Generated Events are %d and Succeeded Events are %d",fNGenerated,fNSucceded)); |
368 | } | |
369 | ||
370 |