Added RuleChecker jars in trunk and updated CMake .pkg files
[u/mrichter/AliRoot.git] / PYTHIA8 / AliPythia8.cxx
1
2 /**************************************************************************
3  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4  *                                                                        *
5  * Author: The ALICE Off-line Project.                                    *
6  * Contributors are mentioned in the code where appropriate.              *
7  *                                                                        *
8  * Permission to use, copy, modify and distribute this software and its   *
9  * documentation strictly for non-commercial purposes is hereby granted   *
10  * without fee, provided that the above copyright notice appears in all   *
11  * copies and that both the copyright notice and this permission notice   *
12  * appear in the supporting documentation. The authors make no claims     *
13  * about the suitability of this software for any purpose. It is          *
14  * provided "as is" without express or implied warranty.                  *
15  **************************************************************************/
16
17 /* $Id$ */
18 #include <TString.h>
19 #include <TVector3.h>
20 #include <TMath.h>
21
22 #include "AliPythia8.h"
23 #include "AliLog.h"
24 #include "AliStack.h"
25 #include "AliPythiaRndm.h"
26
27
28 ClassImp(AliPythia8)
29     
30 // Particles produced in string fragmentation point directly to either of the two endpoints
31 // of the string (depending in the side they were generated from).
32 //    SetMSTU(16,2); // ????
33 //  String drawing almost completely minimizes string length.
34 //  Probability that an additional interaction gives two gluons
35 //  ... with color connection to nearest neighbours
36 //    SetPARP(85,0.9);
37 //  ... as closed gluon loop
38 //    SetPARP(86,0.95);
39 // Lambda_FSR scale.
40 //      SetPARJ(81, 0.29);
41 // Baryon production model
42 //      SetMSTJ(12,3); 
43 // String fragmentation
44 //      SetMSTJ(1,1);  
45 // sea quarks can be used for baryon formation
46 //      SetMSTP(88,2); 
47 // choice of max. virtuality for ISR
48 //      SetMSTP(68,1);
49 // regularisation scheme of ISR 
50 //      SetMSTP(70,2);
51 // all resonance decays switched on
52 //    SetMSTP(41,1);   
53 AliPythia8* AliPythia8::fgAliPythia8=NULL;
54
55 AliPythia8::AliPythia8():
56     AliTPythia8(),
57     AliPythiaBase(),
58     fProcess(kPyMb),
59     fEcms(0.),
60     fStrucFunc(kCTEQ5L),
61     fEtSeed(0.),
62     fMinEtJet(0.),
63     fRJet(0.),
64     fYScale(0.),
65     fPtScale(0.),
66     fNJetMin(0),
67     fNJetMax(0)
68 {
69 // Default Constructor
70 //
71 //  Set random number
72     if (!AliPythiaRndm::GetPythiaRandom()) 
73         AliPythiaRndm::SetPythiaRandom(GetRandom());
74 }
75
76 AliPythia8::AliPythia8(const AliPythia8& pythia):
77     AliTPythia8(), 
78     AliPythiaBase(),
79     fProcess(kPyMb),
80     fEcms(0.),
81     fStrucFunc(kCTEQ5L),
82     fEtSeed(0.),
83     fMinEtJet(0.),
84     fRJet(0.),
85     fYScale(0.),
86     fPtScale(0.),
87     fNJetMin(0),
88     fNJetMax(0)
89 {
90     // Copy Constructor
91     pythia.Copy(*this);
92 }
93
94 void AliPythia8::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfunc)
95 {
96 // Initialise the process to generate 
97     if (!AliPythiaRndm::GetPythiaRandom()) 
98       AliPythiaRndm::SetPythiaRandom(GetRandom());
99     
100     fProcess   = process;
101     fEcms      = energy;
102     fStrucFunc = strucfunc;
103 //...Switch off decay of pi0, K0S, Lambda, Sigma+-, Xi0-, Omega-.
104     ReadString("111:mayDecay  = off");
105     ReadString("310:mayDecay  = off");
106     ReadString("3122:mayDecay = off");
107     ReadString("3112:mayDecay = off");
108     ReadString("3212:mayDecay = off");
109     ReadString("3222:mayDecay = off");
110     ReadString("3312:mayDecay = off");
111     ReadString("3322:mayDecay = off");
112     ReadString("3334:mayDecay = off");
113     // Select structure function 
114
115     ReadString("PDF:useLHAPDF = on");
116     ReadString(Form("PDF:LHAPDFset = %s", AliStructFuncType::PDFsetName(fStrucFunc).Data()));
117
118     // Particles produced in string fragmentation point directly to either of the two endpoints
119     // of the string (depending in the side they were generated from).
120
121 //    SetMSTU(16,2); // ????
122
123 //
124 // Pythia initialisation for selected processes//
125 //
126     switch (process) 
127     {
128     case kPyOldUEQ2ordered:  //Old underlying events with Q2 ordered QCD processes
129 //        Multiple interactions on.
130         ReadString("PartonLevel:MI = on");
131 // Double Gaussian matter distribution.
132         ReadString("MultipleInteractions:bProfile = 2");
133         ReadString("MultipleInteractions:coreFraction = 0.5");
134         ReadString("MultipleInteractions:coreRadius = 0.4");
135 //  pT0.
136         ReadString("MultipleInteractions:pTmin = 2.0");
137 //  Reference energy for pT0 and energy rescaling pace.
138         ReadString("MultipleInteractions:ecmRef = 1800.");
139         ReadString("MultipleInteractions:ecmPow = 0.25");
140 //  String drawing almost completely minimizes string length.
141 //      SetPARP(85,0.9);
142 //      SetPARP(86,0.95);
143 // ISR and FSR activity.
144 // Q^2 scale of the hard scattering
145         ReadString("SigmaProcess:factorMultFac = 4.");
146 // Lambda_FSR scale.
147 //      SetPARJ(81, 0.29);
148         break;
149     case kPyOldUEQ2ordered2:   
150 // Old underlying events with Q2 ordered QCD processes
151 // Multiple interactions on.
152         ReadString("PartonLevel:MI = on");
153 // Double Gaussian matter distribution.
154         ReadString("MultipleInteractions:bProfile = 2");
155         ReadString("MultipleInteractions:coreFraction = 0.5");
156         ReadString("MultipleInteractions:coreRadius = 0.4");
157 // pT0.
158         ReadString("MultipleInteractions:pTmin = 2.0");
159 // Reference energy for pT0 and energy rescaling pace.
160         ReadString("MultipleInteractions:ecmRef = 1800.");
161         ReadString("MultipleInteractions:ecmPow = 0.16");
162 // String drawing almost completely minimizes string length.
163 //      SetPARP(85,0.9);
164 //      SetPARP(86,0.95);
165 // ISR and FSR activity.
166         ReadString("SigmaProcess:factorMultFac = 4.");
167 // Lambda_FSR scale.
168 //      SetPARJ(81,0.29);       
169         break;
170     case kPyOldPopcorn:  
171 // Old production mechanism: Old Popcorn
172         ReadString("HardQCD:all = on");
173 //      SetMSTJ(12,3); 
174 // (D=2) Like MSTJ(12)=2 but added prod ofthe 1er rank baryon
175 //      SetMSTP(88,2); 
176 // (D=1)see can be used to form  baryons (BARYON JUNCTION)
177 //      SetMSTJ(1,1);  
178         AtlasTuning();
179         break;
180     case kPyCharm:
181         ReadString("HardQCD:gg2ccbar = on");
182         ReadString("HardQCD:qqbar2ccbar = on");
183 //  heavy quark masses
184         ReadString("ParticleData:mcRun = 1.2");
185 //
186 //    primordial pT
187         ReadString("BeamRemnants:primordialKT = on");
188         ReadString("BeamRemnants:primordialKTsoft = 0.");
189         ReadString("BeamRemnants:primordialKThard = 1.");
190         ReadString("BeamRemnants:halfScaleForKT = 0.");
191         ReadString("BeamRemnants:halfMassForKT = 0.");
192         break;
193     case kPyBeauty:
194         ReadString("HardQCD:gg2bbbar = on");
195         ReadString("HardQCD:qqbar2bbbar = on");
196         ReadString("ParticleData:mbRun = 4.75");
197         break;
198     case kPyJpsi:
199 // gg->J/Psi g
200         ReadString("Charmonium:gg2QQbar[3S1(1)]g = on");
201         break;
202     case kPyJpsiChi:
203         ReadString("Charmonium:all = on");
204         break;
205     case kPyCharmUnforced:
206 // gq->qg   
207         ReadString("HardQCD:gq2qg = on");
208 // gg->qq
209         ReadString("HardQCD:gg2qq = on");
210 // gg->gg
211         ReadString("HardQCD:gg2gg = on");
212         break;
213     case kPyBeautyUnforced:
214 // gq->qg   
215         ReadString("HardQCD:gq2qg = on");
216 // gg->qq
217         ReadString("HardQCD:gg2qq = on");
218 // gg->gg
219         ReadString("HardQCD:gg2gg = on");
220         break;
221     case kPyMb:
222 // Minimum Bias pp-Collisions
223 //
224 //   
225 //      select Pythia min. bias model
226 // single diffraction AB-->XB
227         ReadString("SoftQCD:minBias = on");
228         ReadString("SoftQCD:singleDiffractive = on");
229         ReadString("SoftQCD:doubleDiffractive = on");
230         AtlasTuning();
231         break;
232     case kPyMbDefault:
233 // Minimum Bias pp-Collisions
234 //
235 //   
236 //      select Pythia min. bias model
237         ReadString("SoftQCD:minBias = on");
238         ReadString("SoftQCD:singleDiffractive = on");
239         ReadString("SoftQCD:doubleDiffractive = on");
240         break;
241     case kPyLhwgMb:
242 // Les Houches Working Group 05 Minimum Bias pp-Collisions: hep-ph/0604120
243 //  -> Pythia 6.3 or above is needed
244 //   
245         ReadString("SoftQCD:minBias = on");
246         ReadString("SoftQCD:singleDiffractive = on");
247         ReadString("SoftQCD:doubleDiffractive = on");
248         ReadString(Form("PDF:LHAPDFset = %s", AliStructFuncType::PDFsetName(kCTEQ6ll).Data()));
249
250 //      SetMSTP(68,1);
251 //      SetMSTP(70,2);
252 //      ReadString("PartonLevel:MI = on");
253 // Double Gaussian matter distribution.
254         ReadString("MultipleInteractions:bProfile = 2");
255         ReadString("MultipleInteractions:coreFraction = 0.5");
256         ReadString("MultipleInteractions:coreRadius = 0.5");
257         ReadString("MultipleInteractions:expPow = 0.16");
258         ReadString("MultipleInteractions:pTmin = 2.3");
259 //      SetMSTP(88,1);
260 //      SetPARP(85,0.9);           // Regulates gluon prod. mechanism
261         break;
262     case kPyMbNonDiffr:
263 // Minimum Bias pp-Collisions
264 //
265 //   
266 //      select Pythia min. bias model
267         ReadString("SoftQCD:minBias = on");
268         AtlasTuning();
269         break;
270     case kPyMbMSEL1:
271         ConfigHeavyFlavor();
272 // Intrinsic <kT^2>
273         ReadString("BeamRemnants:primordialKT = on");
274         ReadString("BeamRemnants:primordialKTsoft = 0.");
275         ReadString("BeamRemnants:primordialKThard = 1.");
276         ReadString("BeamRemnants:halfScaleForKT = 0.");
277         ReadString("BeamRemnants:halfMassForKT = 0.");
278 // Set Q-quark mass
279         ReadString("ParticleData:mcRun = 1.20");
280         ReadString("ParticleData:mbRun = 4.78");
281 // Atlas Tuning
282         AtlasTuning();
283         break;
284     case kPyJets:
285 //
286 //  QCD Jets
287 //
288         ReadString("HardQCD:all = on");
289 //
290 //  Pythia Tune A (CDF)
291 //
292         ReadString("PartonLevel:MI = on");
293         ReadString("MultipleInteractions:pTmin = 2.0");
294         ReadString("MultipleInteractions:pT0Ref = 2.8");
295         ReadString("MultipleInteractions:ecmRef = 1800.");
296         ReadString("MultipleInteractions:expPow = 0.25");
297         ReadString("MultipleInteractions:bProfile = 2");
298         ReadString("MultipleInteractions:coreFraction = 0.16");
299         ReadString("MultipleInteractions:coreRadius = 0.4");
300         ReadString("SigmaProcess:factorMultFac = 2.5");
301 //      SetPARP(85,0.90) ;         // Regulates gluon prod. mechanism
302 //      SetPARP(86,0.95);          // Regulates gluon prod. mechanism
303        break;
304     case kPyDirectGamma:
305         ReadString("PromptPhoton:all = on");
306         break;
307     case kPyCharmPbPbMNR:
308     case kPyD0PbPbMNR:
309     case kPyDPlusPbPbMNR:
310     case kPyDPlusStrangePbPbMNR:
311       // Tuning of Pythia parameters aimed to get a resonable agreement
312       // between with the NLO calculation by Mangano, Nason, Ridolfi for the
313       // c-cbar single inclusive and double differential distributions.
314       // This parameter settings are meant to work with Pb-Pb collisions
315       // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
316       // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
317       // has to be set to 2.1GeV. Example in ConfigCharmPPR.C.
318         ConfigHeavyFlavor();
319       // Intrinsic <kT>
320         ReadString("BeamRemnants:primordialKT = on");
321         ReadString("BeamRemnants:primordialKTsoft = 0.");
322         ReadString("BeamRemnants:primordialKThard = 1.304");
323         ReadString("BeamRemnants:halfScaleForKT = 0.");
324         ReadString("BeamRemnants:halfMassForKT = 0.");
325       // Set c-quark mass
326         ReadString("ParticleData:mcRun = 1.20");
327         break;
328     case kPyCharmpPbMNR:
329     case kPyD0pPbMNR:
330     case kPyDPluspPbMNR:
331     case kPyDPlusStrangepPbMNR:
332       // Tuning of Pythia parameters aimed to get a resonable agreement
333       // between with the NLO calculation by Mangano, Nason, Ridolfi for the
334       // c-cbar single inclusive and double differential distributions.
335       // This parameter settings are meant to work with p-Pb collisions
336       // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
337       // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
338       // has to be set to 2.1GeV. Example in ConfigCharmPPR.C.
339         ConfigHeavyFlavor();
340       // Intrinsic <kT>
341         ReadString("BeamRemnants:primordialKT = on");
342         ReadString("BeamRemnants:primordialKTsoft = 0.");
343         ReadString("BeamRemnants:primordialKThard = 1.16");
344         ReadString("BeamRemnants:halfScaleForKT = 0.");
345         ReadString("BeamRemnants:halfMassForKT = 0.");
346       // Set c-quark mass
347         ReadString("ParticleData:mcRun = 1.20");
348       break;
349     case kPyCharmppMNR:
350     case kPyD0ppMNR:
351     case kPyDPlusppMNR:
352     case kPyDPlusStrangeppMNR:
353     case kPyLambdacppMNR:
354       // Tuning of Pythia parameters aimed to get a resonable agreement
355       // between with the NLO calculation by Mangano, Nason, Ridolfi for the
356       // c-cbar single inclusive and double differential distributions.
357       // This parameter settings are meant to work with pp collisions
358       // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
359       // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
360       // has to be set to 2.1GeV. Example in ConfigCharmPPR.C.
361         ConfigHeavyFlavor();
362       // Intrinsic <kT^2>
363         ReadString("BeamRemnants:primordialKT = on");
364         ReadString("BeamRemnants:primordialKTsoft = 0.");
365         ReadString("BeamRemnants:primordialKThard = 1.");
366         ReadString("BeamRemnants:halfScaleForKT = 0.");
367         ReadString("BeamRemnants:halfMassForKT = 0.");
368       // Set c-quark mass
369         ReadString("ParticleData:mcRun = 1.20");
370       break;
371     case kPyCharmppMNRwmi:
372       // Tuning of Pythia parameters aimed to get a resonable agreement
373       // between with the NLO calculation by Mangano, Nason, Ridolfi for the
374       // c-cbar single inclusive and double differential distributions.
375       // This parameter settings are meant to work with pp collisions
376       // and with kCTEQ5L PDFs.
377       // Added multiple interactions according to ATLAS tune settings.
378       // To get a "reasonable" agreement with MNR results, events have to be 
379       // generated with the minimum ptHard (AliGenPythia::SetPtHard)
380       // set to 2.76 GeV.
381       // To get a "perfect" agreement with MNR results, events have to be 
382       // generated in four ptHard bins with the following relative 
383       // normalizations:
384       // 2.76-3 GeV: 25%
385       //    3-4 GeV: 40%
386       //    4-8 GeV: 29%
387       //     >8 GeV:  6%
388         ConfigHeavyFlavor();
389       // Intrinsic <kT^2>
390         ReadString("BeamRemnants:primordialKT = on");
391         ReadString("BeamRemnants:primordialKTsoft = 0.");
392         ReadString("BeamRemnants:primordialKThard = 1.");
393         ReadString("BeamRemnants:halfScaleForKT = 0.");
394         ReadString("BeamRemnants:halfMassForKT = 0.");
395       // Set c-quark mass
396         ReadString("ParticleData:mcRun = 1.20");
397         AtlasTuning();
398         break;
399     case kPyBeautyPbPbMNR:
400       // Tuning of Pythia parameters aimed to get a resonable agreement
401       // between with the NLO calculation by Mangano, Nason, Ridolfi for the
402       // b-bbar single inclusive and double differential distributions.
403       // This parameter settings are meant to work with Pb-Pb collisions
404       // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
405       // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
406       // has to be set to 2.75GeV. Example in ConfigBeautyPPR.C.
407         ConfigHeavyFlavor();
408       // QCD scales
409         ReadString("SigmaProcess:factorMultFac = 1.");
410       // Intrinsic <kT>
411         ReadString("BeamRemnants:primordialKT = on");
412         ReadString("BeamRemnants:primordialKTsoft = 0.");
413         ReadString("BeamRemnants:primordialKThard = 2.035");
414         ReadString("BeamRemnants:halfScaleForKT = 0.");
415         ReadString("BeamRemnants:halfMassForKT = 0.");
416       // Set b-quark mass
417         ReadString("ParticleData:mbRun = 4.75");
418       break;
419     case kPyBeautypPbMNR:
420       // Tuning of Pythia parameters aimed to get a resonable agreement
421       // between with the NLO calculation by Mangano, Nason, Ridolfi for the
422       // b-bbar single inclusive and double differential distributions.
423       // This parameter settings are meant to work with p-Pb collisions
424       // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
425       // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
426       // has to be set to 2.75GeV. Example in ConfigBeautyPPR.C.
427         ConfigHeavyFlavor();
428       // QCD scales
429         ReadString("SigmaProcess:factorMultFac = 1.");
430       // Intrinsic <kT>
431         ReadString("BeamRemnants:primordialKT = on");
432         ReadString("BeamRemnants:primordialKTsoft = 0.");
433         ReadString("BeamRemnants:primordialKThard = 1.6");
434         ReadString("BeamRemnants:halfScaleForKT = 0.");
435         ReadString("BeamRemnants:halfMassForKT = 0.");
436       // Set b-quark mass
437         ReadString("ParticleData:mbRun = 4.75");
438       break;
439     case kPyBeautyppMNR:
440       // Tuning of Pythia parameters aimed to get a resonable agreement
441       // between with the NLO calculation by Mangano, Nason, Ridolfi for the
442       // b-bbar single inclusive and double differential distributions.
443       // This parameter settings are meant to work with pp collisions
444       // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
445       // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
446       // has to be set to 2.75GeV. Example in ConfigBeautyPPR.C.
447         ConfigHeavyFlavor();
448         // QCD scales
449         ReadString("SigmaProcess:factorMultFac = 1.");
450         // Intrinsic <kT>
451         ReadString("BeamRemnants:primordialKT = on");
452         ReadString("BeamRemnants:primordialKTsoft = 0.");
453         ReadString("BeamRemnants:primordialKThard = 1.0");
454         ReadString("BeamRemnants:halfScaleForKT = 0.");
455         ReadString("BeamRemnants:halfMassForKT = 0.");
456         // Set b-quark mass
457         ReadString("ParticleData:mbRun = 4.75");
458       break;
459      case kPyBeautyppMNRwmi:
460       // Tuning of Pythia parameters aimed to get a resonable agreement
461       // between with the NLO calculation by Mangano, Nason, Ridolfi for the
462       // b-bbar single inclusive and double differential distributions.
463       // This parameter settings are meant to work with pp collisions
464       // and with kCTEQ5L PDFs.
465       // Added multiple interactions according to ATLAS tune settings.
466       // To get a "reasonable" agreement with MNR results, events have to be 
467       // generated with the minimum ptHard (AliGenPythia::SetPtHard)
468       // set to 2.76 GeV.
469       // To get a "perfect" agreement with MNR results, events have to be 
470       // generated in four ptHard bins with the following relative 
471       // normalizations:
472       // 2.76-4 GeV:  5% 
473       //    4-6 GeV: 31%
474       //    6-8 GeV: 28%
475       //     >8 GeV: 36%
476          ConfigHeavyFlavor();
477          // QCD scales
478          ReadString("SigmaProcess:factorMultFac = 1.");
479          // Intrinsic <kT>
480         ReadString("BeamRemnants:primordialKT = on");
481         ReadString("BeamRemnants:primordialKTsoft = 0.");
482         ReadString("BeamRemnants:primordialKThard = 1.0");
483         ReadString("BeamRemnants:halfScaleForKT = 0.");
484         ReadString("BeamRemnants:halfMassForKT = 0.");
485         // Set b-quark mass
486         ReadString("ParticleData:mbRun = 4.75");
487         AtlasTuning();
488         break; 
489     case kPyW:
490         //Inclusive production of W+/-
491         //f fbar -> W+ 
492         ReadString("WeakSingleBoson:ffbar2W = on");
493         // Initial/final parton shower on (Pythia default)
494         // With parton showers on we are generating "W inclusive process"
495         ReadString("PartonLevel:ISR = on");
496         ReadString("PartonLevel:FSR = on");
497         break;  
498     case kPyZ:
499         //Inclusive production of Z
500         //f fbar -> Z/gamma
501         ReadString("WeakSingleBoson:ffbar2gmZ = on");
502         //only Z included, not gamma 
503         ReadString("WeakZ0:gmZmode = 2");
504         // Initial/final parton shower on (Pythia default)
505         // With parton showers on we are generating "Z inclusive process"
506         ReadString("PartonLevel:ISR = on");
507         ReadString("PartonLevel:FSR = on");
508         break;  
509     }
510 //
511 //  Initialize PYTHIA
512 //    SetMSTP(41,1);   // all resonance decays switched on
513     Initialize(2212, 2212, fEcms);
514 }
515
516 void AliPythia8::SetNuclei(Int_t /*a1*/, Int_t /*a2*/)
517 {
518 // Treat protons as inside nuclei with mass numbers a1 and a2  
519 //    The MSTP array in the PYPARS common block is used to enable and 
520 //    select the nuclear structure functions. 
521 //    MSTP(52)  : (D=1) choice of proton and nuclear structure-function library
522 //            =1: internal PYTHIA acording to MSTP(51) 
523 //            =2: PDFLIB proton  s.f., with MSTP(51)  = 1000xNGROUP+NSET
524 //    If the following mass number both not equal zero, nuclear corrections of the stf are used.
525 //    MSTP(192) : Mass number of nucleus side 1
526 //    MSTP(193) : Mass number of nucleus side 2
527 //    SetMSTP(52,2);
528 //    SetMSTP(192, a1);
529 //    SetMSTP(193, a2);  
530 }
531         
532
533 AliPythia8* AliPythia8::Instance()
534
535 // Set random number generator 
536     if (fgAliPythia8) {
537         return fgAliPythia8;
538     } else {
539         fgAliPythia8 = new AliPythia8();
540         return fgAliPythia8;
541     }
542 }
543
544 void AliPythia8::PrintParticles()
545
546 // Print list of particl properties
547     ReadString("Main:showAllParticleData");
548 }
549
550 void  AliPythia8::ResetDecayTable()
551 {
552 //  Set default values for pythia decay switches
553 //    Int_t i;
554 //    for (i = 1; i <  501; i++) SetMDCY(i,1,fDefMDCY[i]);
555 //    for (i = 1; i < 2001; i++) SetMDME(i,1,fDefMDME[i]);
556 }
557
558 void  AliPythia8::SetDecayTable()
559 {
560 //  Set default values for pythia decay switches
561 //
562 //    Int_t i;
563 //    for (i = 1; i <  501; i++) fDefMDCY[i] = GetMDCY(i,1);
564 //    for (i = 1; i < 2001; i++) fDefMDME[i] = GetMDME(i,1);
565 }
566
567 void  AliPythia8::Pyclus(Int_t& njet)
568 {
569 //  Call Pythia clustering algorithm
570 //
571     Bool_t ok = fClusterJet.analyze(Pythia8()->event, fYScale, fPtScale, fNJetMin, fNJetMax);
572     njet = 0;
573     if (ok) njet = fClusterJet.size();
574 }
575
576 void  AliPythia8::Pycell(Int_t& njet)
577 {
578 //  Call Pythia jet reconstruction algorithm
579 //
580     Bool_t ok = fCellJet.analyze(Pythia8()->event, fMinEtJet, fRJet, fEtSeed);
581     njet = 0;
582     if (ok) njet = fCellJet.size();
583 }
584
585 void AliPythia8::GetJet(Int_t i, Float_t& px, Float_t& py, Float_t& pz, Float_t& e)
586 {
587     // Get jet number i
588     Float_t et = fCellJet.eT(i);
589     px = et * TMath::Cos(fCellJet.phiWeighted(i));
590     py = et * TMath::Sin(fCellJet.phiWeighted(i));
591     pz = et * TMath::SinH(fCellJet.etaWeighted(i));
592     e  = et * TMath::CosH(fCellJet.etaWeighted(i));
593 }
594
595 void AliPythia8::GenerateEvent()
596 {
597     // Generate one event
598     AliTPythia8::GenerateEvent();
599 }
600
601 void AliPythia8::GenerateMIEvent()
602 {
603     // New multiple interaction scenario
604     AliWarning("Not implemented. No event will be generated");
605 }
606
607 void AliPythia8::PrintStatistics()
608 {
609     // End of run statistics
610     AliTPythia8::PrintStatistics();
611 }
612
613 void AliPythia8::EventListing()
614 {
615     // End of run statistics
616     AliTPythia8::EventListing();
617 }
618
619 Int_t AliPythia8::ProcessCode()
620 {
621     // Returns the subprocess code for the current event
622     return Pythia8()->info.code();
623 }
624
625 void AliPythia8::ConfigHeavyFlavor()
626 {
627     //
628     // Default configuration for Heavy Flavor production
629     //
630     // All QCD processes
631     //
632     ReadString("HardQCD:all = on");
633
634     // No multiple interactions
635     ReadString("PartonLevel:MI = off");
636     ReadString("MultipleInteractions:pTmin = 0.0");
637     ReadString("MultipleInteractions:pT0Ref = 0.0");
638
639     // Initial/final parton shower on (Pythia default)
640     ReadString("PartonLevel:ISR = on");
641     ReadString("PartonLevel:FSR = on");
642     
643     // 2nd order alpha_s
644     ReadString("SigmaProcess:alphaSorder = 2");
645
646     // QCD scales 
647     ReadString("SigmaProcess:renormScale2 = 2");
648     ReadString("SigmaProcess:renormMultFac = 1.");
649 }
650
651 void AliPythia8::AtlasTuning()
652 {
653     //
654     // Configuration for the ATLAS tuning
655     ReadString(Form("PDF:LHAPDFset = %s", AliStructFuncType::PDFsetName(kCTEQ5L).Data()));
656     ReadString("PartonLevel:MI = on");
657     ReadString("MultipleInteractions:pTmin = 1.9");
658     ReadString("MultipleInteractions:pT0Ref = 1.8");
659     ReadString("MultipleInteractions:ecmRef = 1000.");
660     ReadString("MultipleInteractions:expPow = 0.16");
661     ReadString("MultipleInteractions:bProfile = 2");
662     ReadString("MultipleInteractions:coreFraction = 0.16");
663     ReadString("MultipleInteractions:coreRadius = 0.5");
664 //      SetPARP(85,0.33);          // Regulates gluon prod. mechanism
665 //      SetPARP(86,0.66);          // Regulates gluon prod. mechanism
666     ReadString("SigmaProcess:factorMultFac = 1.");
667 }
668
669 void AliPythia8::SetPtHardRange(Float_t ptmin, Float_t ptmax)
670 {
671     // Set the pt hard range
672     ReadString(Form("PhaseSpace:pTHatMin = %13.3f", ptmin));
673     ReadString(Form("PhaseSpace:pTHatMax = %13.3f", ptmax));
674 }
675
676 void AliPythia8::SetYHardRange(Float_t /*ymin*/, Float_t /*ymax*/)
677 {
678     // Set the y hard range
679     printf("YHardRange not implemented in Pythia8 !!!\n");
680     
681 }
682
683
684 void AliPythia8::SetFragmentation(Int_t flag)
685 {
686     // Switch fragmentation on/off
687     if (flag) {
688         ReadString("HadronLevel:Hadronize = on");
689     } else {
690         ReadString("HadronLevel:Hadronize = off");
691     }
692 }
693
694 void AliPythia8::SetInitialAndFinalStateRadiation(Int_t flag1, Int_t flag2)
695 {
696 //  initial state radiation    
697     if (flag1) {
698         ReadString("PartonLevel:ISR = on");
699     } else {
700         ReadString("PartonLevel:ISR = off");
701     }
702 // final state radiation    
703     if (flag2) {
704         ReadString("PartonLevel:FSR = on");
705     } else {
706         ReadString("PartonLevel:FSR = off");
707     }
708 }
709
710 void AliPythia8::SetIntrinsicKt(Float_t kt)
711 {
712         ReadString("BeamRemnants:primordialKT = on");
713         ReadString("BeamRemnants:primordialKTsoft = 0.");
714         ReadString(Form("BeamRemnants:primordialKThard = %13.3f", kt));
715         ReadString("BeamRemnants:halfScaleForKT = 0.");
716         ReadString("BeamRemnants:halfMassForKT = 0.");
717 }
718
719 void AliPythia8::SwitchHFOff()
720 {
721     // Switch off heavy flavor
722     // Maximum number of quark flavours used in pdf 
723     ReadString("PDFinProcess:nQuarkIn = 3");
724     // Maximum number of flavors that can be used in showers
725     ReadString("TimeShower:nGluonToQuark = 3");
726     ReadString("SpaceShower:nQuarkIn = 3");
727     
728
729 }
730
731 void AliPythia8::SetPycellParameters(Float_t etaMax, Int_t nEta, Int_t nPhi,
732                                        Float_t thresh, Float_t etseed, Float_t minet, Float_t r)
733 {
734 // Set pycell parameters
735     fCellJet  = Pythia8::CellJet( etaMax, nEta, nPhi, 2, 0, 0., 0., thresh);
736     fEtSeed   = etseed;
737     fMinEtJet = minet;
738     fRJet     = r;
739 }
740
741 void AliPythia8::ModifiedSplitting()
742 {
743 //
744 // We have to see how to implement this in Pythia8 !!!
745 //
746     // Modified splitting probability as a model for quenching
747 //    SetPARJ(200, 0.8);
748 //    SetMSTJ(41, 1);  // QCD radiation only
749 //    SetMSTJ(42, 2);  // angular ordering
750 //    SetMSTJ(44, 2);  // option to run alpha_s
751 //    SetMSTJ(47, 0);  // No correction back to hard scattering element
752 //    SetMSTJ(50, 0);  // No coherence in first branching
753 //    SetPARJ(82, 1.); // Cut off for parton showers
754 }
755
756     
757 void AliPythia8::InitQuenching(Float_t /*cMin*/, Float_t /*cMax*/, Float_t /*k*/, Int_t /*iECMethod*/, Float_t /*zmax*/, Int_t /*ngmax*/)
758 {
759     //
760     //
761     AliWarning("Not implemented !");
762 }
763
764 void AliPythia8::SwitchHadronisationOff()
765 {
766     // Switch off hadronisation
767     ReadString("HadronLevel:Hadronize = off");
768 }
769
770 void AliPythia8::SwitchHadronisationOn()
771 {
772     // Switch on hadronisarion
773     ReadString("HadronLevel:Hadronize = on");
774 }
775
776
777 void AliPythia8::GetXandQ(Float_t& x1, Float_t& x2, Float_t& q)
778 {
779     // Get x1, x2 and Q for this event
780     
781     q  = Pythia8()->info.QFac();
782     x1 = Pythia8()->info.x1();
783     x2 = Pythia8()->info.x2();
784     
785 }
786
787 Float_t AliPythia8::GetXSection()
788 {
789     // Get the total cross-section
790     return Pythia8()->info.sigmaGen();
791 }
792
793 Float_t AliPythia8::GetPtHard()
794 {
795     // Get the pT hard for this event
796     return Pythia8()->info.pTHat();
797 }
798
799
800
801
802 AliPythia8& AliPythia8::operator=(const  AliPythia8& rhs)
803 {
804 // Assignment operator
805     rhs.Copy(*this);
806     return *this;
807 }
808
809  void AliPythia8::Copy(TObject&) const
810 {
811     //
812     // Copy 
813     //
814     Fatal("Copy","Not implemented!\n");
815 }
816
817 //
818 // To be implemented
819 //
820 void AliPythia8::SetNumberOfParticles(Int_t /*i*/)
821 {
822     AliWarning("Not implemented");
823 }
824
825 void AliPythia8::EditEventList(Int_t /*i*/)
826 {
827     AliWarning("Not implemented");
828 }
829
830 void AliPythia8::Pyquen(Double_t /*a*/, Int_t /*b*/, Double_t /*c*/)
831 {
832     AliWarning("Cannot be used with Pythia8");
833 }
834
835 void AliPythia8::HadronizeEvent()
836 {
837     // Needs access to HadronLevel ?
838     AliWarning("Not yet implemented");
839 }
840
841 void AliPythia8::GetQuenchingParameters(Double_t& /*xp*/, Double_t& /*yp*/, Double_t* /*z[4]*/)
842 {
843     AliWarning("Not yet implemented");
844 }
845
846 void AliPythia8::LoadEvent(AliStack* /*stack*/, Int_t /*flag*/, Int_t /*reHadr*/)
847 {
848     AliWarning("Not yet implemented");
849 }