]> git.uio.no Git - u/mrichter/AliRoot.git/blame - test/genkine/gen/fastMcProduction.C
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Corrected library names and paths to macros
[u/mrichter/AliRoot.git] / test / genkine / gen / fastMcProduction.C
CommitLineData
b41cb7c3 1/***************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: ALICE OFFLICE. *
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 $Satyajit Jena || alien:sjena Sun Apr 21 14:05:19 CEST 2013$
17
18 Sterring macro for fast production of MC events - Followed from
19 the original macro of fastGenAmpt.C
20
21 Implemented Generators: (Version 1.0: Sun Apr 21 14:05:19 CEST 2013)
22 ----------------------------------------------------------------
23 kPythia6, kPythia8, kPythia6D6T,
24 kPythiaPerugia0, kPythia6ATLAS,
25 kPythiaJets,
26 kPhojet,
27 kDPMjet, kDPMjet_pA,
28 kAmptDefault, kAmptStringMelting, kAmptStringMeltingNoART,
29 kAmptpA, kAmptFlow,
30 kAmptReducedFlow,
31
32 FIXME:
33 kPythia6ATLAS_Flat,
34 kHijing, kHijing2000, kHijing_pA,
35
36
37 ***************************************************************************/
38
39#if !defined(__CINT__) || defined(__MAKECINT__)
40#include <Riostream.h>
41#include <TStopwatch.h>
42#include <TDatime.h>
43#include <TRandom.h>
44#include <TFile.h>
45#include <TTree.h>
46#include <TMath.h>
47#include <TParticle.h>
48#include <TPDGCode.h>
49#include <TDatabasePDG.h>
50#include <TRandom3.h>
51#include <TChain.h>
52#include "AliRun.h"
53#include "AliRunLoader.h"
54#include "AliHeader.h"
55#include "AliStack.h"
56#include "AliPDG.h"
57#include "AliGenAmpt.h"
58#include "TAmpt.h"
59#include "PYTHIA6/AliDecayerPythia.h"
60#include <TVirtualMC.h>
61#include <TF1.h>
62#include "STEER/AliConfig.h"
63#include "STEER/AliGenerator.h"
64#include "STEER/AliLog.h"
65#include "EVGEN/AliGenHIJINGpara.h"
66#include "THijing/AliGenHijing.h"
67#include "EVGEN/AliGenCocktail.h"
68#include "EVGEN/AliGenSlowNucleons.h"
69#include "EVGEN/AliSlowNucleonModelExp.h"
70#include "EVGEN/AliGenParam.h"
71#include "EVGEN/AliGenMUONlib.h"
72#include "EVGEN/AliGenSTRANGElib.h"
73#include "EVGEN/AliGenMUONCocktail.h"
74#include "EVGEN/AliGenCocktail.h"
75#include "EVGEN/AliGenGeVSim.h"
76#include "EVGEN/AliGeVSimParticle.h"
77#include "PYTHIA6/AliGenPythia.h"
78#endif
79
80//___________________________________________________________________________
81enum PDC06Proc_t {
82 kPythia6,
83 kPythia8,
84 kPythia6D6T,
85 kPythiaPerugia0,
86 kPythia6ATLAS,
87 kPythia6ATLAS_Flat,
88 kPythiaJets,
89 kPhojet,
90 kHijing,
91 kHijing2000,
92 kHijing_pA,
93 kDPMjet,
94 kDPMjet_pA,
95 kAmptDefault,
96 kAmpt,
97 kAmptpA,
98 kAmptFlowStringMelting,
99 kAmptStringMeltingNoART,
100 kAmptFlow,
101 kAmptReducedFlow,
102 kRunMax
103};
104
105//___________________________________________________________________________
106const char * pprRunName[] = {
107 "kPythia6",
108 "kPythia8",
109 "kPythia6D6T",
110 "kPythiaPerugia0",
111 "kPythia6ATLAS",
112 "kPythia6ATLAS_Flat",
113 "kPythiaJets",
114 "kPhojet",
115 "kHijing",
116 "kHijing2000",
117 "kHijing_pA",
118 "kDPMjet",
119 "kDPMjet_pA",
120 "kAmpt",
121 "kAmptpA",
122 "kAmptFlow",
123 "kAmptReducedFlow"
124};
125
126enum PprTrigConf_t {kDefaultPPTrig, kDefaultPbPbTrig };
127const char * pprTrigConfName[] = {"p-p","Pb-Pb"};
128
129//___________________________________________________________________________
130void ProcessEnvironmentVars();
131class AliGenPythia;
132AliGenerator *MbPythia();
133AliGenerator *Pythia8();
134AliGenerator *MbPythiaTuneD6T();
135AliGenerator *MbPythiaTunePerugia0();
136AliGenerator *MbPythiaTuneATLAS();
137AliGenerator *MbPythiaTuneATLAS_Flat();
138AliGenerator *PythiaJets();
139AliGenerator *MbPhojet();
140AliGenerator *Hijing();
141AliGenerator *Hijing2000();
142AliGenerator *Hijing_pA(Bool_t kSlowN);
143AliGenerator *DPMjet();
144AliGenerator *DPMjet_pA(Bool_t fragments);
145AliGenerator *Ampt();
146AliGenerator *AmptpA();
147AliGenerator* AmptFlow();
148AliGenerator *AmptReducedFlow();
149
150
151//_________________________________________________________________________
152//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
153// Geterator, field, beam energy
154
155//static Double_t pBeamEnergy = 4000.0; // Used during pA runs
156//static Double_t energy = 2.*pBeamEnergy*TMath::Sqrt(82./208.); //energy in CMS
157
158static Double_t energy = 7000.;
159static PDC06Proc_t proc = kPythia6;
160static Float_t bMin = 0.;
161static Float_t bMax = 100.;
162static PprTrigConf_t strig = kDefaultPPTrig; // default pp trigger configuration
163
164
165static Double_t JpsiPol = 0; // Jpsi polarisation
166static Bool_t JpsiHarderPt = kFALSE; // Jpsi harder pt spectrum (8.8 TeV)
167static TString comment;
168//static PprTrigConf_t strig = kDefaultPbPbTrig; // default pp trigger configuration
169TDatime dt;
170static UInt_t seed = dt.Get();
171
172//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
173//_________________________________________________________________________
174void fastMcProduction(Int_t nev = 300) {
175
176 ProcessEnvironmentVars();
177
178 gRandom->SetSeed(seed);
179 cerr<<" Seed for random number generation = "<<seed<<endl;
180
181#if defined(__CINT__)
182 gSystem->Load("liblhapdf");
183 gSystem->Load("libEGPythia6");
184
185 if (proc == kPythia6 || proc == kPhojet || proc == kDPMjet || proc==kDPMjet_pA) {
186 gSystem->Load("libpythia6"); // Pythia 6.2
187 gSystem->Load("libAliPythia6"); // ALICE specific implementations
188 }
189
190 if (proc == kHijing || proc == kHijing2000 || proc == kHijing_pA ) {
257e1b2b 191 gSystem->Load("libHIJING");
b41cb7c3 192 gSystem->Load("libTHijing");
193 }
194
195 else if ( proc == kDPMjet || proc== kDPMjet_pA ) {
257e1b2b 196 gSystem->Load("libDPMJET");
b41cb7c3 197 gSystem->Load("libTDPMjet");
198 }
199
200 else if (proc == kAmptDefault || kAmptFlowStringMelting || proc == kAmptStringMeltingNoART || proc == kAmptpA || proc == kAmptReducedFlow) {
201 gSystem->Load("libampt");
202 gSystem->Load("libTAmpt");
203 gSystem->Load("libpythia6.so");
204 gSystem->Load("libAliPythia6.so");
205 }
206
207 if (proc == kPythia8) {
208 gSystem->Load("libpythia8.so");
209 gSystem->Load("libAliPythia8.so");
210 gSystem->Setenv("PYTHIA8DATA", gSystem->ExpandPathName("$ALICE_ROOT/PYTHIA8/pythia8145/xmldoc"));
211 gSystem->Setenv("LHAPDF", gSystem->ExpandPathName("$ALICE_ROOT/LHAPDF"));
212 gSystem->Setenv("LHAPATH", gSystem->ExpandPathName("$ALICE_ROOT/LHAPDF/PDFsets"));
213 }
214#endif
215
216
217 AliGenerator* gener = 0x0;
218
219 cout<<"Run type set to ------------- "<<pprRunName[proc]<<" " << proc << " " << kDPMjet_pA<< endl;
220
221 if (proc == kPythia6) {
222 Printf("<<<<<<<<<<<<<<<<<<< Processing Mb. PYTHIA >>>>>>>>>>>>>>>>>>>>");
223 gener = MbPythia();
224 }
225
226 else if (proc == kPythia8) {
227 Printf("<<<<<<<<<<<<<<<<<<< Processing Mb. Pythia 8 >>>>>>>>>>>>>>>>>>>>");
228 gener = Pythia8();
229 }
230
231 else if (proc == kPythia6D6T) {
232 Printf("<<<<<<<<<<<<<<<<<<< Processing Mb. PYTHIA D6T >>>>>>>>>>>>>>>>>>>>");
233 gener = MbPythiaTuneD6T();
234 }
235
236 else if (proc == kPythiaPerugia0) {
237 Printf("<<<<<<<<<<<<<<<<<<< Processing Mb. PYTHIA Perugia0 >>>>>>>>>>>>>>>>>>>>");
238 gener = MbPythiaTunePerugia0();
239 }
240
241 else if (proc == kPythia6ATLAS) {
242 Printf("<<<<<<<<<<<<<<<<<<< Processing Mb. PYTHIA ATLAS>>>>>>>>>>>>>>>>>>>>");
243 gener = MbPythiaTuneATLAS();
244 }
245
246 else if (proc == kPythia6ATLAS_Flat) {
247 Printf("<<<<<<<<<<<<<<<<<<< Processing Mb. PYTHIA ATLAS_FLAT >>>>>>>>>>>>>>>>>>>>");
248 gener = MbPythiaTuneATLAS_Flat();
249 }
250
251 else if (proc == kPythiaJets ) {
252 Printf("<<<<<<<<<<<<<<<<<<< Processing Pythia Jets >>>>>>>>>>>>>>>>>>>>");
253 gener = PythiaJets();
254 }
255
256 else if (proc == kPhojet) {
257 Printf("<<<<<<<<<<<<<<<<<<< Processing Mb. PHOJET >>>>>>>>>>>>>>>>>>>>");
258 gener = MbPhojet();
259 }
260
261 else if (proc == kHijing) {
262 Printf("<<<<<<<<<<<<<<<<<<< Processing Mb. HIJING >>>>>>>>>>>>>>>>>>>>");
263 gener = Hijing();
264 }
265
266 else if (proc == kHijing2000) {
267 Printf("<<<<<<<<<<<<<<<<<<< Processing Mb. HIJING 2000 >>>>>>>>>>>>>>>>>>>>");
268 gener = Hijing2000();
269 }
270
271 else if (proc ==kHijing_pA) {
272 Printf("<<<<<<<<<<<<<<<<<<< Processing Mb. pA Hijing >>>>>>>>>>>>>>>>>>>>");
273 gener = Hijing_pA(kTRUE);
274 }
275
276 else if (proc == kDPMjet) {
277 Printf("<<<<<<<<<<<<<<<<<<< Processing DMPJet >>>>>>>>>>>>>>>>>>>>");
278 gener = DPMjet();
279 }
280
281 else if (proc == kDPMjet_pA) {
282 Printf("<<<<<<<<<<<<<<<<<<< Processing DMPJet pA >>>>>>>>>>>>>>>>>>>>");
283 gener = DPMjet_pA(kFALSE);
284 }
285
286 else if (proc == kAmptDefault) {
287 Printf("<<<<<<<<<<<<<<<<<<< Processing AMPT Default >>>>>>>>>>>>>>>>>>>>");
288 gener = AmptDefault();
289 }
290
291 else if (proc == kAmptStringMelting) {
292 Printf("<<<<<<<<<<<<<<<<<<< Processing AMPT With Flow >>>>>>>>>>>>>>>>>>>>");
293 gener = AmptStringMelting();
294 }
295
296 else if (proc == kAmptStringMeltingNoART) {
297 Printf("<<<<<<<<<<<<<<<<<<< Processing AMPT With Flow >>>>>>>>>>>>>>>>>>>>");
298 gener = AmptStringMeltingNoART();
299 }
300
301 else if (proc == kAmptpA) {
302 Printf("<<<<<<<<<<<<<<<<<<< Processing AMPT pA >>>>>>>>>>>>>>>>>>>>");
303 gener = AmptpA();
304 }
305
306
307 else if (proc == kAmptReducedFlow) {
308 // Specific Fastgen
309 Printf("<<<<<<<<<<<<<<<<<<< Processing AMPT With Reduced Flow >>>>>>>>>>>>>>>>>>>>");
310 gener = AmptReducedFlow();
311 }
312
313 else {
314 cout << "ERROR : Wrong Procss Selcted !!!" << endl;
315 return;
316 }
317
318
319 AliPDG::AddParticlesToPdgDataBase();
320 TDatabasePDG::Instance();
321
322 const char* filename = "galice.root";
323 AliRunLoader* rl = AliRunLoader::Open(filename,"FASTRUN","recreate");
324
325 rl->SetCompressionLevel(2);
326 rl->SetNumberOfEventsPerFile(nev);
327 rl->LoadKinematics("RECREATE");
328 rl->MakeTree("E");
329 gAlice->SetRunLoader(rl);
330 rl->MakeStack();
331 AliStack* stack = rl->Stack();
332
333 AliHeader* header = rl->GetHeader();
334
335 /*
336 Float_t sigmaz = 5.4 / TMath::Sqrt(2.); // [cm]
337 Float_t betast = 3.5; // beta* [m]
338 Float_t eps = 3.75e-6; // emittance [m]
339 Float_t gamma = energy / 2.0 / 0.938272; // relativistic gamma [1]
340 Float_t sigmaxy = TMath::Sqrt(eps * betast / gamma) / TMath::Sqrt(2.) * 100.; // [cm]
341 printf("\n \n Diamond size x-y: %10.3e z: %10.3e\n \n", sigmaxy, sigmaz);
342 gener->SetSigma(sigmaxy, sigmaxy, sigmaz); // Sigma in (X,Y,Z) (cm) on IP position
343 gener->SetVertexSmear(kPerEvent);
344 */
345
346
347 gener->Init();
348 gener->SetStack(stack);
349
350 rl->CdGAFile();
351
352 TStopwatch timer;
353 timer.Start();
354 for (Int_t iev = 0; iev < nev; iev++) {
355 cout <<"============================================= Event number "<< iev << endl;
356 // Initialize event
357 header->Reset(0,iev);
358 rl->SetEventNumber(iev);
359 stack->Reset();
360 rl->MakeTree("K");
361 Int_t nprim = 0;
362 Int_t ntrial = 0;
363 // Int_t ndstar = 0;
364 stack->Reset();
365 stack->ConnectTree(rl->TreeK());
366 gener->Generate();
367 ntrial++;
368 nprim = stack->GetNprimary();
369 cout << "Number of particles " << nprim << endl;
370 cout << "Number of trials " << ntrial << endl;
371 header->SetNprimary(stack->GetNprimary());
372 header->SetNtrack(stack->GetNtrack());
373 stack->FinishEvent();
374 header->SetStack(stack);
375 rl->TreeE()->Fill();
376 rl->WriteKinematics("OVERWRITE");
377 } // event loop
378 timer.Stop();
379 timer.Print();
380 gener->FinishRun();
381 rl->WriteHeader("OVERWRITE");
382 gener->Write();
383 rl->Write();
384}
385
386//___________________________________________________//
387void ProcessEnvironmentVars() {
388 // Run type
389 if (gSystem->Getenv("CONFIG_RUN_TYPE")) {
390 for (Int_t iRun = 0; iRun < kRunMax; iRun++) {
391 if (strcmp(gSystem->Getenv("CONFIG_RUN_TYPE"), pprRunName[iRun]) == 0) {
392 proc = (PDC06Proc_t)iRun;
393 cout<<"Run type set to "<<pprRunName[iRun]<<endl;
394 }
395 }
396 }
397
398
399 // Energy
400 if (gSystem->Getenv("CONFIG_ENERGY")) {
401 energy = atoi(gSystem->Getenv("CONFIG_ENERGY"));
402 cout<<"Energy set to "<<energy<<" GeV"<<endl;
403 }
404
405 // Random Number seed
406 if (gSystem->Getenv("CONFIG_SEED")) {
407 seed = atoi(gSystem->Getenv("CONFIG_SEED"));
408 }
409
410 // Impact param
411 if (gSystem->Getenv("CONFIG_BMIN")) {
412 bMin = atof(gSystem->Getenv("CONFIG_BMIN"));
413 }
414
415 if (gSystem->Getenv("CONFIG_BMAX")) {
416 bMax = atof(gSystem->Getenv("CONFIG_BMAX"));
417 }
418 cout<<"Impact parameter in ["<<bMin<<","<<bMax<<"]"<<endl;
419}
420
421
422
423//______________________________________________________________________
424AliGenerator* MbPythia() // Mb Pythia
425{
426 comment = comment.Append(" pp: Pythia low-pt");
427 AliGenPythia* pythia = new AliGenPythia(-1);
428 /* pythia->SetMomentumRange(0, 999999.);
429 pythia->SetThetaRange(0., 180.);
430 pythia->SetYRange(-12.,12.);
431 pythia->SetPtRange(0,1000.);*/
432 pythia->SetProcess(kPyMb);
433 pythia->SetEnergyCMS(energy);
434
435 return pythia;
436}
437
438
439//______________________________________________________________________
440AliGenerator* Pythia8()
441{
442 AliGenPythiaPlus* gener = new AliGenPythiaPlus(AliPythia8::Instance());
443 gener->SetProcess(kPyMbDefault);
444 gener->SetEnergyCMS(energy);
445 gener->SetEventListRange(-1, 2);
446 return gener;
447}
448
449
450
451//______________________________________________________________________
452AliGenerator* MbPythiaTuneD6T()
453{
454 comment = comment.Append(" pp: Pythia low-pt");
455 AliGenPythia* pythia = new AliGenPythia(-1);
456 pythia->SetMomentumRange(0, 999999.);
457 pythia->SetThetaRange(0., 180.);
458 pythia->SetYRange(-12.,12.);
459 pythia->SetPtRange(0,1000.);
460 pythia->SetProcess(kPyMb);
461 pythia->SetEnergyCMS(energy);
462// Tune
463// 109 D6T : Rick Field's CDF Tune D6T (NB: needs CTEQ6L pdfs externally)
464 pythia->SetTune(109); // F I X
465 pythia->SetStrucFunc(kCTEQ6l);
466//
467 return pythia;
468}
469
470//______________________________________________________________________
471AliGenerator* MbPythiaTunePerugia0()
472{
473 comment = comment.Append(" pp: Pythia low-pt (Perugia0)");
474//
475// Pythia
476 AliGenPythia* pythia = new AliGenPythia(-1);
477 /* pythia->SetMomentumRange(0, 999999.);
478 pythia->SetThetaRange(0., 180.);
479 pythia->SetYRange(-12.,12.);
480 pythia->SetPtRange(0,1000.);*/
481 pythia->SetProcess(kPyMb);
482 pythia->SetEnergyCMS(energy);
483// Tune
484// 320 Perugia 0
485 pythia->SetTune(320);
486 pythia->UseNewMultipleInteractionsScenario();
487//
488 return pythia;
489}
490
491//______________________________________________________________________
492AliGenerator* MbPythiaTuneATLAS()
493{
494 comment = comment.Append(" pp: Pythia low-pt");
495//
496// Pythia
497 AliGenPythia* pythia = new AliGenPythia(-1);
498 /* pythia->SetMomentumRange(0, 999999.);
499 pythia->SetThetaRange(0., 180.);
500 pythia->SetYRange(-12.,12.);
501 pythia->SetPtRange(0,1000.);*/
502 pythia->SetProcess(kPyMb);
503 pythia->SetEnergyCMS(energy);
504// Tune
505// C 306 ATLAS-CSC: Arthur Moraes' (new) ATLAS tune (needs CTEQ6L externally)
506 pythia->SetTune(306);
507 pythia->SetStrucFunc(kCTEQ6l);
508//
509 return pythia;
510}
511
512
513//______________________________________________________________________
514AliGenerator* MbPythiaTuneATLAS_Flat()
515{
516 AliGenPythia* pythia = MbPythiaTuneATLAS();
517
518 comment = comment.Append("; flat multiplicity distribution");
519
520 // set high multiplicity trigger
521 // this weight achieves a flat multiplicity distribution
522 TH1 *weight = new TH1D("weight","weight",201,-0.5,200.5);
523 weight->SetBinContent(1,5.49443);
524 weight->SetBinContent(2,8.770816);
525 weight->SetBinContent(6,0.4568624);
526 weight->SetBinContent(7,0.2919915);
527 weight->SetBinContent(8,0.6674189);
528 weight->SetBinContent(9,0.364737);
529 weight->SetBinContent(10,0.8818444);
530 weight->SetBinContent(11,0.531885);
531 weight->SetBinContent(12,1.035197);
532 weight->SetBinContent(13,0.9394057);
533 weight->SetBinContent(14,0.9643193);
534 weight->SetBinContent(15,0.94543);
535 weight->SetBinContent(16,0.9426507);
536 weight->SetBinContent(17,0.9423649);
537 weight->SetBinContent(18,0.789456);
538 weight->SetBinContent(19,1.149026);
539 weight->SetBinContent(20,1.100491);
540 weight->SetBinContent(21,0.6350525);
541 weight->SetBinContent(22,1.351941);
542 weight->SetBinContent(23,0.03233504);
543 weight->SetBinContent(24,0.9574557);
544 weight->SetBinContent(25,0.868133);
545 weight->SetBinContent(26,1.030998);
546 weight->SetBinContent(27,1.08897);
547 weight->SetBinContent(28,1.251382);
548 weight->SetBinContent(29,0.1391099);
549 weight->SetBinContent(30,1.192876);
550 weight->SetBinContent(31,0.448944);
551 weight->SetBinContent(32,1);
552 weight->SetBinContent(33,1);
553 weight->SetBinContent(34,1);
554 weight->SetBinContent(35,1);
555 weight->SetBinContent(36,0.9999997);
556 weight->SetBinContent(37,0.9999997);
557 weight->SetBinContent(38,0.9999996);
558 weight->SetBinContent(39,0.9999996);
559 weight->SetBinContent(40,0.9999995);
560 weight->SetBinContent(41,0.9999993);
561 weight->SetBinContent(42,1);
562 weight->SetBinContent(43,1);
563 weight->SetBinContent(44,1);
564 weight->SetBinContent(45,1);
565 weight->SetBinContent(46,1);
566 weight->SetBinContent(47,0.9999999);
567 weight->SetBinContent(48,0.9999998);
568 weight->SetBinContent(49,0.9999998);
569 weight->SetBinContent(50,0.9999999);
570 weight->SetBinContent(51,0.9999999);
571 weight->SetBinContent(52,0.9999999);
572 weight->SetBinContent(53,0.9999999);
573 weight->SetBinContent(54,0.9999998);
574 weight->SetBinContent(55,0.9999998);
575 weight->SetBinContent(56,0.9999998);
576 weight->SetBinContent(57,0.9999997);
577 weight->SetBinContent(58,0.9999996);
578 weight->SetBinContent(59,0.9999995);
579 weight->SetBinContent(60,1);
580 weight->SetBinContent(61,1);
581 weight->SetBinContent(62,1);
582 weight->SetBinContent(63,1);
583 weight->SetBinContent(64,1);
584 weight->SetBinContent(65,0.9999999);
585 weight->SetBinContent(66,0.9999998);
586 weight->SetBinContent(67,0.9999998);
587 weight->SetBinContent(68,0.9999999);
588 weight->SetBinContent(69,1);
589 weight->SetBinContent(70,1);
590 weight->SetBinContent(71,0.9999997);
591 weight->SetBinContent(72,0.9999995);
592 weight->SetBinContent(73,0.9999994);
593 weight->SetBinContent(74,1);
594 weight->SetBinContent(75,1);
595 weight->SetBinContent(76,1);
596 weight->SetBinContent(77,1);
597 weight->SetBinContent(78,0.9999999);
598 weight->SetBinContent(79,1);
599 weight->SetBinContent(80,1);
600 weight->SetEntries(526);
601 Int_t limit = weight->GetRandom();
602 pythia->SetTriggerChargedMultiplicity(limit, 1.4);
603 comment = comment.Append(Form("; multiplicity threshold set to %d in |eta| < 1.4", limit));
604 return pythia;
605}
606
607//______________________________________________________________________
608AliGenerator* PythiaJets()
609{
610 comment = comment.Append(" pp: Pythia low-pt");
611//
612// Pythia
613 AliGenPythia* pythia = new AliGenPythia(-1);
614 /* pythia->SetMomentumRange(0, 999999.);
615 pythia->SetThetaRange(0., 180.);
616 pythia->SetYRange(-12.,12.);
617 pythia->SetPtRange(0,1000.);*/
618 pythia->SetProcess(kPyJets);
619 pythia->SetEnergyCMS(energy);
620 pythia->SetStrucFunc(kCTEQ6l);
621 // pythia->SetPtHard(50, 1000);
622//
623 return pythia;
624}
625
626
627
628//___________________________________________________________________
629AliGenerator* MbPhojet()
630{
631 comment = comment.Append(" pp: Pythia low-pt");
632#if defined(__CINT__)
257e1b2b 633 gSystem->Load("libDPMJET"); // Parton density functions
b41cb7c3 634 gSystem->Load("libTDPMjet"); // Parton density functions
635#endif
636 AliGenDPMjet* dpmjet = new AliGenDPMjet(-1);
637
638 /* dpmjet->SetMomentumRange(0, 999999.);
639 dpmjet->SetThetaRange(0., 180.);
640 dpmjet->SetYRange(-12.,12.);
641 dpmjet->SetPtRange(0,1000.);*/
642
643 dpmjet->SetProcess(kDpmMb);
644 dpmjet->SetEnergyCMS(energy);
645 return dpmjet;
646}
647
648//__________________________________________________________________
649AliGenerator* Hijing()
650{
651 AliGenHijing *gener = new AliGenHijing(-1);
652 // centre of mass energy
653 gener->SetEnergyCMS(energy);
654 gener->SetImpactParameterRange(bMin, bMax);
655 // reference frame
656 gener->SetReferenceFrame("CMS");
657 // projectile
658 gener->SetProjectile("A", 208, 82);
659 gener->SetTarget ("A", 208, 82);
660 // tell hijing to keep the full parent child chain
661 gener->KeepFullEvent();
662 // enable jet quenching
663 gener->SetJetQuenching(1);
664 // enable shadowing
665 gener->SetDecaysOff(0);
666 gener->SetShadowing(1);
667 // Don't track spectators
668 gener->SetSpectators(0);
669 // kinematic selection
670 gener->SetSelectAll(0);
671 return gener;
672}
673
674
675//__________________________________________________________________
676AliGenerator* Hijing2000()
677{
678 AliGenHijing *gener = (AliGenHijing*) Hijing();
679 gener->SetJetQuenching(0);
680 gener->SetPtHardMin (2.3);
681 return gener;
682}
683
684
685
686//_____________________________________________________________
687AliGenerator* Hijing_pA(Bool_t kSlowN) {
688 AliGenCocktail *gener = 0x0;
689 if (kSlowN) {
690 gener = new AliGenCocktail();
691 gener->SetProjectile("A", 208, 82);
692 gener->SetTarget ("P", 1, 1);
693 gener->SetEnergyCMS(energy);
694 }
695
696 AliGenHijing *hijing = new AliGenHijing(-1);
697 // centre of mass energy
698 hijing->SetEnergyCMS(energy);
699 // impact parameter range
700 hijing->SetImpactParameterRange(0., 100.);
701 // reference frame
702 hijing->SetReferenceFrame("CMS");
703 hijing->SetBoostLHC(1);
704 // projectile
705 hijing->SetProjectile("A", 208, 82);
706 hijing->SetTarget ("P", 1, 1);
707 // tell hijing to keep the full parent child chain
708 hijing->KeepFullEvent();
709 // enable jet quenching
710 hijing->SetJetQuenching(0);
711 // enable shadowing
712 hijing->SetShadowing(1);
713 // kinematic selection
714 hijing->SetSelectAll(0);
715
716 if (!kSlowN) {
717 // DO track spectators
718 hijing->SetSpectators(1);
719 return hijing;
720 }
721 else {
722 // Cocktail with slow nucleons generator
723 // DO NOT track spectators
724 hijing->SetSpectators(0);
725 AliGenSlowNucleons* gray = new AliGenSlowNucleons(1);
726 AliCollisionGeometry* coll = hijing->CollisionGeometry();
727 AliSlowNucleonModelExp* model = new AliSlowNucleonModelExp();
728 // Not yet in the release...
729 // model->SetSaturation(kTRUE);
730 gray->SetSlowNucleonModel(model);
731 gray->SetTarget(208,82);
732 gray->SetThetaDist(1);
733 gray->SetProtonDirection(1);
734 // gray->SetDebug(1);
735 gray->SetNominalCmsEnergy(2.*pBeamEnergy);
736 gray->NeedsCollisionGeometry();
737 gray->SetCollisionGeometry(coll);
738
739 gener->AddGenerator(hijing, "Hijing pPb", 1);
740 gener->AddGenerator(gray, "Gray Particles", 1);
741
742 return gener;
743 }
744}
745
746
747
748//__________________________________________________________________
749AliGenerator* DPMjet()
750{
751 AliGenDPMjet* dpmjet = new AliGenDPMjet(-1);
752 dpmjet->SetEnergyCMS(energy);
753 dpmjet->SetProjectile("A", 208, 82);
754 dpmjet->SetTarget ("A", 208, 82);
755 dpmjet->SetImpactParameterRange(bMin, bMax);
756 dpmjet->SetPi0Decay(0);
757 return dpmjet;
758}
759
760
761//____________________________________________
762AliGenerator* DPMjet_pA(Bool_t fragments)
763{
764 AliGenDPMjet *gener = new AliGenDPMjet(-1);
765 // A-p
766 gener->SetProjectile("A", 208, 82);
767 gener->SetTarget("P", 1, 1);
768 //
769 gener->SetEnergyCMS(energy);
770 gener->SetProjectileBeamEnergy(82.*pBeamEnergy/208.);
771
772 //
773 gener->SetProcess(kDpmMb);
774 gener->SetImpactParameterRange(0., 100.);
775 // DO NOT BOOST !
776 // gener->SetBoostLHC(1);
777 //
778 gener->SetFragmentProd(fragments);
779
780 return gener;
781
782}
783
784//__________________________________________________________________
785AliGenerator* AmptDefault()
786{
787 AliGenAmpt *genHi = new AliGenAmpt(-1);
788 //=========================================================================
789 // THE DECAYER
790 AliDecayer *decayer = new AliDecayerPythia();
791 cout << "*****************************************" << endl;
792 genHi->SetForceDecay( kHadronicD );
793 genHi->SetDecayer( decayer );
794 //=========================================================================
795 genHi->SetEnergyCMS(energy);
796 genHi->SetReferenceFrame("CMS");
797 genHi->SetProjectile("A",208,82);
798 genHi->SetTarget("A",208,82);
799
800 genHi->SetIsoft(1); //1: defaul - 4: string melting
801 genHi->SetStringFrag(0.5,0.9); //Lund string frgmentation parameters
802 genHi->SetPtHardMin (3);
803 //genHi->SetImpactParameterRange(9.,9.5);
804 genHi->SetImpactParameterRange(bMin,bMax);
805
806 // Xmu = 3.2 fm^-1 and as = 0.33 ==> sigma_{partonic} = 1.5mb
807 // Ntmax = 150
808 // v_2{2} = 0.102105 +/- 0.000266894
809 // v_2{4} = 0.0829477 +/- 0.00106158
810
811 genHi->SetNtMax(150); // NTMAX: number of timesteps (D=150)
812 genHi->SetXmu(3.2264); // parton screening mass in fm^(-1) (D=3.2264d0)
813
814 genHi->SetJetQuenching(0); // enable jet quenching
815 genHi->SetShadowing(1); // enable shadowing
816 genHi->SetDecaysOff(1); // neutral pion and heavy particle decays switched off
817 genHi->SetSpectators(0); // track spectators
818 //Boost into LHC lab frame
819 genHi->SetBoostLHC(1);
820 // genHi->Init();
821 genHi->SetRandomReactionPlane(kTRUE);
822
823 return genHi;
824}
825
826//__________________________________________________________________
827AliGenerator* AmptStringMelting()
828{
829 AliGenAmpt *genHi = new AliGenAmpt(-1);
830 //=========================================================================
831 // THE DECAYER
832 AliDecayer *decayer = new AliDecayerPythia();
833 cout << "*****************************************" << endl;
834 genHi->SetForceDecay( kHadronicD );
835 genHi->SetDecayer( decayer );
836 //=========================================================================
837 genHi->SetEnergyCMS(energy);
838 genHi->SetReferenceFrame("CMS");
839 genHi->SetProjectile("A",208,82);
840 genHi->SetTarget("A",208,82);
841
842 genHi->SetIsoft(4); //1: defaul - 4: string melting
843 genHi->SetStringFrag(0.5,0.9); //Lund string frgmentation parameters
844 genHi->SetPtHardMin (3);
845 //genHi->SetImpactParameterRange(9.,9.5);
846 genHi->SetImpactParameterRange(bMin,bMax);
847
848 // Xmu = 3.2 fm^-1 and as = 0.33 ==> sigma_{partonic} = 1.5mb
849 // Ntmax = 150
850 // v_2{2} = 0.102105 +/- 0.000266894
851 // v_2{4} = 0.0829477 +/- 0.00106158
852
853 genHi->SetNtMax(150); // NTMAX: number of timesteps (D=150)
854 genHi->SetXmu(3.2264); // parton screening mass in fm^(-1) (D=3.2264d0)
855
856 genHi->SetJetQuenching(0); // enable jet quenching
857 genHi->SetShadowing(1); // enable shadowing
858 genHi->SetDecaysOff(1); // neutral pion and heavy particle decays switched off
859 genHi->SetSpectators(0); // track spectators
860 //Boost into LHC lab frame
861 genHi->SetBoostLHC(1);
862// genHi->Init();
863 genHi->SetRandomReactionPlane(kTRUE);
864 return genHi;
865
866}
867
868//__________________________________________________________________
869AliGenerator* AmptStringMeltingNoART()
870{
871 AliGenAmpt *genHi = new AliGenAmpt(-1);
872 //=========================================================================
873 // THE DECAYER
874 AliDecayer *decayer = new AliDecayerPythia();
875 cout << "*****************************************" << endl;
876 genHi->SetForceDecay( kHadronicD );
877 genHi->SetDecayer( decayer );
878 //=========================================================================
879 genHi->SetEnergyCMS(energy);
880 genHi->SetReferenceFrame("CMS");
881 genHi->SetProjectile("A",208,82);
882 genHi->SetTarget("A",208,82);
883
884 genHi->SetIsoft(4); //1: defaul - 4: string melting
885 genHi->SetStringFrag(0.5,0.9); //Lund string frgmentation parameters
886 genHi->SetPtHardMin (3);
887 //genHi->SetImpactParameterRange(9.,9.5);
888 genHi->SetImpactParameterRange(bMin,bMax);
889
890 // Xmu = 3.2 fm^-1 and as = 0.33 ==> sigma_{partonic} = 1.5mb
891 // Ntmax = 150
892 // v_2{2} = 0.102105 +/- 0.000266894
893 // v_2{4} = 0.0829477 +/- 0.00106158
894
895 genHi->SetNtMax(3); // NTMAX: number of timesteps (D=150)
896 genHi->SetXmu(3.2264); // parton screening mass in fm^(-1) (D=3.2264d0)
897
898 genHi->SetJetQuenching(0); // enable jet quenching
899 genHi->SetShadowing(1); // enable shadowing
900 genHi->SetDecaysOff(1); // neutral pion and heavy particle decays switched off
901 genHi->SetSpectators(0); // track spectators
902 //Boost into LHC lab frame
903 genHi->SetBoostLHC(1);
904 // genHi->Init();
905 genHi->SetRandomReactionPlane(kTRUE);
906 return genHi;
907
908}
909
910
911//__________________________________________________________________
912AliGenerator* AmptReducedFlow()
913{
914 AliGenAmpt *genHi = new AliGenAmpt(-1);
915 //=========================================================================
916 // THE DECAYER
917 AliDecayer *decayer = new AliDecayerPythia();
918 cout << "*****************************************" << endl;
919 genHi->SetForceDecay( kHadronicD );
920 genHi->SetDecayer( decayer );
921 //=========================================================================
922 genHi->SetEnergyCMS(energy);
923 genHi->SetReferenceFrame("CMS");
924 genHi->SetProjectile("A",208,82);
925 genHi->SetTarget("A",208,82);
926
927 genHi->SetIsoft(4); //1: defaul - 4: string melting
928 genHi->SetStringFrag(0.5,0.9); //Lund string frgmentation parameters
929 genHi->SetPtHardMin (3);
930 //genHi->SetImpactParameterRange(9.,9.5);
931 genHi->SetImpactParameterRange(bMin,bMax);
932
933 // Xmu = 12.4 fm^-1 and as = 0.33 ==> sigma_{partonic} = 0.1mb
934 // Ntmax = 20
935 // flow estimates from Q-cumulants
936 // (POI, without weights)
937 // v_2{2} = 0.0549735 +/- 0.000270249
938 // v_2{4} = 0.0421905 +/- 0.00189449
939
940 genHi->SetNtMax(20); // NTMAX: number of timesteps (D=150)
941 genHi->SetXmu(12.4); // parton screening mass in fm^(-1) (D=3.2264d0)
942
943 genHi->SetJetQuenching(0); // enable jet quenching
944 genHi->SetShadowing(1); // enable shadowing
945 genHi->SetDecaysOff(1); // neutral pion and heavy particle decays switched off
946 genHi->SetSpectators(0); // track spectators
947 //Boost into LHC lab frame
948 genHi->SetBoostLHC(1);
949 // genHi->Init();
950 genHi->SetRandomReactionPlane(kTRUE);
951 return genHi;
952
953}
954
955//__________________________________________________________________
956AliGenerator* AmptpA()
957{
958 AliGenAmpt *genHi = new AliGenAmpt(-1);
959 //=========================================================================
960 // THE DECAYER
961 AliDecayer *decayer = new AliDecayerPythia();
962 cout << "*****************************************" << endl;
963 genHi->SetForceDecay( kHadronicD );
964 genHi->SetDecayer( decayer );
965 //=========================================================================
966 genHi->SetEnergyCMS(energy);
967 genHi->SetReferenceFrame("CMS");
968 genHi->SetProjectile("A", 208, 82);
969 genHi->SetTarget ("P", 1, 1);
970 genHi->SetIsoft(4); //1: defaul - 4: string melting
971 genHi->SetStringFrag(0.5,0.9); //Lund string frgmentation parameters
972 genHi->SetPtHardMin (3);
973 //genHi->SetImpactParameterRange(9.,9.5);
974 genHi->SetImpactParameterRange(bMin,bMax);
975 genHi->SetNtMax(1500); //NTMAX: number of timesteps (D=150)
976 genHi->SetXmu(3.2264); //parton screening mass in fm^(-1) (D=3.2264d0)
977 genHi->SetJetQuenching(0); // enable jet quenching
978 genHi->SetShadowing(1); // enable shadowing
979 genHi->SetDecaysOff(1); // neutral pion and heavy particle decays switched off
980 genHi->SetSpectators(0); // track spectators
981 //Boost into LHC lab frame
982 genHi->SetBoostLHC(1);
983 // genHi->Init();
984 genHi->SetRandomReactionPlane(kTRUE);
985 return genHi;
986
987}