]> git.uio.no Git - u/mrichter/AliRoot.git/blame - FMD/scripts/pdc06_config.C
Fixed up this script
[u/mrichter/AliRoot.git] / FMD / scripts / pdc06_config.C
CommitLineData
cf806c6a 1//
2// Configuration for the Physics Data Challenge 2006
3//
4// Physics Run composition (Presented by A.Dainese, 21/03/06
5//
6// Assuming inel = 70 mb (PPR v.1, p.64)
7//
8// 84.98% of MSEL=0 events (including diffractive) with
9// QQbar switched off (these events will include injected J/psi) => kPyMbNoHvq
10//
11// 14.14% of MSEL=1 events with ccbar (in 4 subsamples) => kCharmpp14000wmi
12// bin 1 25% (3.535%): 2.76 < pthard < 3 GeV/c
13// bin 2 40% (5.656%): 3 < pthard < 4 GeV/c
14// bin 3 29% (4.101%): 4 < pthard < 8 GeV/c
15// bin 4 6% (0.848%): pthard > 8 GeV/c
16//
17// 0.73% of MSEL=1 events with bbbar (in 4 subsamples) => kBeautypp14000wmi
18// bin 1 5% (0.037%): 2.76 < pthard < 4 GeV/c
19// bin 2 31% (0.226%): 4 < pthard < 6 GeV/c
20// bin 3 28% (0.204%): 6 < pthard < 8 GeV/c
21// bin 4 36% (0.263%): pthard >8 GeV/c
22//
23// 0.075% of MSEL=0 events with QQbar switched off and 1 Omega- => kPyOmegaMinus
24// 0.075% of MSEL=0 events with QQbar switched off and 1 OmegaBar+ => kPyOmegaPlus
25
26// One can use the configuration macro in compiled mode by
27// root [0] gSystem->Load("libgeant321");
28// root [0] gSystem->SetIncludePath("-I$ROOTSYS/include -I$ALICE_ROOT/include\
29// -I$ALICE_ROOT -I$ALICE/geant3/TGeant3");
30// root [0] .x grun.C(1,"Config_PDC06.C++")
31
32#if !defined(__CINT__) || defined(__MAKECINT__)
33#include <Riostream.h>
34#include <TRandom.h>
35#include <TDatime.h>
36#include <TSystem.h>
37#include <TVirtualMC.h>
38#include <TGeant3TGeo.h>
39#include "EVGEN/AliGenCocktail.h"
40#include "EVGEN/AliGenParam.h"
41#include "EVGEN/AliGenMUONlib.h"
42#include "STEER/AliRunLoader.h"
43#include "STEER/AliRun.h"
44#include "STEER/AliConfig.h"
45#include "PYTHIA6/AliDecayerPythia.h"
46#include "PYTHIA6/AliGenPythia.h"
f7a1cc68 47#include "STEER/AliMagF.h"
cf806c6a 48#include "STRUCT/AliBODY.h"
49#include "STRUCT/AliMAG.h"
50#include "STRUCT/AliABSOv0.h"
51#include "STRUCT/AliDIPOv2.h"
52#include "STRUCT/AliHALL.h"
53#include "STRUCT/AliFRAMEv2.h"
54#include "STRUCT/AliSHILv2.h"
55#include "STRUCT/AliPIPEv0.h"
56#include "ITS/AliITSgeom.h"
57#include "ITS/AliITSvPPRasymmFMD.h"
58#include "TPC/AliTPCv2.h"
59#include "TOF/AliTOFv5T0.h"
60#include "RICH/AliRICHv1.h"
61#include "ZDC/AliZDCv2.h"
62#include "TRD/AliTRDv1.h"
63#include "FMD/AliFMDv1.h"
64#include "MUON/AliMUONv1.h"
65#include "PHOS/AliPHOSv1.h"
66#include "PMD/AliPMDv1.h"
67#include "START/AliSTARTv1.h"
68#include "EMCAL/AliEMCALv2.h"
69#include "CRT/AliCRTv0.h"
70#include "VZERO/AliVZEROv7.h"
71#endif
72
73
74enum PDC06Proc_t
75{
76//--- Heavy Flavour Production ---
77 kCharmPbPb5500, kCharmpPb8800, kCharmpp14000, kCharmpp14000wmi,
78 kD0PbPb5500, kD0pPb8800, kD0pp14000,
79 kDPlusPbPb5500, kDPluspPb8800, kDPluspp14000,
80 kBeautyPbPb5500, kBeautypPb8800, kBeautypp14000, kBeautypp14000wmi,
81// -- Pythia Mb
82 kPyMbNoHvq, kPyOmegaPlus, kPyOmegaMinus, kRunMax
83};
84
85const char * pprRunName[] = {
86 "kCharmPbPb5500", "kCharmpPb8800", "kCharmpp14000", "kCharmpp14000wmi",
87 "kD0PbPb5500", "kD0pPb8800", "kD0pp14000",
88 "kDPlusPbPb5500", "kDPluspPb8800", "kDPluspp14000",
89 "kBeautyPbPb5500", "kBeautypPb8800", "kBeautypp14000", "kBeautypp14000wmi",
90 "kPyMbNoHvq", "kPyOmegaPlus", "kPyOmegaMinus"
91};
92
93
94//--- Decay Mode ---
95enum DecayHvFl_t
96{
97 kNature, kHadr, kSemiEl, kSemiMu
98};
99//--- Rapidity Cut ---
100enum YCut_t
101{
102 kFull, kBarrel, kMuonArm
103};
cf806c6a 104//--- Trigger config ---
105enum TrigConf_t
106{
107 kDefaultPPTrig, kDefaultPbPbTrig
108};
109
110const char * TrigConfName[] = {
111 "p-p","Pb-Pb"
112};
113
114//--- Functions ---
115AliGenPythia *PythiaHVQ(PDC06Proc_t proc);
116AliGenerator *MbCocktail();
117AliGenerator *PyMbTriggered(Int_t pdg);
118void ProcessEnvironmentVars();
119
120// This part for configuration
121static PDC06Proc_t proc = kPyOmegaPlus;
122static DecayHvFl_t decHvFl = kNature;
123static YCut_t ycut = kFull;
f7a1cc68 124static AliMagF::BMap_t mag = AliMagF::k5kG;
cf806c6a 125static TrigConf_t trig = kDefaultPPTrig; // default pp trigger configuration
126static Int_t runNumber= 0;
127//========================//
128// Set Random Number seed //
129//========================//
130TDatime dt;
131static UInt_t seed = dt.Get();
132
133// nEvts = -1 : you get 1 QQbar pair and all the fragmentation and
134// decay chain
135// nEvts = N>0 : you get N charm / beauty Hadrons
136Int_t nEvts = -1;
137// stars = kTRUE : all heavy resonances and their decay stored
138// = kFALSE: only final heavy hadrons and their decays stored
139Bool_t stars = kTRUE;
140
141// To be used only with kCharmpp14000wmi and kBeautypp14000wmi
142// To get a "reasonable" agreement with MNR results, events have to be
143// generated with the minimum ptHard set to 2.76 GeV.
144// To get a "perfect" agreement with MNR results, events have to be
145// generated in four ptHard bins with the following relative
146// normalizations:
147// CHARM
148// 2.76-3 GeV: 25%
149// 3-4 GeV: 40%
150// 4-8 GeV: 29%
151// >8 GeV: 6%
152// BEAUTY
153// 2.76-4 GeV: 5%
154// 4-6 GeV: 31%
155// 6-8 GeV: 28%
156// >8 GeV: 36%
157Float_t ptHardMin = 2.76;
158Float_t ptHardMax = -1.;
159
160
161// Comment line
162static TString comment;
163
164void Config()
165{
166
167
168 // Get settings from environment variables
169 ProcessEnvironmentVars();
170
171 // libraries required by geant321
172#if defined(__CINT__)
173 gSystem->Load("libgeant321");
174#endif
175
176 new TGeant3TGeo("C++ Interface to Geant3");
177
178 // Output every 100 tracks
179 ((TGeant3*)gMC)->SetSWIT(4,100);
180
181 //=======================================================================
182
183 // Run loader
184 AliRunLoader* rl=0x0;
185 rl = AliRunLoader::Open("galice.root",
186 AliConfig::GetDefaultEventFolderName(),
187 "recreate");
188 if (rl == 0x0)
189 {
190 gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
191 return;
192 }
193 rl->SetCompressionLevel(2);
194 rl->SetNumberOfEventsPerFile(1000);
195 gAlice->SetRunLoader(rl);
196
197 // Run number
198 gAlice->SetRunNumber(runNumber);
199
200 // Set the trigger configuration
28da60d3 201 AliSimulation::Instance()->SetTriggerConfig(TrigConfName[trig]);
cf806c6a 202 cout<<"Trigger configuration is set to "<<TrigConfName[trig]<<endl;
203
204 //
205 //=======================================================================
206 // ************* STEERING parameters FOR ALICE SIMULATION **************
207 // --- Specify event type to be tracked through the ALICE setup
208 // --- All positions are in cm, angles in degrees, and P and E in GeV
209
210
211 gMC->SetProcess("DCAY",1);
212 gMC->SetProcess("PAIR",1);
213 gMC->SetProcess("COMP",1);
214 gMC->SetProcess("PHOT",1);
215 gMC->SetProcess("PFIS",0);
216 gMC->SetProcess("DRAY",0);
217 gMC->SetProcess("ANNI",1);
218 gMC->SetProcess("BREM",1);
219 gMC->SetProcess("MUNU",1);
220 gMC->SetProcess("CKOV",1);
221 gMC->SetProcess("HADR",1);
222 gMC->SetProcess("LOSS",2);
223 gMC->SetProcess("MULS",1);
224 gMC->SetProcess("RAYL",1);
225
226 Float_t cut = 1.e-3; // 1MeV cut by default
227 Float_t tofmax = 1.e10;
228
229 gMC->SetCut("CUTGAM", cut);
230 gMC->SetCut("CUTELE", cut);
231 gMC->SetCut("CUTNEU", cut);
232 gMC->SetCut("CUTHAD", cut);
233 gMC->SetCut("CUTMUO", cut);
234 gMC->SetCut("BCUTE", cut);
235 gMC->SetCut("BCUTM", cut);
236 gMC->SetCut("DCUTE", cut);
237 gMC->SetCut("DCUTM", cut);
238 gMC->SetCut("PPCUTM", cut);
239 gMC->SetCut("TOFMAX", tofmax);
240
241 //======================//
242 // Set External decayer //
243 //======================//
244 TVirtualMCDecayer* decayer = new AliDecayerPythia();
245 // DECAYS
246 //
247 switch(decHvFl) {
248 case kNature:
249 decayer->SetForceDecay(kAll);
250 break;
251 case kHadr:
252 decayer->SetForceDecay(kHadronicD);
253 break;
254 case kSemiEl:
255 decayer->SetForceDecay(kSemiElectronic);
256 break;
257 case kSemiMu:
258 decayer->SetForceDecay(kSemiMuonic);
259 break;
260 }
261 decayer->Init();
262 gMC->SetExternalDecayer(decayer);
263
264 //=========================//
265 // Generator Configuration //
266 //=========================//
267 AliGenerator* gener = 0x0;
268
269 if (proc <= kBeautypp14000wmi) {
270 AliGenPythia *pythia = PythiaHVQ(proc);
271 // FeedDown option
272 pythia->SetFeedDownHigherFamily(kFALSE);
273 // Stack filling option
274 if(!stars) pythia->SetStackFillOpt(AliGenPythia::kParentSelection);
275 // Set Count mode
276 if(nEvts>0) pythia->SetCountMode(AliGenPythia::kCountParents);
277 //
278 // DECAYS
279 //
280 switch(decHvFl) {
281 case kNature:
282 pythia->SetForceDecay(kAll);
283 break;
284 case kHadr:
285 pythia->SetForceDecay(kHadronicD);
286 break;
287 case kSemiEl:
288 pythia->SetForceDecay(kSemiElectronic);
289 break;
290 case kSemiMu:
291 pythia->SetForceDecay(kSemiMuonic);
292 break;
293 }
294 //
295 // GEOM & KINE CUTS
296 //
297 pythia->SetMomentumRange(0,99999999);
298 pythia->SetPhiRange(0., 360.);
299 pythia->SetThetaRange(0,180);
300 switch(ycut) {
301 case kFull:
302 pythia->SetYRange(-999,999);
303 break;
304 case kBarrel:
305 pythia->SetYRange(-2,2);
306 break;
307 case kMuonArm:
308 pythia->SetYRange(1,6);
309 break;
310 }
311 gener = pythia;
312 } else if (proc == kPyMbNoHvq) {
313 gener = MbCocktail();
314 } else if (proc == kPyOmegaMinus) {
315 gener = PyMbTriggered(3334);
316 } else if (proc == kPyOmegaPlus) {
317 gener = PyMbTriggered(-3334);
318 }
319
320
321
322 // PRIMARY VERTEX
323
324 gener->SetOrigin(0., 0., 0.); // vertex position
325
326 // Size of the interaction diamond
327 // Longitudinal
328 Float_t sigmaz = 7.55 / TMath::Sqrt(2.); // [cm]
329
330 // Transverse
331 Float_t betast = 10; // beta* [m]
332 Float_t eps = 3.75e-6; // emittance [m]
333 Float_t gamma = 7000. / 0.938272; // relativistic gamma [1]
334 Float_t sigmaxy = TMath::Sqrt(eps * betast / gamma) / TMath::Sqrt(2.) * 100.; // [cm]
335 printf("\n \n Diamond size x-y: %10.3e z: %10.3e\n \n", sigmaxy, sigmaz);
336
337 gener->SetSigma(sigmaxy, sigmaxy, sigmaz); // Sigma in (X,Y,Z) (cm) on IP position
338 gener->SetCutVertexZ(3.); // Truncate at 3 sigma
339 gener->SetVertexSmear(kPerEvent);
340
341 gener->Init();
342
343 // FIELD
344
345 if (mag == k2kG) {
346 comment = comment.Append(" | L3 field 0.2 T");
cf806c6a 347 } else if (mag == k5kG) {
348 comment = comment.Append(" | L3 field 0.5 T");
349 }
350 printf("\n \n Comment: %s \n \n", comment.Data());
351
4642ac4b 352 AliMagF* field = new AliMagF("Maps","Maps",1.,1.,mag);
cf806c6a 353 field->SetL3ConstField(0); //Using const. field in the barrel
f7a1cc68 354 TGeoGlobalMagField::Instance()->SetField(field);
cf806c6a 355
f7a1cc68 356 rl->CdGAFile();
cf806c6a 357
358 Int_t iABSO = 1;
359 Int_t iCRT = 0;
360 Int_t iDIPO = 1;
361 Int_t iEMCAL = 1;
362 Int_t iFMD = 1;
363 Int_t iFRAME = 1;
364 Int_t iHALL = 1;
365 Int_t iITS = 1;
366 Int_t iMAG = 1;
367 Int_t iMUON = 1;
368 Int_t iPHOS = 1;
369 Int_t iPIPE = 1;
370 Int_t iPMD = 1;
371 Int_t iRICH = 1;
372 Int_t iSHIL = 1;
373 Int_t iSTART = 1;
374 Int_t iTOF = 1;
375 Int_t iTPC = 1;
376 Int_t iTRD = 1;
377 Int_t iVZERO = 1;
378 Int_t iZDC = 1;
379
380
381 //=================== Alice BODY parameters =============================
382 AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
383
384
385 if (iMAG)
386 {
387 //=================== MAG parameters ============================
388 // --- Start with Magnet since detector layouts may be depending ---
389 // --- on the selected Magnet dimensions ---
390 AliMAG *MAG = new AliMAG("MAG", "Magnet");
391 }
392
393
394 if (iABSO)
395 {
396 //=================== ABSO parameters ============================
397 AliABSO *ABSO = new AliABSOv0("ABSO", "Muon Absorber");
398 }
399
400 if (iDIPO)
401 {
402 //=================== DIPO parameters ============================
403
404 AliDIPO *DIPO = new AliDIPOv2("DIPO", "Dipole version 2");
405 }
406
407 if (iHALL)
408 {
409 //=================== HALL parameters ============================
410
411 AliHALL *HALL = new AliHALL("HALL", "Alice Hall");
412 }
413
414
415 if (iFRAME)
416 {
417 //=================== FRAME parameters ============================
418
419 AliFRAMEv2 *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
420 }
421
422 if (iSHIL)
423 {
424 //=================== SHIL parameters ============================
425
426 AliSHIL *SHIL = new AliSHILv2("SHIL", "Shielding Version 2");
427 }
428
429
430 if (iPIPE)
431 {
432 //=================== PIPE parameters ============================
433
434 AliPIPE *PIPE = new AliPIPEv0("PIPE", "Beam Pipe");
435 }
436
437 if(iITS) {
438
439 //=================== ITS parameters ============================
440 //
441 // As the innermost detector in ALICE, the Inner Tracking System "impacts" on
442 // almost all other detectors. This involves the fact that the ITS geometry
443 // still has several options to be followed in parallel in order to determine
444 // the best set-up which minimizes the induced background. All the geometries
445 // available to date are described in the following. Read carefully the comments
446 // and use the default version (the only one uncommented) unless you are making
447 // comparisons and you know what you are doing. In this case just uncomment the
448 // ITS geometry you want to use and run Aliroot.
449 //
450 // Detailed geometries:
451 //
452 //
453 //AliITS *ITS = new AliITSv5symm("ITS","Updated ITS TDR detailed version with symmetric services");
454 //
455 //AliITS *ITS = new AliITSv5asymm("ITS","Updates ITS TDR detailed version with asymmetric services");
456 //
457 AliITSvPPRasymmFMD *ITS = new AliITSvPPRasymmFMD("ITS","New ITS PPR detailed version with asymmetric services");
458 ITS->SetMinorVersion(2); // don't touch this parameter if you're not an ITS developer
459 ITS->SetReadDet(kFALSE); // don't touch this parameter if you're not an ITS developer
460 // ITS->SetWriteDet("$ALICE_ROOT/ITS/ITSgeometry_vPPRasymm2.det"); // don't touch this parameter if you're not an ITS developer
461 ITS->SetThicknessDet1(200.); // detector thickness on layer 1 must be in the range [100,300]
462 ITS->SetThicknessDet2(200.); // detector thickness on layer 2 must be in the range [100,300]
463 ITS->SetThicknessChip1(150.); // chip thickness on layer 1 must be in the range [150,300]
464 ITS->SetThicknessChip2(150.); // chip thickness on layer 2 must be in the range [150,300]
465 ITS->SetRails(0); // 1 --> rails in ; 0 --> rails out
466 ITS->SetCoolingFluid(1); // 1 --> water ; 0 --> freon
467
468 // Coarse geometries (warning: no hits are produced with these coarse geometries and they unuseful
469 // for reconstruction !):
470 //
471 //
472 //AliITSvPPRcoarseasymm *ITS = new AliITSvPPRcoarseasymm("ITS","New ITS PPR coarse version with asymmetric services");
473 //ITS->SetRails(0); // 1 --> rails in ; 0 --> rails out
474 //ITS->SetSupportMaterial(0); // 0 --> Copper ; 1 --> Aluminum ; 2 --> Carbon
475 //
476 //AliITS *ITS = new AliITSvPPRcoarsesymm("ITS","New ITS PPR coarse version with symmetric services");
477 //ITS->SetRails(0); // 1 --> rails in ; 0 --> rails out
478 //ITS->SetSupportMaterial(0); // 0 --> Copper ; 1 --> Aluminum ; 2 --> Carbon
479 //
480 //
481 //
482 // Geant3 <-> EUCLID conversion
483 // ============================
484 //
485 // SetEUCLID is a flag to output (=1) or not to output (=0) both geometry and
486 // media to two ASCII files (called by default ITSgeometry.euc and
487 // ITSgeometry.tme) in a format understandable to the CAD system EUCLID.
488 // The default (=0) means that you dont want to use this facility.
489 //
490 ITS->SetEUCLID(0);
491 }
492
493 if (iTPC)
494 {
495 //============================ TPC parameters =====================
496 AliTPC *TPC = new AliTPCv2("TPC", "Default");
497 }
498
499
500 if (iTOF) {
501 //=================== TOF parameters ============================
502 AliTOF *TOF = new AliTOFv5T0("TOF", "normal TOF");
503 // Partial geometry: modules at 2,3,4,6,7,11,12,14,15,16,17
504 // starting at 6h in positive direction
505 Int_t TOFSectors[18]={-1,-1,0,0,0,-1,0,0,-1,-1,-1,0,0,-1,0,0,0,0};
506 TOF->SetTOFSectors(TOFSectors);
507 }
508
509
510 if (iRICH)
511 {
512 //=================== RICH parameters ===========================
513 AliRICH *RICH = new AliRICHv1("RICH", "normal RICH");
514
515 }
516
517
518 if (iZDC)
519 {
520 //=================== ZDC parameters ============================
521
522 AliZDC *ZDC = new AliZDCv2("ZDC", "normal ZDC");
523 }
524
525 if (iTRD)
526 {
527 //=================== TRD parameters ============================
528
529 AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
530 AliTRDgeometry *geoTRD = TRD->GetGeometry();
531 // Partial geometry: modules at 2,3,4,6,11,12,14,15
532 // starting at 6h in positive direction
533 geoTRD->SetSMstatus(0,0);
534 geoTRD->SetSMstatus(1,0);
535 geoTRD->SetSMstatus(5,0);
536 geoTRD->SetSMstatus(7,0);
537 geoTRD->SetSMstatus(8,0);
538 geoTRD->SetSMstatus(9,0);
539 geoTRD->SetSMstatus(10,0);
540 geoTRD->SetSMstatus(13,0);
541 geoTRD->SetSMstatus(16,0);
542 geoTRD->SetSMstatus(17,0);
543 }
544
545 if (iFMD)
546 {
547 //=================== FMD parameters ============================
548 AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
549 }
550
551 if (iMUON)
552 {
553 //=================== MUON parameters ===========================
554 // New MUONv1 version (geometry defined via builders)
555 AliMUON *MUON = new AliMUONv1("MUON", "default");
556 }
557 //=================== PHOS parameters ===========================
558
559 if (iPHOS)
560 {
561 AliPHOS *PHOS = new AliPHOSv1("PHOS", "IHEP");
562 }
563
564
565 if (iPMD)
566 {
567 //=================== PMD parameters ============================
568 AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
569 }
570
571 if (iSTART)
572 {
573 //=================== START parameters ============================
574 AliSTART *START = new AliSTARTv1("START", "START Detector");
575 }
576
577 if (iEMCAL)
578 {
579 //=================== EMCAL parameters ============================
580 AliEMCAL *EMCAL = new AliEMCALv2("EMCAL", "SHISH_77_TRD1_2X2_FINAL_110DEG");
581 }
582
583 if (iCRT)
584 {
585 //=================== CRT parameters ============================
586 AliCRT *CRT = new AliCRTv0("CRT", "normal ACORDE");
587 }
588
589 if (iVZERO)
590 {
591 //=================== CRT parameters ============================
592 AliVZERO *VZERO = new AliVZEROv7("VZERO", "normal VZERO");
593 }
594}
595
596// PYTHIA
597
598AliGenPythia *PythiaHVQ(PDC06Proc_t proc) {
599//*******************************************************************//
600// Configuration file for charm / beauty generation with PYTHIA //
601// //
602// The parameters have been tuned in order to reproduce the inclusive//
603// heavy quark pt distribution given by the NLO pQCD calculation by //
604// Mangano, Nason and Ridolfi. //
605// //
606// For details and for the NORMALIZATION of the yields see: //
607// N.Carrer and A.Dainese, //
608// "Charm and beauty production at the LHC", //
609// ALICE-INT-2003-019, [arXiv:hep-ph/0311225]; //
610// PPR Chapter 6.6, CERN/LHCC 2005-030 (2005). //
611//*******************************************************************//
612 AliGenPythia * gener = 0x0;
613
614 switch(proc) {
615 case kCharmPbPb5500:
616 comment = comment.Append(" Charm in Pb-Pb at 5.5 TeV");
617 gener = new AliGenPythia(nEvts);
618 gener->SetProcess(kPyCharmPbPbMNR);
619 gener->SetStrucFunc(kCTEQ4L);
620 gener->SetPtHard(2.1,-1.0);
621 gener->SetEnergyCMS(5500.);
622 gener->SetNuclei(208,208);
623 break;
624 case kCharmpPb8800:
625 comment = comment.Append(" Charm in p-Pb at 8.8 TeV");
626 gener = new AliGenPythia(nEvts);
627 gener->SetProcess(kPyCharmpPbMNR);
628 gener->SetStrucFunc(kCTEQ4L);
629 gener->SetPtHard(2.1,-1.0);
630 gener->SetEnergyCMS(8800.);
631 gener->SetProjectile("P",1,1);
632 gener->SetTarget("Pb",208,82);
633 break;
634 case kCharmpp14000:
635 comment = comment.Append(" Charm in pp at 14 TeV");
636 gener = new AliGenPythia(nEvts);
637 gener->SetProcess(kPyCharmppMNR);
638 gener->SetStrucFunc(kCTEQ4L);
639 gener->SetPtHard(2.1,-1.0);
640 gener->SetEnergyCMS(14000.);
641 break;
642 case kCharmpp14000wmi:
643 comment = comment.Append(" Charm in pp at 14 TeV with mult. interactions");
644 gener = new AliGenPythia(-1);
645 gener->SetProcess(kPyCharmppMNRwmi);
646 gener->SetStrucFunc(kCTEQ5L);
647 gener->SetPtHard(ptHardMin,ptHardMax);
648 gener->SetEnergyCMS(14000.);
649 break;
650 case kD0PbPb5500:
651 comment = comment.Append(" D0 in Pb-Pb at 5.5 TeV");
652 gener = new AliGenPythia(nEvts);
653 gener->SetProcess(kPyD0PbPbMNR);
654 gener->SetStrucFunc(kCTEQ4L);
655 gener->SetPtHard(2.1,-1.0);
656 gener->SetEnergyCMS(5500.);
657 gener->SetNuclei(208,208);
658 break;
659 case kD0pPb8800:
660 comment = comment.Append(" D0 in p-Pb at 8.8 TeV");
661 gener = new AliGenPythia(nEvts);
662 gener->SetProcess(kPyD0pPbMNR);
663 gener->SetStrucFunc(kCTEQ4L);
664 gener->SetPtHard(2.1,-1.0);
665 gener->SetEnergyCMS(8800.);
666 gener->SetProjectile("P",1,1);
667 gener->SetTarget("Pb",208,82);
668 break;
669 case kD0pp14000:
670 comment = comment.Append(" D0 in pp at 14 TeV");
671 gener = new AliGenPythia(nEvts);
672 gener->SetProcess(kPyD0ppMNR);
673 gener->SetStrucFunc(kCTEQ4L);
674 gener->SetPtHard(2.1,-1.0);
675 gener->SetEnergyCMS(14000.);
676 break;
677 case kDPlusPbPb5500:
678 comment = comment.Append(" DPlus in Pb-Pb at 5.5 TeV");
679 gener = new AliGenPythia(nEvts);
680 gener->SetProcess(kPyDPlusPbPbMNR);
681 gener->SetStrucFunc(kCTEQ4L);
682 gener->SetPtHard(2.1,-1.0);
683 gener->SetEnergyCMS(5500.);
684 gener->SetNuclei(208,208);
685 break;
686 case kDPluspPb8800:
687 comment = comment.Append(" DPlus in p-Pb at 8.8 TeV");
688 gener = new AliGenPythia(nEvts);
689 gener->SetProcess(kPyDPluspPbMNR);
690 gener->SetStrucFunc(kCTEQ4L);
691 gener->SetPtHard(2.1,-1.0);
692 gener->SetEnergyCMS(8800.);
693 gener->SetProjectile("P",1,1);
694 gener->SetTarget("Pb",208,82);
695 break;
696 case kDPluspp14000:
697 comment = comment.Append(" DPlus in pp at 14 TeV");
698 gener = new AliGenPythia(nEvts);
699 gener->SetProcess(kPyDPlusppMNR);
700 gener->SetStrucFunc(kCTEQ4L);
701 gener->SetPtHard(2.1,-1.0);
702 gener->SetEnergyCMS(14000.);
703 break;
704 case kBeautyPbPb5500:
705 comment = comment.Append(" Beauty in Pb-Pb at 5.5 TeV");
706 gener = new AliGenPythia(nEvts);
707 gener->SetProcess(kPyBeautyPbPbMNR);
708 gener->SetStrucFunc(kCTEQ4L);
709 gener->SetPtHard(2.75,-1.0);
710 gener->SetEnergyCMS(5500.);
711 gener->SetNuclei(208,208);
712 break;
713 case kBeautypPb8800:
714 comment = comment.Append(" Beauty in p-Pb at 8.8 TeV");
715 gener = new AliGenPythia(nEvts);
716 gener->SetProcess(kPyBeautypPbMNR);
717 gener->SetStrucFunc(kCTEQ4L);
718 gener->SetPtHard(2.75,-1.0);
719 gener->SetEnergyCMS(8800.);
720 gener->SetProjectile("P",1,1);
721 gener->SetTarget("Pb",208,82);
722 break;
723 case kBeautypp14000:
724 comment = comment.Append(" Beauty in pp at 14 TeV");
725 gener = new AliGenPythia(nEvts);
726 gener->SetProcess(kPyBeautyppMNR);
727 gener->SetStrucFunc(kCTEQ4L);
728 gener->SetPtHard(2.75,-1.0);
729 gener->SetEnergyCMS(14000.);
730 break;
731 case kBeautypp14000wmi:
732 comment = comment.Append(" Beauty in pp at 14 TeV with mult. interactions");
733 gener = new AliGenPythia(-1);
734 gener->SetProcess(kPyBeautyppMNRwmi);
735 gener->SetStrucFunc(kCTEQ5L);
736 gener->SetPtHard(ptHardMin,ptHardMax);
737 gener->SetEnergyCMS(14000.);
738 break;
739 }
740
741 return gener;
742}
743
744AliGenerator* MbCocktail()
745{
746 comment = comment.Append(" pp at 14 TeV: Pythia low-pt, no heavy quarks + J/Psi from parameterisation");
747 AliGenCocktail * gener = new AliGenCocktail();
748 gener->UsePerEventRates();
749
750// Pythia
751
752 AliGenPythia* pythia = new AliGenPythia(-1);
753 pythia->SetMomentumRange(0, 999999.);
754 pythia->SetThetaRange(0., 180.);
755 pythia->SetYRange(-12.,12.);
756 pythia->SetPtRange(0,1000.);
757 pythia->SetProcess(kPyMb);
758 pythia->SetEnergyCMS(14000.);
759 pythia->SwitchHFOff();
760
761// J/Psi parameterisation
762
763 AliGenParam* jpsi = new AliGenParam(1, AliGenMUONlib::kJpsi, "CDF scaled", "Jpsi");
764 jpsi->SetPtRange(0.,100.);
765 jpsi->SetYRange(-8., 8.);
766 jpsi->SetPhiRange(0., 360.);
767 jpsi->SetForceDecay(kAll);
768
769 gener->AddGenerator(pythia, "Pythia", 1.);
770 gener->AddGenerator(jpsi, "J/Psi", 8.e-4);
771
772 return gener;
773}
774
775AliGenerator* PyMbTriggered(Int_t pdg)
776{
777 AliGenPythia* pythia = new AliGenPythia(-1);
778 pythia->SetMomentumRange(0, 999999.);
779 pythia->SetThetaRange(0., 180.);
780 pythia->SetYRange(-12.,12.);
781 pythia->SetPtRange(0,1000.);
782 pythia->SetProcess(kPyMb);
783 pythia->SetEnergyCMS(14000.);
784 pythia->SetTriggerParticle(pdg, 0.9);
785 return pythia;
786}
787
788void ProcessEnvironmentVars()
789{
790 cout << "Processing environment variables" << endl;
791 // Random Number seed
792 if (gSystem->Getenv("CONFIG_SEED")) {
793 seed = atoi(gSystem->Getenv("CONFIG_SEED"));
794 }
795
796 gRandom->SetSeed(seed);
797 cout<<"Seed for random number generation= "<<seed<<endl;
798
799 // Run Number
800 if (gSystem->Getenv("DC_RUN")) {
801 runNumber = atoi(gSystem->Getenv("DC_RUN"));
802 }
803 cout<<"Run number "<<runNumber<<endl;
804
805 // Run type
806 if (gSystem->Getenv("DC_RUN_TYPE")) {
807 for (Int_t iRun = 0; iRun < kRunMax; iRun++) {
808 if (strcmp(gSystem->Getenv("DC_RUN_TYPE"), pprRunName[iRun])==0) {
809 proc = (PDC06Proc_t)iRun;
810 cout<<"Run type set to "<<pprRunName[iRun]<<endl;
811 }
812 }
813 } else {
814 // Define the run type randomly
815
816 // The array below contains the cumulative probability
817 // for the following cases:
818 // kPyMbNoHvq, kCharmpp14000wmi,kBeautypp14000wmi,kPyOmegaMinus,kPyOmegaPlus
819 Double_t probType[] = {0.0,0.8498,0.9912,0.9985,0.99925,1.0};
820 Int_t iType = TMath::BinarySearch(6,probType,gRandom->Rndm());
821
822 switch (iType) {
823 case 0:
824 proc = kPyMbNoHvq;
825 break;
826 case 1:
827 proc = kCharmpp14000wmi;
828 {
829 // Define ptHardMin,ptHardMax
830 // The array below contains the cumulative probability
831 // for the different pthard bins
832 Double_t probPtCharm[] = {0.0,0.25,0.65,0.94,1.0};
833 Int_t iPt = TMath::BinarySearch(5,probPtCharm,gRandom->Rndm());
834 switch (iPt) {
835 case 0:
836 ptHardMin = 2.76;
837 ptHardMax = 3.0;
838 break;
839 case 1:
840 ptHardMin = 3.0;
841 ptHardMax = 4.0;
842 break;
843 case 2:
844 ptHardMin = 4.0;
845 ptHardMax = 8.0;
846 break;
847 case 3:
848 ptHardMin = 8.0;
849 ptHardMax = -1.0;
850 break;
851 default:
852 cout << "ProcessEnvironmentVars: Wrong pthard bin" << endl;
853 }
854 }
855 break;
856 case 2:
857 proc = kBeautypp14000wmi;
858 {
859 // Define ptHardMin,ptHardMax
860 // The array below contains the cumulative probability
861 // for the different pthard bins
862 Double_t probPtBeauty[] = {0.0,0.05,0.36,0.64,1.0};
863 Int_t iPt = TMath::BinarySearch(5,probPtBeauty,gRandom->Rndm());
864 switch (iPt) {
865 case 0:
866 ptHardMin = 2.76;
867 ptHardMax = 4.0;
868 break;
869 case 1:
870 ptHardMin = 4.0;
871 ptHardMax = 6.0;
872 break;
873 case 2:
874 ptHardMin = 6.0;
875 ptHardMax = 8.0;
876 break;
877 case 3:
878 ptHardMin = 8.0;
879 ptHardMax = -1.0;
880 break;
881 default:
882 cout << "ProcessEnvironmentVars: Wrong pthard bin" << endl;
883 }
884 }
885 break;
886 case 3:
887 proc = kPyOmegaMinus;
888 break;
889 case 4:
890 proc = kPyOmegaPlus;
891 break;
892 default:
893 cout << "ProcessEnvironmentVars: Wrong run type" << endl;
894 }
895 cout<<"Run type set to "<<pprRunName[proc]<<endl;
896 cout<<"ptHard limits: "<<ptHardMin<<" to " <<ptHardMax<<endl;
897 }
898
899}
900
901
902