A bug concerning the calculation of the track length and time-of-flight hypotheses...
[u/mrichter/AliRoot.git] / PYTHIA6 / AliPythia.cxx
CommitLineData
8d2cd130 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/*
17$Log$
adf4d898 18Revision 1.1 2003/03/15 15:00:48 morsch
19Classed imported from EVGEN.
20
8d2cd130 21Revision 1.28 2002/12/09 08:22:56 morsch
22UA1 jet finder (Pycell) for software triggering added.
23
24Revision 1.27 2002/11/15 00:39:37 morsch
25- Correct initialisation of sRandom.
26- QCD Jets with initial and final state gluon radiation is default
27- pt kick for jets default
28- Interface to Pyclus added.
29
30Revision 1.26 2002/11/14 00:37:32 morsch
31Warning message for kPyJets added.
32
33Revision 1.25 2002/10/14 14:55:35 hristov
34Merging the VirtualMC branch to the main development branch (HEAD)
35
36Revision 1.20.6.1 2002/06/10 14:57:41 hristov
37Merged with v3-08-02
38
39Revision 1.24 2002/05/22 13:22:53 morsch
40Process kPyMbNonDiffr added.
41
42Revision 1.23 2002/05/06 07:17:29 morsch
43Pyr gives random number r in interval 0 < r < 1.
44
45Revision 1.22 2002/04/26 10:28:48 morsch
46Option kPyBeautyPbMNR added (N. Carrer).
47
48Revision 1.21 2002/03/25 14:46:16 morsch
49Case kPyD0PbMNR added (N. Carrer).
50
51Revision 1.20 2002/03/03 13:48:50 morsch
52Option kPyCharmPbMNR added. Produce charm pairs in agreement with MNR
53NLO calculations (Nicola Carrer).
54
55Revision 1.19 2002/02/20 08:52:20 morsch
56Correct documentation of SetNuclei method.
57
58Revision 1.18 2002/02/07 10:43:06 morsch
59Tuned pp-min.bias settings (M.Monteno, R.Ugoccioni and N.Carrer)
60
61Revision 1.17 2001/12/19 15:40:43 morsch
62For kPyJets enforce simple jet topology, i.e no initial or final state
63gluon radiation and no primordial pT.
64
65Revision 1.16 2001/10/12 11:13:59 morsch
66Missing break statements added (thanks to Nicola Carrer)
67
68Revision 1.15 2001/03/27 10:54:50 morsch
69Add ResetDecayTable() and SsetDecayTable() methods.
70
71Revision 1.14 2001/03/09 13:03:40 morsch
72Process_t and Struc_Func_t moved to AliPythia.h
73
74Revision 1.13 2000/12/18 08:55:35 morsch
75Make AliPythia dependent generartors work with new scheme of random number generation
76
77Revision 1.12 2000/11/30 07:12:50 alibrary
78Introducing new Rndm and QA classes
79
80Revision 1.11 2000/10/20 06:30:06 fca
81Use version 0 to avoid streamer generation
82
83Revision 1.10 2000/10/06 14:18:44 morsch
84Upper cut of prim. pT distribution set to 5. GeV
85
86Revision 1.9 2000/09/18 10:41:35 morsch
87Add possibility to use nuclear structure functions from PDF library V8.
88
89Revision 1.8 2000/09/06 14:26:24 morsch
90Decayer functionality of AliPythia has been moved to AliDecayerPythia.
91Class is now a singleton.
92
93Revision 1.7 2000/06/09 20:34:50 morsch
94All coding rule violations except RS3 corrected
95
96Revision 1.6 1999/11/09 07:38:48 fca
97Changes for compatibility with version 2.23 of ROOT
98
99Revision 1.5 1999/11/03 17:43:20 fca
100New version from G.Martinez & A.Morsch
101
102Revision 1.4 1999/09/29 09:24:14 fca
103Introduction of the Copyright and cvs Log
104
105*/
106
107
108#include "AliPythia.h"
109
110ClassImp(AliPythia)
111
112#ifndef WIN32
113# define pyclus pyclus_
114# define pycell pycell_
115# define type_of_call
116#else
117# define pyclus PYCLUS
118# define pycell PYCELL
119# define type_of_call _stdcall
120#endif
121
122extern "C" void type_of_call pyclus(Int_t & );
123extern "C" void type_of_call pycell(Int_t & );
124
125//_____________________________________________________________________________
126
127AliPythia* AliPythia::fgAliPythia=NULL;
128
129AliPythia::AliPythia()
130{
131// Default Constructor
132//
133// Set random number
134 if (!sRandom) sRandom=fRandom;
135
136}
137
138void AliPythia::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfunc)
139{
140// Initialise the process to generate
141 if (!sRandom) sRandom = gRandom;
142
143 fProcess = process;
144 fEcms = energy;
145 fStrucFunc = strucfunc;
146// don't decay p0
147 SetMDCY(Pycomp(111),1,0);
148// select structure function
149 SetMSTP(52,2);
150 SetMSTP(51,strucfunc);
151//
152// Pythia initialisation for selected processes//
153//
154// Make MSEL clean
155//
156 for (Int_t i=1; i<= 200; i++) {
157 SetMSUB(i,0);
158 }
159// select charm production
160 switch (process)
161 {
162 case kPyCharm:
163 SetMSEL(4);
164//
165// heavy quark masses
166
167 SetPMAS(4,1,1.2);
168 SetMSTU(16,2);
169//
170// primordial pT
171 SetMSTP(91,1);
172 SetPARP(91,1.);
173 SetPARP(93,5.);
174//
175 break;
176 case kPyBeauty:
177 SetMSEL(5);
178 SetPMAS(5,1,4.75);
179 SetMSTU(16,2);
180 break;
181 case kPyJpsi:
182 SetMSEL(0);
183// gg->J/Psi g
184 SetMSUB(86,1);
185 break;
186 case kPyJpsiChi:
187 SetMSEL(0);
188// gg->J/Psi g
189 SetMSUB(86,1);
190// gg-> chi_0c g
191 SetMSUB(87,1);
192// gg-> chi_1c g
193 SetMSUB(88,1);
194// gg-> chi_2c g
195 SetMSUB(89,1);
196 break;
197 case kPyCharmUnforced:
198 SetMSEL(0);
199// gq->qg
200 SetMSUB(28,1);
201// gg->qq
202 SetMSUB(53,1);
203// gg->gg
204 SetMSUB(68,1);
205 break;
206 case kPyBeautyUnforced:
207 SetMSEL(0);
208// gq->qg
209 SetMSUB(28,1);
210// gg->qq
211 SetMSUB(53,1);
212// gg->gg
213 SetMSUB(68,1);
214 break;
215 case kPyMb:
216// Minimum Bias pp-Collisions
217//
218//
219// select Pythia min. bias model
220 SetMSEL(0);
221 SetMSUB(92,1); // single diffraction AB-->XB
222 SetMSUB(93,1); // single diffraction AB-->AX
223 SetMSUB(94,1); // double diffraction
224 SetMSUB(95,1); // low pt production
225 SetMSTP(81,1); // multiple interactions switched on
226 SetMSTP(82,3); // model with varying impact param. & a single Gaussian
227 SetPARP(82,3.47); // set value pT_0 for turn-off of the cross section of
228 // multiple interaction at a reference energy = 14000 GeV
229 SetPARP(89,14000.); // reference energy for the above parameter
230 SetPARP(90,0.174); // set exponent for energy dependence of pT_0
231 case kPyMbNonDiffr:
232// Minimum Bias pp-Collisions
233//
234//
235// select Pythia min. bias model
236 SetMSEL(0);
237 SetMSUB(95,1); // low pt production
238 SetMSTP(81,1); // multiple interactions switched on
239 SetMSTP(82,3); // model with varying impact param. & a single Gaussian
240 SetPARP(82,3.47); // set value pT_0 for turn-off of the cross section of
241 // multiple interaction at a reference energy = 14000 GeV
242 SetPARP(89,14000.); // reference energy for the above parameter
243 SetPARP(90,0.174); // set exponent for energy dependence of pT_0
244
245 break;
246 case kPyJets:
247//
248// QCD Jets
249//
250 SetMSEL(1);
251 break;
252 case kPyDirectGamma:
253 SetMSEL(10);
254 break;
adf4d898 255 case kPyCharmPbPbMNR:
256 case kPyD0PbPbMNR:
8d2cd130 257 // Tuning of Pythia parameters aimed to get a resonable agreement
258 // between with the NLO calculation by Mangano, Nason, Ridolfi for the
259 // c-cbar single inclusive and double differential distributions.
260 // This parameter settings are meant to work with Pb-Pb collisions
adf4d898 261 // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
8d2cd130 262 // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
263 // has to be set to 2.1GeV. Example in ConfigCharmPPR.C.
264
265 // All QCD processes
266 SetMSEL(1);
267
268 // No multiple interactions
269 SetMSTP(81,0);
270 SetPARP(81,0.0);
271 SetPARP(82,0.0);
272
273 // Initial/final parton shower on (Pythia default)
274 SetMSTP(61,1);
275 SetMSTP(71,1);
276
277 // 2nd order alpha_s
278 SetMSTP(2,2);
279
280 // QCD scales
281 SetMSTP(32,2);
282 SetPARP(34,1.0);
283
adf4d898 284 // Intrinsic <kT>
8d2cd130 285 SetMSTP(91,1);
286 SetPARP(91,1.304);
287 SetPARP(93,6.52);
288
289 // Set c-quark mass
290 SetPMAS(4,1,1.2);
291
292 break;
adf4d898 293 case kPyCharmpPbMNR:
294 case kPyD0pPbMNR:
295 // Tuning of Pythia parameters aimed to get a resonable agreement
296 // between with the NLO calculation by Mangano, Nason, Ridolfi for the
297 // c-cbar single inclusive and double differential distributions.
298 // This parameter settings are meant to work with p-Pb collisions
299 // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
300 // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
301 // has to be set to 2.1GeV. Example in ConfigCharmPPR.C.
302
303 // All QCD processes
304 SetMSEL(1);
305
306 // No multiple interactions
307 SetMSTP(81,0);
308 SetPARP(81,0.0);
309 SetPARP(82,0.0);
310
311 // Initial/final parton shower on (Pythia default)
312 SetMSTP(61,1);
313 SetMSTP(71,1);
314
315 // 2nd order alpha_s
316 SetMSTP(2,2);
317
318 // QCD scales
319 SetMSTP(32,2);
320 SetPARP(34,1.0);
321
322 // Intrinsic <kT>
323 SetMSTP(91,1);
324 SetPARP(91,1.16);
325 SetPARP(93,5.8);
326
327 // Set c-quark mass
328 SetPMAS(4,1,1.2);
329
330 break;
331 case kPyCharmppMNR:
332 case kPyD0ppMNR:
333 // Tuning of Pythia parameters aimed to get a resonable agreement
334 // between with the NLO calculation by Mangano, Nason, Ridolfi for the
335 // c-cbar single inclusive and double differential distributions.
336 // This parameter settings are meant to work with pp collisions
337 // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
338 // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
339 // has to be set to 2.1GeV. Example in ConfigCharmPPR.C.
340
341 // All QCD processes
342 SetMSEL(1);
343
344 // No multiple interactions
345 SetMSTP(81,0);
346 SetPARP(81,0.0);
347 SetPARP(82,0.0);
348
349 // Initial/final parton shower on (Pythia default)
350 SetMSTP(61,1);
351 SetMSTP(71,1);
352
353 // 2nd order alpha_s
354 SetMSTP(2,2);
355
356 // QCD scales
357 SetMSTP(32,2);
358 SetPARP(34,1.0);
359
360 // Intrinsic <kT^2>
361 SetMSTP(91,1);
362 SetPARP(91,1.);
363 SetPARP(93,5.);
364
365 // Set c-quark mass
366 SetPMAS(4,1,1.2);
367
368 break;
369 case kPyBeautyPbPbMNR:
8d2cd130 370 // Tuning of Pythia parameters aimed to get a resonable agreement
371 // between with the NLO calculation by Mangano, Nason, Ridolfi for the
372 // b-bbar single inclusive and double differential distributions.
373 // This parameter settings are meant to work with Pb-Pb collisions
374 // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
375 // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
376 // has to be set to 2.75GeV. Example in ConfigBeautyPPR.C.
377
378 // All QCD processes
379 SetMSEL(1);
380
381 // No multiple interactions
382 SetMSTP(81,0);
383 SetPARP(81,0.0);
384 SetPARP(82,0.0);
385
386 // Initial/final parton shower on (Pythia default)
387 SetMSTP(61,1);
388 SetMSTP(71,1);
389
390 // 2nd order alpha_s
391 SetMSTP(2,2);
392
393 // QCD scales
394 SetMSTP(32,2);
395 SetPARP(34,1.0);
396 SetPARP(67,1.0);
397 SetPARP(71,1.0);
398
adf4d898 399 // Intrinsic <kT>
8d2cd130 400 SetMSTP(91,1);
401 SetPARP(91,2.035);
402 SetPARP(93,10.17);
403
404 // Set b-quark mass
405 SetPMAS(5,1,4.75);
406
adf4d898 407 break;
408 case kPyBeautypPbMNR:
409 // Tuning of Pythia parameters aimed to get a resonable agreement
410 // between with the NLO calculation by Mangano, Nason, Ridolfi for the
411 // b-bbar single inclusive and double differential distributions.
412 // This parameter settings are meant to work with p-Pb collisions
413 // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
414 // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
415 // has to be set to 2.75GeV. Example in ConfigBeautyPPR.C.
416
417 // All QCD processes
418 SetMSEL(1);
419
420 // No multiple interactions
421 SetMSTP(81,0);
422 SetPARP(81,0.0);
423 SetPARP(82,0.0);
424
425 // Initial/final parton shower on (Pythia default)
426 SetMSTP(61,1);
427 SetMSTP(71,1);
428
429 // 2nd order alpha_s
430 SetMSTP(2,2);
431
432 // QCD scales
433 SetMSTP(32,2);
434 SetPARP(34,1.0);
435 SetPARP(67,1.0);
436 SetPARP(71,1.0);
437
438 // Intrinsic <kT>
439 SetMSTP(91,1);
440 SetPARP(91,1.60);
441 SetPARP(93,8.00);
442
443 // Set b-quark mass
444 SetPMAS(5,1,4.75);
445
446 break;
447 case kPyBeautyppMNR:
448 // Tuning of Pythia parameters aimed to get a resonable agreement
449 // between with the NLO calculation by Mangano, Nason, Ridolfi for the
450 // b-bbar single inclusive and double differential distributions.
451 // This parameter settings are meant to work with pp collisions
452 // (AliGenPythia::SetNuclei) and with kCTEQ4L PDFs.
453 // To get a good agreement the minimum ptHard (AliGenPythia::SetPtHard)
454 // has to be set to 2.75GeV. Example in ConfigBeautyPPR.C.
455
456 // All QCD processes
457 SetMSEL(1);
458
459 // No multiple interactions
460 SetMSTP(81,0);
461 SetPARP(81,0.0);
462 SetPARP(82,0.0);
463
464 // Initial/final parton shower on (Pythia default)
465 SetMSTP(61,1);
466 SetMSTP(71,1);
467
468 // 2nd order alpha_s
469 SetMSTP(2,2);
470
471 // QCD scales
472 SetMSTP(32,2);
473 SetPARP(34,1.0);
474 SetPARP(67,1.0);
475 SetPARP(71,1.0);
476
477 // Intrinsic <kT>
478 SetMSTP(91,1);
479 SetPARP(91,1.);
480 SetPARP(93,5.);
481
482 // Set b-quark mass
483 SetPMAS(5,1,4.75);
484
8d2cd130 485 break;
486 }
487//
488// Initialize PYTHIA
489 SetMSTP(41,1); // all resonance decays switched on
490
491 Initialize("CMS","p","p",fEcms);
492
493}
494
495Int_t AliPythia::CheckedLuComp(Int_t kf)
496{
497// Check Lund particle code (for debugging)
498 Int_t kc=Pycomp(kf);
499 printf("\n Lucomp kf,kc %d %d",kf,kc);
500 return kc;
501}
502
503void AliPythia::SetNuclei(Int_t a1, Int_t a2)
504{
505// Treat protons as inside nuclei with mass numbers a1 and a2
506// The MSTP array in the PYPARS common block is used to enable and
507// select the nuclear structure functions.
508// MSTP(52) : (D=1) choice of proton and nuclear structure-function library
509// =1: internal PYTHIA acording to MSTP(51)
510// =2: PDFLIB proton s.f., with MSTP(51) = 1000xNGROUP+NSET
511// If the following mass number both not equal zero, nuclear corrections of the stf are used.
512// MSTP(192) : Mass number of nucleus side 1
513// MSTP(193) : Mass number of nucleus side 2
514 SetMSTP(52,2);
515 SetMSTP(192, a1);
516 SetMSTP(193, a2);
517}
518
519
520AliPythia* AliPythia::Instance()
521{
522// Set random number generator
523 if (fgAliPythia) {
524 return fgAliPythia;
525 } else {
526 fgAliPythia = new AliPythia();
527 return fgAliPythia;
528 }
529}
530
531void AliPythia::PrintParticles()
532{
533// Print list of particl properties
534 Int_t np = 0;
535
536 for (Int_t kf=0; kf<1000000; kf++) {
537 for (Int_t c = 1; c > -2; c-=2) {
538
539 Int_t kc = Pycomp(c*kf);
540 if (kc) {
541 Float_t mass = GetPMAS(kc,1);
542 Float_t width = GetPMAS(kc,2);
543 Float_t tau = GetPMAS(kc,4);
544
545 char* name = new char[8];
546 Pyname(kf,name);
547
548 np++;
549
550 printf("\n mass, width, tau: %6d %s %10.3f %10.3e %10.3e",
551 c*kf, name, mass, width, tau);
552 }
553 }
554 }
555 printf("\n Number of particles %d \n \n", np);
556}
557
558void AliPythia::ResetDecayTable()
559{
560// Set default values for pythia decay switches
561 Int_t i;
562 for (i = 1; i < 501; i++) SetMDCY(i,1,fDefMDCY[i]);
563 for (i = 1; i < 2001; i++) SetMDME(i,1,fDefMDME[i]);
564}
565
566void AliPythia::SetDecayTable()
567{
568// Set default values for pythia decay switches
569//
570 Int_t i;
571 for (i = 1; i < 501; i++) fDefMDCY[i] = GetMDCY(i,1);
572 for (i = 1; i < 2001; i++) fDefMDME[i] = GetMDME(i,1);
573}
574
575void AliPythia::Pyclus(Int_t& njet)
576{
577// Call Pythia clustering algorithm
578//
579 pyclus(njet);
580}
581
582void AliPythia::Pycell(Int_t& njet)
583{
584// Call Pythia jet reconstruction algorithm
585//
586 pycell(njet);
587}
588
589
590
591#ifndef WIN32
592#define pyr pyr_
593#define pyrset pyrset_
594#define pyrget pyrget_
595#define pyclus pyclus_
596#define pycell pycell_
597#else
598#define pyr PYR
599#define pyrset PYRSET
600#define pyrget PYRGET
601#define pyclus PYCLUS
602#define pycell PYCELL
603#endif
604
605extern "C" {
606 Double_t pyr(Int_t*)
607{
608 Float_t r;
609 do r=sRandom->Rndm(); while(0 >= r || r >= 1);
610 return r;
611}
612 void pyrset(Int_t*,Int_t*) {}
613 void pyrget(Int_t*,Int_t*) {}
614}
615
616
617
618
619