1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
19 // Generator using the TPythia interface (via AliPythia)
20 // to generate pp collisions.
21 // Using SetNuclei() also nuclear modifications to the structure functions
22 // can be taken into account. This makes, of course, only sense for the
23 // generation of the products of hard processes (heavy flavor, jets ...)
25 // andreas.morsch@cern.ch
28 #include <TClonesArray.h>
29 #include <TDatabasePDG.h>
30 #include <TParticle.h>
32 #include <TObjArray.h>
36 #include "AliDecayerPythia.h"
37 #include "AliGenPythia.h"
38 #include "AliFastGlauber.h"
39 #include "AliHeader.h"
40 #include "AliGenPythiaEventHeader.h"
41 #include "AliPythia.h"
42 #include "AliPythiaRndm.h"
45 #include "AliRunLoader.h"
48 #include "PyquenCommon.h"
51 ClassImp(AliGenPythia)
54 AliGenPythia::AliGenPythia():
95 fDecayer(new AliDecayerPythia()),
103 fPhiMaxJet(2.* TMath::Pi()),
104 fJetReconstruction(kCell),
108 fPhiMaxGamma(2. * TMath::Pi()),
112 fPycellThreshold(0.),
114 fPycellMinEtJet(10.),
115 fPycellMaxRadius(1.),
116 fStackFillOpt(kFlavorSelection),
118 fFragmentation(kTRUE),
127 fTriggerMultiplicity(0),
128 fTriggerMultiplicityEta(0),
129 fTriggerMultiplicityPtMin(0),
130 fCountMode(kCountAll),
134 fFragPhotonInCalo(kFALSE),
136 fPhotonInCalo(kFALSE),
140 fCheckPHOSeta(kFALSE),
141 fFragPhotonOrPi0MinPt(0),
153 // Default Constructor
155 if (!AliPythiaRndm::GetPythiaRandom())
156 AliPythiaRndm::SetPythiaRandom(GetRandom());
159 AliGenPythia::AliGenPythia(Int_t npart)
171 fInteractionRate(0.),
185 fHadronisation(kTRUE),
186 fPatchOmegaDalitz(0),
188 fReadFromFile(kFALSE),
200 fDecayer(new AliDecayerPythia()),
201 fDebugEventFirst(-1),
208 fPhiMaxJet(2.* TMath::Pi()),
209 fJetReconstruction(kCell),
213 fPhiMaxGamma(2. * TMath::Pi()),
217 fPycellThreshold(0.),
219 fPycellMinEtJet(10.),
220 fPycellMaxRadius(1.),
221 fStackFillOpt(kFlavorSelection),
223 fFragmentation(kTRUE),
232 fTriggerMultiplicity(0),
233 fTriggerMultiplicityEta(0),
234 fTriggerMultiplicityPtMin(0),
235 fCountMode(kCountAll),
239 fFragPhotonInCalo(kFALSE),
241 fPhotonInCalo(kFALSE),
245 fCheckPHOSeta(kFALSE),
246 fFragPhotonOrPi0MinPt(0),
258 // default charm production at 5. 5 TeV
260 // structure function GRVHO
264 fTitle= "Particle Generator using PYTHIA";
266 // Set random number generator
267 if (!AliPythiaRndm::GetPythiaRandom())
268 AliPythiaRndm::SetPythiaRandom(GetRandom());
271 AliGenPythia::~AliGenPythia()
274 if(fEventsTime) delete fEventsTime;
277 void AliGenPythia::SetInteractionRate(Float_t rate,Float_t timewindow)
279 // Generate pileup using user specified rate
280 fInteractionRate = rate;
281 fTimeWindow = timewindow;
285 void AliGenPythia::GeneratePileup()
287 // Generate sub events time for pileup
289 if(fInteractionRate == 0.) {
290 Warning("GeneratePileup","Zero interaction specified. Skipping pileup generation.\n");
294 Int_t npart = NumberParticles();
296 Warning("GeneratePileup","Negative number of particles. Skipping pileup generation.\n");
300 if(fEventsTime) delete fEventsTime;
301 fEventsTime = new TArrayF(npart);
302 TArrayF &array = *fEventsTime;
303 for(Int_t ipart = 0; ipart < npart; ipart++)
306 Float_t eventtime = 0.;
309 eventtime += (AliPythiaRndm::GetPythiaRandom())->Exp(1./fInteractionRate);
310 if(eventtime > fTimeWindow) break;
311 array.Set(array.GetSize()+1);
312 array[array.GetSize()-1] = eventtime;
318 eventtime -= (AliPythiaRndm::GetPythiaRandom())->Exp(1./fInteractionRate);
319 if(TMath::Abs(eventtime) > fTimeWindow) break;
320 array.Set(array.GetSize()+1);
321 array[array.GetSize()-1] = eventtime;
324 SetNumberParticles(fEventsTime->GetSize());
327 void AliGenPythia::SetPycellParameters(Float_t etamax, Int_t neta, Int_t nphi,
328 Float_t thresh, Float_t etseed, Float_t minet, Float_t r)
330 // Set pycell parameters
331 fPycellEtaMax = etamax;
334 fPycellThreshold = thresh;
335 fPycellEtSeed = etseed;
336 fPycellMinEtJet = minet;
337 fPycellMaxRadius = r;
342 void AliGenPythia::SetEventListRange(Int_t eventFirst, Int_t eventLast)
344 // Set a range of event numbers, for which a table
345 // of generated particle will be printed
346 fDebugEventFirst = eventFirst;
347 fDebugEventLast = eventLast;
348 if (fDebugEventLast==-1) fDebugEventLast=fDebugEventFirst;
351 void AliGenPythia::Init()
355 SetMC(AliPythia::Instance());
356 fPythia=(AliPythia*) fMCEvGen;
359 fParentWeight=1./Float_t(fNpart);
363 fPythia->SetCKIN(3,fPtHardMin);
364 fPythia->SetCKIN(4,fPtHardMax);
365 fPythia->SetCKIN(7,fYHardMin);
366 fPythia->SetCKIN(8,fYHardMax);
368 if (fAProjectile > 0 && fATarget > 0) fPythia->SetNuclei(fAProjectile, fATarget, fNucPdf);
370 if (fFragmentation) {
371 fPythia->SetMSTP(111,1);
373 fPythia->SetMSTP(111,0);
377 // initial state radiation
378 fPythia->SetMSTP(61,fGinit);
379 // final state radiation
380 fPythia->SetMSTP(71,fGfinal);
383 fPythia->SetMSTP(91,1);
384 fPythia->SetPARP(91,fPtKick);
385 fPythia->SetPARP(93, 4. * fPtKick);
387 fPythia->SetMSTP(91,0);
392 fRL = AliRunLoader::Open(fkFileName, "Partons");
393 fRL->LoadKinematics();
399 fPythia->ProcInit(fProcess,fEnergyCMS,fStrucFunc, fItune);
400 // Forward Paramters to the AliPythia object
401 fDecayer->SetForceDecay(fForceDecay);
402 // Switch off Heavy Flavors on request
404 // Maximum number of quark flavours used in pdf
405 fPythia->SetMSTP(58, 3);
406 // Maximum number of flavors that can be used in showers
407 fPythia->SetMSTJ(45, 3);
408 // Switch off g->QQbar splitting in decay table
409 ((AliDecayerPythia*) fDecayer)->HeavyFlavourOff();
415 // Parent and Children Selection
418 case kPyOldUEQ2ordered:
419 case kPyOldUEQ2ordered2:
423 case kPyCharmUnforced:
424 case kPyCharmPbPbMNR:
427 case kPyCharmppMNRwmi:
428 fParentSelect[0] = 411;
429 fParentSelect[1] = 421;
430 fParentSelect[2] = 431;
431 fParentSelect[3] = 4122;
432 fParentSelect[4] = 4232;
433 fParentSelect[5] = 4132;
434 fParentSelect[6] = 4332;
440 fParentSelect[0] = 421;
443 case kPyDPlusPbPbMNR:
446 fParentSelect[0] = 411;
449 case kPyDPlusStrangePbPbMNR:
450 case kPyDPlusStrangepPbMNR:
451 case kPyDPlusStrangeppMNR:
452 fParentSelect[0] = 431;
455 case kPyLambdacppMNR:
456 fParentSelect[0] = 4122;
461 case kPyBeautyPbPbMNR:
462 case kPyBeautypPbMNR:
464 case kPyBeautyppMNRwmi:
465 fParentSelect[0]= 511;
466 fParentSelect[1]= 521;
467 fParentSelect[2]= 531;
468 fParentSelect[3]= 5122;
469 fParentSelect[4]= 5132;
470 fParentSelect[5]= 5232;
471 fParentSelect[6]= 5332;
474 case kPyBeautyUnforced:
475 fParentSelect[0] = 511;
476 fParentSelect[1] = 521;
477 fParentSelect[2] = 531;
478 fParentSelect[3] = 5122;
479 fParentSelect[4] = 5132;
480 fParentSelect[5] = 5232;
481 fParentSelect[6] = 5332;
486 fParentSelect[0] = 443;
489 case kPyMbAtlasTuneMC09:
491 case kPyMbWithDirectPhoton:
504 // JetFinder for Trigger
506 // Configure detector (EMCAL like)
508 fPythia->SetPARU(51, fPycellEtaMax);
509 fPythia->SetMSTU(51, fPycellNEta);
510 fPythia->SetMSTU(52, fPycellNPhi);
512 // Configure Jet Finder
514 fPythia->SetPARU(58, fPycellThreshold);
515 fPythia->SetPARU(52, fPycellEtSeed);
516 fPythia->SetPARU(53, fPycellMinEtJet);
517 fPythia->SetPARU(54, fPycellMaxRadius);
518 fPythia->SetMSTU(54, 2);
520 // This counts the total number of calls to Pyevnt() per run.
535 Warning("Init","SetNuclei used. Use SetProjectile + SetTarget instead. fDyBoost has been reset to 0\n");
538 fPythia->SetPARJ(200, 0.0);
539 fPythia->SetPARJ(199, 0.0);
540 fPythia->SetPARJ(198, 0.0);
541 fPythia->SetPARJ(197, 0.0);
544 fPythia->InitQuenching(0., 0.1, 0.6e6, 0);
547 if(fQuench ==2){fPythia->SetPyquenParameters(fpyquenT,fpyquenTau,fpyquenNf,fpyquenEloss,fpyquenAngle);}
550 // Nestor's change of the splittings
551 fPythia->SetPARJ(200, 0.8);
552 fPythia->SetMSTJ(41, 1); // QCD radiation only
553 fPythia->SetMSTJ(42, 2); // angular ordering
554 fPythia->SetMSTJ(44, 2); // option to run alpha_s
555 fPythia->SetMSTJ(47, 0); // No correction back to hard scattering element
556 fPythia->SetMSTJ(50, 0); // No coherence in first branching
557 fPythia->SetPARJ(82, 1.); // Cut off for parton showers
558 } else if (fQuench == 4) {
559 // Armesto-Cunqueiro-Salgado change of the splittings.
560 AliFastGlauber* glauber = AliFastGlauber::Instance();
562 //read and store transverse almonds corresponding to differnt
564 glauber->SetCentralityClass(0.,0.1);
565 fPythia->SetPARJ(200, 1.);
566 fPythia->SetPARJ(198, fQhat);
567 fPythia->SetPARJ(199, fLength);
568 fPythia->SetMSTJ(42, 2); // angular ordering
569 fPythia->SetMSTJ(44, 2); // option to run alpha_s
570 fPythia->SetPARJ(82, 1.); // Cut off for parton showers
574 void AliGenPythia::Generate()
576 // Generate one event
577 if (!fPythia) fPythia=(AliPythia*) fMCEvGen;
578 fDecayer->ForceDecay();
580 Double_t polar[3] = {0,0,0};
581 Double_t origin[3] = {0,0,0};
583 // converts from mm/c to s
584 const Float_t kconv=0.001/2.999792458e8;
594 // Set collision vertex position
595 if (fVertexSmear == kPerEvent) Vertex();
604 // Switch hadronisation off
606 fPythia->SetMSTJ(1, 0);
610 // Quenching comes through medium-modified splitting functions.
611 AliFastGlauber::Instance()->GetRandomBHard(bimp);
612 fPythia->SetPARJ(197, bimp);
617 // Either produce new event or read partons from file
619 if (!fReadFromFile) {
625 fNpartons = fPythia->GetN();
627 printf("Loading Event %d\n",AliRunLoader::Instance()->GetEventNumber());
628 fRL->GetEvent(AliRunLoader::Instance()->GetEventNumber());
630 LoadEvent(fRL->Stack(), 0 , 1);
635 // Run quenching routine
639 } else if (fQuench == 2){
640 fPythia->Pyquen(208., 0, 0.);
641 } else if (fQuench == 3) {
642 // Quenching is via multiplicative correction of the splittings
646 // Switch hadronisation on
648 if (fHadronisation) {
649 fPythia->SetMSTJ(1, 1);
651 // .. and perform hadronisation
652 // printf("Calling hadronisation %d\n", fPythia->GetN());
654 if (fPatchOmegaDalitz) {
655 fPythia->SetMDCY(fPythia->Pycomp(111) ,1, 0);
657 fPythia->DalitzDecays();
658 fPythia->SetMDCY(fPythia->Pycomp(111) ,1, 1);
663 fPythia->ImportParticles(&fParticles,"All");
665 if (TMath::Abs(fDyBoost) > 1.e-4) Boost();
666 if(TMath::Abs(fXingAngleY) > 1.e-10) BeamCrossAngle();
674 Int_t np = fParticles.GetEntriesFast();
676 if (np == 0) continue;
680 Int_t* pParent = new Int_t[np];
681 Int_t* pSelected = new Int_t[np];
682 Int_t* trackIt = new Int_t[np];
683 for (i = 0; i < np; i++) {
689 Int_t nc = 0; // Total n. of selected particles
690 Int_t nParents = 0; // Selected parents
691 Int_t nTkbles = 0; // Trackable particles
692 if (fProcess != kPyMbDefault &&
694 fProcess != kPyMbAtlasTuneMC09 &&
695 fProcess != kPyMbWithDirectPhoton &&
696 fProcess != kPyJets &&
697 fProcess != kPyDirectGamma &&
698 fProcess != kPyMbNonDiffr &&
699 fProcess != kPyMbMSEL1 &&
702 fProcess != kPyCharmppMNRwmi &&
703 fProcess != kPyBeautyppMNRwmi &&
704 fProcess != kPyBeautyJets) {
706 for (i = 0; i < np; i++) {
707 TParticle* iparticle = (TParticle *) fParticles.At(i);
708 Int_t ks = iparticle->GetStatusCode();
709 kf = CheckPDGCode(iparticle->GetPdgCode());
710 // No initial state partons
711 if (ks==21) continue;
713 // Heavy Flavor Selection
720 if (kfl > 100000) kfl %= 100000;
721 if (kfl > 10000) kfl %= 10000;
723 if (kfl > 10) kfl/=100;
725 if (kfl > 10) kfl/=10;
726 Int_t ipa = iparticle->GetFirstMother()-1;
729 // Establish mother daughter relation between heavy quarks and mesons
731 if (kf >= fFlavorSelect && kf <= 6) {
732 Int_t idau = iparticle->GetFirstDaughter() - 1;
734 TParticle* daughter = (TParticle *) fParticles.At(idau);
735 Int_t pdgD = daughter->GetPdgCode();
736 if (pdgD == 91 || pdgD == 92) {
737 Int_t jmin = daughter->GetFirstDaughter() - 1;
738 Int_t jmax = daughter->GetLastDaughter() - 1;
739 for (Int_t jp = jmin; jp <= jmax; jp++)
740 ((TParticle *) fParticles.At(jp))->SetFirstMother(i+1);
741 } // is string or cluster
747 TParticle * mother = (TParticle *) fParticles.At(ipa);
748 kfMo = TMath::Abs(mother->GetPdgCode());
751 // What to keep in Stack?
752 Bool_t flavorOK = kFALSE;
753 Bool_t selectOK = kFALSE;
755 if (kfl >= fFlavorSelect) flavorOK = kTRUE;
757 if (kfl > fFlavorSelect) {
761 if (kfl == fFlavorSelect) flavorOK = kTRUE;
763 switch (fStackFillOpt) {
765 case kFlavorSelection:
768 case kParentSelection:
769 if (ParentSelected(kf) || kf <= 10) selectOK = kTRUE;
772 if (flavorOK && selectOK) {
774 // Heavy flavor hadron or quark
776 // Kinematic seletion on final state heavy flavor mesons
777 if (ParentSelected(kf) && !KinematicSelection(iparticle, 0))
782 if (ParentSelected(kf)) ++nParents; // Update parent count
783 // printf("\n particle (HF) %d %d %d", i, pSelected[i], kf);
785 // Kinematic seletion on decay products
786 if (fCutOnChild && ParentSelected(kfMo) && ChildSelected(kf)
787 && !KinematicSelection(iparticle, 1))
793 // Select if mother was selected and is not tracked
795 if (pSelected[ipa] &&
796 !trackIt[ipa] && // mother will be tracked ?
797 kfMo != 5 && // mother is b-quark, don't store fragments
798 kfMo != 4 && // mother is c-quark, don't store fragments
799 kf != 92) // don't store string
802 // Semi-stable or de-selected: diselect decay products:
805 if (pSelected[i] == -1 || fDecayer->GetLifetime(kf) > fMaxLifeTime)
807 Int_t ipF = iparticle->GetFirstDaughter();
808 Int_t ipL = iparticle->GetLastDaughter();
809 if (ipF > 0) for (j = ipF-1; j < ipL; j++) pSelected[j] = -1;
811 // printf("\n particle (decay) %d %d %d", i, pSelected[i], kf);
812 pSelected[i] = (pSelected[i] == -1) ? 0 : 1;
815 if (pSelected[i] == -1) pSelected[i] = 0;
816 if (!pSelected[i]) continue;
817 // Count quarks only if you did not include fragmentation
818 if (fFragmentation && kf <= 10) continue;
821 // Decision on tracking
824 // Track final state particle
825 if (ks == 1) trackIt[i] = 1;
826 // Track semi-stable particles
827 if ((ks == 1) || (fDecayer->GetLifetime(kf) > fMaxLifeTime)) trackIt[i] = 1;
828 // Track particles selected by process if undecayed.
829 if (fForceDecay == kNoDecay) {
830 if (ParentSelected(kf)) trackIt[i] = 1;
832 if (ParentSelected(kf)) trackIt[i] = 0;
834 if (trackIt[i] == 1) ++nTkbles; // Update trackable counter
838 } // particle selection loop
840 for (i = 0; i<np; i++) {
841 if (!pSelected[i]) continue;
842 TParticle * iparticle = (TParticle *) fParticles.At(i);
843 kf = CheckPDGCode(iparticle->GetPdgCode());
844 Int_t ks = iparticle->GetStatusCode();
845 p[0] = iparticle->Px();
846 p[1] = iparticle->Py();
847 p[2] = iparticle->Pz();
848 p[3] = iparticle->Energy();
850 origin[0] = fVertex[0]+iparticle->Vx()/10; // [cm]
851 origin[1] = fVertex[1]+iparticle->Vy()/10; // [cm]
852 origin[2] = fVertex[2]+iparticle->Vz()/10; // [cm]
854 Float_t tof = fTime + kconv*iparticle->T();
855 Int_t ipa = iparticle->GetFirstMother()-1;
856 Int_t iparent = (ipa > -1) ? pParent[ipa] : -1;
858 PushTrack(fTrackIt*trackIt[i], iparent, kf,
859 p[0], p[1], p[2], p[3],
860 origin[0], origin[1], origin[2], tof,
861 polar[0], polar[1], polar[2],
862 kPPrimary, nt, 1., ks);
879 switch (fCountMode) {
881 // printf(" Count all \n");
885 // printf(" Count parents \n");
888 case kCountTrackables:
889 // printf(" Count trackable \n");
893 if (jev >= fNpart || fNpart == -1) {
894 fKineBias=Float_t(fNpart)/Float_t(fTrials);
896 fQ += fPythia->GetVINT(51);
897 fX1 += fPythia->GetVINT(41);
898 fX2 += fPythia->GetVINT(42);
899 fTrialsRun += fTrials;
906 SetHighWaterMark(nt);
907 // adjust weight due to kinematic selection
910 fXsection=fPythia->GetPARI(1);
913 Int_t AliGenPythia::GenerateMB()
916 // Min Bias selection and other global selections
918 Int_t i, kf, nt, iparent;
921 Double_t polar[3] = {0,0,0};
922 Double_t origin[3] = {0,0,0};
923 // converts from mm/c to s
924 const Float_t kconv=0.001/2.999792458e8;
928 Int_t np = (fHadronisation) ? fParticles.GetEntriesFast() : fNpartons;
932 Int_t* pParent = new Int_t[np];
933 for (i=0; i< np; i++) pParent[i] = -1;
934 if (fProcess == kPyJets || fProcess == kPyDirectGamma || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi) {
935 TParticle* jet1 = (TParticle *) fParticles.At(6);
936 TParticle* jet2 = (TParticle *) fParticles.At(7);
937 if (!CheckTrigger(jet1, jet2)) {
943 // Select jets with fragmentation photon or pi0 going to PHOS or EMCAL
944 if (fProcess == kPyJets && (fFragPhotonInCalo || fPi0InCalo) ) {
949 if (fFragPhotonInCalo) pdg = 22 ; // Photon
950 else if (fPi0InCalo) pdg = 111 ; // Pi0
952 for (i=0; i< np; i++) {
953 TParticle* iparticle = (TParticle *) fParticles.At(i);
954 if(iparticle->GetStatusCode()==1 && iparticle->GetPdgCode()==pdg &&
955 iparticle->Pt() > fFragPhotonOrPi0MinPt){
956 Int_t imother = iparticle->GetFirstMother() - 1;
957 TParticle* pmother = (TParticle *) fParticles.At(imother);
959 (pdg == 22 && pmother->GetStatusCode() != 11)) //No photon from hadron decay
961 Float_t phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees
962 Float_t eta =TMath::Abs(iparticle->Eta()); //in calos etamin=-etamax
963 if((fCheckEMCAL && IsInEMCAL(phi,eta)) ||
964 (fCheckPHOS && IsInPHOS(phi,eta)) )
975 // Select beauty jets with electron in EMCAL
976 if (fProcess == kPyBeautyJets && fEleInEMCAL) {
980 Int_t pdg = 11; //electron
985 for (i=0; i< np; i++) {
986 TParticle* iparticle = (TParticle *) fParticles.At(i);
987 if(iparticle->GetStatusCode()==1 && TMath::Abs(iparticle->GetPdgCode())==pdg &&
988 iparticle->Pt() > fElectronMinPt){
989 pt = iparticle->Pt();
990 phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees
991 eta =TMath::Abs(iparticle->Eta()); //in calos etamin=-etamax
992 if(IsInEMCAL(phi,eta))
1001 AliDebug(5,Form("Found an electron jet (pt,eta,phi) = (%f,%f,%f)",pt,eta,phi));
1003 // Check for diffraction
1005 if(fItune==320 && ( (TMath::Abs(fEnergyCMS - 900) < 1) || (TMath::Abs(fEnergyCMS - 2760) < 1) || (TMath::Abs(fEnergyCMS - 7000) < 1)) ) {
1006 if(!CheckDiffraction()) {
1013 // Check for minimum multiplicity
1014 if (fTriggerMultiplicity > 0) {
1015 Int_t multiplicity = 0;
1016 for (i = 0; i < np; i++) {
1017 TParticle * iparticle = (TParticle *) fParticles.At(i);
1019 Int_t statusCode = iparticle->GetStatusCode();
1021 // Initial state particle
1022 if (statusCode != 1)
1025 if (fTriggerMultiplicityEta > 0 && TMath::Abs(iparticle->Eta()) > fTriggerMultiplicityEta)
1028 if (iparticle->Pt() < fTriggerMultiplicityPtMin)
1031 TParticlePDG* pdgPart = iparticle->GetPDG();
1032 if (pdgPart && pdgPart->Charge() == 0)
1038 if (multiplicity < fTriggerMultiplicity) {
1042 Printf("Triggered on event with multiplicity of %d >= %d", multiplicity, fTriggerMultiplicity);
1045 // Select events with a photon pt > min pt going to PHOS eta acceptance or exactly PHOS eta phi
1046 if ((fProcess == kPyJets || fProcess == kPyDirectGamma) && fPhotonInCalo && (fCheckPHOSeta || fCheckPHOS)){
1048 Bool_t okd = kFALSE;
1052 for (i=0; i< np; i++) {
1053 TParticle* iparticle = (TParticle *) fParticles.At(i);
1054 Float_t phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees
1055 Float_t eta =TMath::Abs(iparticle->Eta());//in calos etamin=-etamax
1057 if(iparticle->GetStatusCode() == 1
1058 && iparticle->GetPdgCode() == pdg
1059 && iparticle->Pt() > fPhotonMinPt
1062 // first check if the photon is in PHOS phi
1063 if(IsInPHOS(phi,eta)){
1067 if(fCheckPHOSeta) iphcand = i; // candiate photon to rotate in phi
1072 if(!okd && iphcand != -1) // execute rotation in phi
1073 RotatePhi(iphcand,okd);
1081 if (fTriggerParticle) {
1082 Bool_t triggered = kFALSE;
1083 for (i = 0; i < np; i++) {
1084 TParticle * iparticle = (TParticle *) fParticles.At(i);
1085 kf = CheckPDGCode(iparticle->GetPdgCode());
1086 if (kf != fTriggerParticle) continue;
1087 if (iparticle->Pt() == 0.) continue;
1088 if (TMath::Abs(iparticle->Eta()) > fTriggerEta) continue;
1089 if ( iparticle->Pt() > fTriggerMaxPt || iparticle->Pt() < fTriggerMinPt ) continue;
1100 // Check if there is a ccbar or bbbar pair with at least one of the two
1101 // in fYMin < y < fYMax
1103 if (fProcess == kPyCharmppMNRwmi || fProcess == kPyBeautyppMNRwmi || fProcess == kPyBeautyJets) {
1104 TParticle *partCheck;
1106 Bool_t theQ=kFALSE,theQbar=kFALSE,inYcut=kFALSE;
1107 Bool_t theChild=kFALSE;
1108 Bool_t triggered=kFALSE;
1110 Int_t pdg,mpdg,mpdgUpperFamily;
1111 for(i=0; i<np; i++) {
1112 partCheck = (TParticle*)fParticles.At(i);
1113 pdg = partCheck->GetPdgCode();
1114 if(TMath::Abs(pdg) == fFlavorSelect) { // quark
1115 if(pdg>0) { theQ=kTRUE; } else { theQbar=kTRUE; }
1116 y = 0.5*TMath::Log((partCheck->Energy()+partCheck->Pz()+1.e-13)/
1117 (partCheck->Energy()-partCheck->Pz()+1.e-13));
1118 if(y>fYMin && y<fYMax) inYcut=kTRUE;
1120 if(fTriggerParticle) {
1121 if(TMath::Abs(pdg)==fTriggerParticle) triggered=kTRUE;
1123 if(fCutOnChild && TMath::Abs(pdg) == fPdgCodeParticleforAcceptanceCut) {
1124 Int_t mi = partCheck->GetFirstMother() - 1;
1126 mother = (TParticle*)fParticles.At(mi);
1127 mpdg=TMath::Abs(mother->GetPdgCode());
1128 mpdgUpperFamily=(mpdg>1000 ? mpdg+1000 : mpdg+100); // keep e from c from b
1129 if ( ParentSelected(mpdg) ||
1130 (fFlavorSelect==5 && ParentSelected(mpdgUpperFamily))) {
1131 if (KinematicSelection(partCheck,1)) {
1137 if (!theQ || !theQbar || !inYcut) { // one of the c/b conditions not satisfied
1141 if (fCutOnChild && !theChild) { // one of the child conditions not satisfied
1145 if(fTriggerParticle && !triggered) { // particle requested is not produced
1152 //Introducing child cuts in case kPyW, kPyZ, kPyMb, and kPyMbNonDiff
1153 if ( (fProcess == kPyW ||
1155 fProcess == kPyMbDefault ||
1156 fProcess == kPyMb ||
1157 fProcess == kPyMbAtlasTuneMC09 ||
1158 fProcess == kPyMbWithDirectPhoton ||
1159 fProcess == kPyMbNonDiffr)
1160 && (fCutOnChild == 1) ) {
1161 if ( !CheckKinematicsOnChild() ) {
1168 for (i = 0; i < np; i++) {
1170 TParticle * iparticle = (TParticle *) fParticles.At(i);
1171 kf = CheckPDGCode(iparticle->GetPdgCode());
1172 Int_t ks = iparticle->GetStatusCode();
1173 Int_t km = iparticle->GetFirstMother();
1175 (((ks == 1 && kf!=0 && KinematicSelection(iparticle, 0)) || (ks !=1)) && ((fStackFillOpt != kHeavyFlavor) || IsFromHeavyFlavor(i))) ||
1176 ((fProcess == kPyJets || fProcess == kPyBeautyJets) && ks == 21 && km == 0 && i>1)
1180 if (ks == 1) trackIt = 1;
1181 Int_t ipa = iparticle->GetFirstMother()-1;
1183 iparent = (ipa > -1) ? pParent[ipa] : -1;
1186 // store track information
1187 p[0] = iparticle->Px();
1188 p[1] = iparticle->Py();
1189 p[2] = iparticle->Pz();
1190 p[3] = iparticle->Energy();
1193 origin[0] = fVertex[0]+iparticle->Vx()/10; // [cm]
1194 origin[1] = fVertex[1]+iparticle->Vy()/10; // [cm]
1195 origin[2] = fVertex[2]+iparticle->Vz()/10; // [cm]
1197 Float_t tof = fTime + fEventTime + kconv * iparticle->T();
1199 PushTrack(fTrackIt*trackIt, iparent, kf,
1200 p[0], p[1], p[2], p[3],
1201 origin[0], origin[1], origin[2], tof,
1202 polar[0], polar[1], polar[2],
1203 kPPrimary, nt, 1., ks);
1207 SetHighWaterMark(nt);
1209 } // select particle
1218 void AliGenPythia::FinishRun()
1220 // Print x-section summary
1229 printf("\nTotal number of Pyevnt() calls %d\n", fTrialsRun);
1230 printf("\nMean Q, x1, x2: %f %f %f\n", fQ, fX1, fX2);
1233 void AliGenPythia::AdjustWeights() const
1235 // Adjust the weights after generation of all events
1239 Int_t ntrack=gAlice->GetMCApp()->GetNtrack();
1240 for (Int_t i=0; i<ntrack; i++) {
1241 part= gAlice->GetMCApp()->Particle(i);
1242 part->SetWeight(part->GetWeight()*fKineBias);
1247 void AliGenPythia::SetNuclei(Int_t a1, Int_t a2, Int_t pdfset)
1249 // Treat protons as inside nuclei with mass numbers a1 and a2
1253 fNucPdf = pdfset; // 0 EKS98 1 EPS08
1258 void AliGenPythia::MakeHeader()
1261 // Make header for the simulated event
1264 if (gAlice->GetEvNumber()>=fDebugEventFirst &&
1265 gAlice->GetEvNumber()<=fDebugEventLast) fPythia->Pylist(2);
1268 // Builds the event header, to be called after each event
1269 if (fHeader) delete fHeader;
1270 fHeader = new AliGenPythiaEventHeader("Pythia");
1274 // if(fProcDiff == 92 || fProcDiff == 93) printf("\n\n\n\n\n");
1275 // printf("fPythia->GetMSTI(1) = %d fProcDiff = %d\n",fPythia->GetMSTI(1), fProcDiff);
1276 ((AliGenPythiaEventHeader*) fHeader)->SetProcessType(fProcDiff);
1279 ((AliGenPythiaEventHeader*) fHeader)->SetProcessType(fPythia->GetMSTI(1));
1282 ((AliGenPythiaEventHeader*) fHeader)->SetTrials(fTrials);
1285 fHeader->SetPrimaryVertex(fVertex);
1286 fHeader->SetInteractionTime(fTime+fEventTime);
1288 // Number of primaries
1289 fHeader->SetNProduced(fNprimaries);
1291 // Jets that have triggered
1293 //Need to store jets for b-jet studies too!
1294 if (fProcess == kPyJets || fProcess == kPyDirectGamma || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi)
1297 Float_t jets[4][10];
1298 GetJets(njet, ntrig, jets);
1301 for (Int_t i = 0; i < ntrig; i++) {
1302 ((AliGenPythiaEventHeader*) fHeader)->AddJet(jets[0][i], jets[1][i], jets[2][i],
1307 // Copy relevant information from external header, if present.
1312 AliGenPythiaEventHeader* exHeader = (AliGenPythiaEventHeader*) (fRL->GetHeader()->GenEventHeader());
1313 for (Int_t i = 0; i < exHeader->NTriggerJets(); i++)
1315 printf("Adding Jet %d %d \n", i, exHeader->NTriggerJets());
1318 exHeader->TriggerJet(i, uqJet);
1319 ((AliGenPythiaEventHeader*) fHeader)->AddUQJet(uqJet[0], uqJet[1], uqJet[2], uqJet[3]);
1323 // Store quenching parameters
1326 Double_t z[4] = {0.};
1332 fPythia->GetQuenchingParameters(xp, yp, z);
1333 } else if (fQuench == 2){
1335 Double_t r1 = PARIMP.rb1;
1336 Double_t r2 = PARIMP.rb2;
1337 Double_t b = PARIMP.b1;
1338 Double_t r = 0.5 * TMath::Sqrt(2. * (r1 * r1 + r2 * r2) - b * b);
1339 Double_t phi = PARIMP.psib1;
1340 xp = r * TMath::Cos(phi);
1341 yp = r * TMath::Sin(phi);
1343 } else if (fQuench == 4) {
1347 AliFastGlauber::Instance()->GetSavedXY(xy);
1348 AliFastGlauber::Instance()->GetSavedI0I1(i0i1);
1351 ((AliGenPythiaEventHeader*) fHeader)->SetImpactParameter(fImpact);
1354 ((AliGenPythiaEventHeader*) fHeader)->SetXYJet(xp, yp);
1355 ((AliGenPythiaEventHeader*) fHeader)->SetZQuench(z);
1359 ((AliGenPythiaEventHeader*) fHeader)->SetPtHard(fPythia->GetVINT(47));
1367 Bool_t AliGenPythia::CheckTrigger(const TParticle* jet1, const TParticle* jet2)
1369 // Check the kinematic trigger condition
1372 eta[0] = jet1->Eta();
1373 eta[1] = jet2->Eta();
1375 phi[0] = jet1->Phi();
1376 phi[1] = jet2->Phi();
1378 pdg[0] = jet1->GetPdgCode();
1379 pdg[1] = jet2->GetPdgCode();
1380 Bool_t triggered = kFALSE;
1382 if (fProcess == kPyJets || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi) {
1385 Float_t jets[4][10];
1387 // Use Pythia clustering on parton level to determine jet axis
1389 GetJets(njets, ntrig, jets);
1391 if (ntrig || fEtMinJet == 0.) triggered = kTRUE;
1396 if (pdg[0] == kGamma) {
1400 //Check eta range first...
1401 if ((eta[ij] < fEtaMaxJet && eta[ij] > fEtaMinJet) &&
1402 (eta[ig] < fEtaMaxGamma && eta[ig] > fEtaMinGamma))
1404 //Eta is okay, now check phi range
1405 if ((phi[ij] < fPhiMaxJet && phi[ij] > fPhiMinJet) &&
1406 (phi[ig] < fPhiMaxGamma && phi[ig] > fPhiMinGamma))
1417 Bool_t AliGenPythia::CheckKinematicsOnChild(){
1419 //Checking Kinematics on Child (status code 1, particle code ?, kin cuts
1421 Bool_t checking = kFALSE;
1422 Int_t j, kcode, ks, km;
1423 Int_t nPartAcc = 0; //number of particles in the acceptance range
1424 Int_t numberOfAcceptedParticles = 1;
1425 if (fNumberOfAcceptedParticles != 0) { numberOfAcceptedParticles = fNumberOfAcceptedParticles; }
1426 Int_t npart = fParticles.GetEntriesFast();
1428 for (j = 0; j<npart; j++) {
1429 TParticle * jparticle = (TParticle *) fParticles.At(j);
1430 kcode = TMath::Abs( CheckPDGCode(jparticle->GetPdgCode()) );
1431 ks = jparticle->GetStatusCode();
1432 km = jparticle->GetFirstMother();
1434 if( (ks == 1) && (kcode == fPdgCodeParticleforAcceptanceCut) && (KinematicSelection(jparticle,1)) ){
1437 if( numberOfAcceptedParticles <= nPartAcc){
1446 void AliGenPythia::LoadEvent(AliStack* stack, Int_t flag, Int_t reHadr)
1449 // Load event into Pythia Common Block
1452 Int_t npart = stack -> GetNprimary();
1456 (fPythia->GetPyjets())->N = npart;
1458 n0 = (fPythia->GetPyjets())->N;
1459 (fPythia->GetPyjets())->N = n0 + npart;
1463 for (Int_t part = 0; part < npart; part++) {
1464 TParticle *mPart = stack->Particle(part);
1466 Int_t kf = mPart->GetPdgCode();
1467 Int_t ks = mPart->GetStatusCode();
1468 Int_t idf = mPart->GetFirstDaughter();
1469 Int_t idl = mPart->GetLastDaughter();
1472 if (ks == 11 || ks == 12) {
1479 Float_t px = mPart->Px();
1480 Float_t py = mPart->Py();
1481 Float_t pz = mPart->Pz();
1482 Float_t e = mPart->Energy();
1483 Float_t m = mPart->GetCalcMass();
1486 (fPythia->GetPyjets())->P[0][part+n0] = px;
1487 (fPythia->GetPyjets())->P[1][part+n0] = py;
1488 (fPythia->GetPyjets())->P[2][part+n0] = pz;
1489 (fPythia->GetPyjets())->P[3][part+n0] = e;
1490 (fPythia->GetPyjets())->P[4][part+n0] = m;
1492 (fPythia->GetPyjets())->K[1][part+n0] = kf;
1493 (fPythia->GetPyjets())->K[0][part+n0] = ks;
1494 (fPythia->GetPyjets())->K[3][part+n0] = idf + 1;
1495 (fPythia->GetPyjets())->K[4][part+n0] = idl + 1;
1496 (fPythia->GetPyjets())->K[2][part+n0] = mPart->GetFirstMother() + 1;
1500 void AliGenPythia::LoadEvent(const TObjArray* stack, Int_t flag, Int_t reHadr)
1503 // Load event into Pythia Common Block
1506 Int_t npart = stack -> GetEntries();
1510 (fPythia->GetPyjets())->N = npart;
1512 n0 = (fPythia->GetPyjets())->N;
1513 (fPythia->GetPyjets())->N = n0 + npart;
1517 for (Int_t part = 0; part < npart; part++) {
1518 TParticle *mPart = dynamic_cast<TParticle *>(stack->At(part));
1519 if (!mPart) continue;
1521 Int_t kf = mPart->GetPdgCode();
1522 Int_t ks = mPart->GetStatusCode();
1523 Int_t idf = mPart->GetFirstDaughter();
1524 Int_t idl = mPart->GetLastDaughter();
1527 if (ks == 11 || ks == 12) {
1534 Float_t px = mPart->Px();
1535 Float_t py = mPart->Py();
1536 Float_t pz = mPart->Pz();
1537 Float_t e = mPart->Energy();
1538 Float_t m = mPart->GetCalcMass();
1541 (fPythia->GetPyjets())->P[0][part+n0] = px;
1542 (fPythia->GetPyjets())->P[1][part+n0] = py;
1543 (fPythia->GetPyjets())->P[2][part+n0] = pz;
1544 (fPythia->GetPyjets())->P[3][part+n0] = e;
1545 (fPythia->GetPyjets())->P[4][part+n0] = m;
1547 (fPythia->GetPyjets())->K[1][part+n0] = kf;
1548 (fPythia->GetPyjets())->K[0][part+n0] = ks;
1549 (fPythia->GetPyjets())->K[3][part+n0] = idf + 1;
1550 (fPythia->GetPyjets())->K[4][part+n0] = idl + 1;
1551 (fPythia->GetPyjets())->K[2][part+n0] = mPart->GetFirstMother() + 1;
1556 void AliGenPythia::RecJetsUA1(Int_t& njets, Float_t jets [4][50])
1559 // Calls the Pythia jet finding algorithm to find jets in the current event
1564 Int_t n = fPythia->GetN();
1568 fPythia->Pycell(njets);
1570 for (i = 0; i < njets; i++) {
1571 Float_t px = (fPythia->GetPyjets())->P[0][n+i];
1572 Float_t py = (fPythia->GetPyjets())->P[1][n+i];
1573 Float_t pz = (fPythia->GetPyjets())->P[2][n+i];
1574 Float_t e = (fPythia->GetPyjets())->P[3][n+i];
1585 void AliGenPythia::GetJets(Int_t& nJets, Int_t& nJetsTrig, Float_t jets[4][10])
1588 // Calls the Pythia clustering algorithm to find jets in the current event
1590 Int_t n = fPythia->GetN();
1593 if (fJetReconstruction == kCluster) {
1595 // Configure cluster algorithm
1597 fPythia->SetPARU(43, 2.);
1598 fPythia->SetMSTU(41, 1);
1600 // Call cluster algorithm
1602 fPythia->Pyclus(nJets);
1604 // Loading jets from common block
1610 fPythia->Pycell(nJets);
1614 for (i = 0; i < nJets; i++) {
1615 Float_t px = (fPythia->GetPyjets())->P[0][n+i];
1616 Float_t py = (fPythia->GetPyjets())->P[1][n+i];
1617 Float_t pz = (fPythia->GetPyjets())->P[2][n+i];
1618 Float_t e = (fPythia->GetPyjets())->P[3][n+i];
1619 Float_t pt = TMath::Sqrt(px * px + py * py);
1620 Float_t phi = TMath::Pi() + TMath::ATan2(-py, -px);
1621 Float_t theta = TMath::ATan2(pt,pz);
1622 Float_t et = e * TMath::Sin(theta);
1623 Float_t eta = -TMath::Log(TMath::Tan(theta / 2.));
1625 eta > fEtaMinJet && eta < fEtaMaxJet &&
1626 phi > fPhiMinJet && phi < fPhiMaxJet &&
1627 et > fEtMinJet && et < fEtMaxJet
1630 jets[0][nJetsTrig] = px;
1631 jets[1][nJetsTrig] = py;
1632 jets[2][nJetsTrig] = pz;
1633 jets[3][nJetsTrig] = e;
1635 // printf("\n........-Jet #%d: %10.3f %10.3f %10.3f %10.3f \n", i, pt, et, eta, phi * kRaddeg);
1637 // printf("\n........-Jet #%d: %10.3f %10.3f %10.3f %10.3f \n", i, pt, et, eta, phi * kRaddeg);
1642 void AliGenPythia::GetSubEventTime()
1644 // Calculates time of the next subevent
1647 TArrayF &array = *fEventsTime;
1648 fEventTime = array[fCurSubEvent++];
1650 // printf(" Event time: %d %f %p",fCurSubEvent,fEventTime,fEventsTime);
1654 Bool_t AliGenPythia::IsInEMCAL(Float_t phi, Float_t eta) const
1656 // Is particle in EMCAL acceptance?
1657 // phi in degrees, etamin=-etamax
1658 if(phi > fEMCALMinPhi && phi < fEMCALMaxPhi &&
1665 Bool_t AliGenPythia::IsInPHOS(Float_t phi, Float_t eta) const
1667 // Is particle in PHOS acceptance?
1668 // Acceptance slightly larger considered.
1669 // phi in degrees, etamin=-etamax
1670 if(phi > fPHOSMinPhi && phi < fPHOSMaxPhi &&
1677 void AliGenPythia::RotatePhi(Int_t iphcand, Bool_t& okdd)
1679 //calculate the new position random between fPHOSMinPhi and fPHOSMaxPhi
1680 Double_t phiPHOSmin = TMath::Pi()*fPHOSMinPhi/180;
1681 Double_t phiPHOSmax = TMath::Pi()*fPHOSMaxPhi/180;
1682 Double_t phiPHOS = (AliPythiaRndm::GetPythiaRandom())->Uniform(phiPHOSmin,phiPHOSmax);
1684 //calculate deltaphi
1685 TParticle* ph = (TParticle *) fParticles.At(iphcand);
1686 Double_t phphi = ph->Phi();
1687 Double_t deltaphi = phiPHOS - phphi;
1691 //loop for all particles and produce the phi rotation
1692 Int_t np = (fHadronisation) ? fParticles.GetEntriesFast() : fNpartons;
1693 Double_t oldphi, newphi;
1694 Double_t newVx, newVy, r, vZ, time;
1695 Double_t newPx, newPy, pt, pz, e;
1696 for(Int_t i=0; i< np; i++) {
1697 TParticle* iparticle = (TParticle *) fParticles.At(i);
1698 oldphi = iparticle->Phi();
1699 newphi = oldphi + deltaphi;
1700 if(newphi < 0) newphi = 2*TMath::Pi() + newphi; // correct angle
1701 if(newphi > 2*TMath::Pi()) newphi = newphi - 2*TMath::Pi(); // correct angle
1704 newVx = r * TMath::Cos(newphi);
1705 newVy = r * TMath::Sin(newphi);
1706 vZ = iparticle->Vz(); // don't transform
1707 time = iparticle->T(); // don't transform
1709 pt = iparticle->Pt();
1710 newPx = pt * TMath::Cos(newphi);
1711 newPy = pt * TMath::Sin(newphi);
1712 pz = iparticle->Pz(); // don't transform
1713 e = iparticle->Energy(); // don't transform
1716 iparticle->SetProductionVertex(newVx, newVy, vZ, time);
1717 iparticle->SetMomentum(newPx, newPy, pz, e);
1719 } //end particle loop
1721 // now let's check that we put correctly the candidate photon in PHOS
1722 Float_t phi = ph->Phi()*180./TMath::Pi(); //Convert to degrees
1723 Float_t eta =TMath::Abs(ph->Eta());//in calos etamin=-etamax
1724 if(IsInPHOS(phi,eta))
1731 Bool_t AliGenPythia::CheckDiffraction()
1733 // use this method only with Perugia-0 tune!
1737 Int_t np = (fHadronisation) ? fParticles.GetEntriesFast() : fNpartons;
1743 Double_t y2 = -1e10;
1745 const Int_t kNstable=20;
1746 const Int_t pdgStable[20] = {
1749 12, // Electron Neutrino
1751 14, // Muon Neutrino
1762 3112, // Sigma Minus
1769 for (Int_t i = 0; i < np; i++) {
1770 TParticle * part = (TParticle *) fParticles.At(i);
1772 Int_t statusCode = part->GetStatusCode();
1774 // Initial state particle
1775 if (statusCode != 1)
1778 Int_t pdg = TMath::Abs(part->GetPdgCode());
1779 Bool_t isStable = kFALSE;
1780 for (Int_t i1 = 0; i1 < kNstable; i1++) {
1781 if (pdg == pdgStable[i1]) {
1789 Double_t y = part->Y();
1803 if(iPart1<0 || iPart2<0) return kFALSE;
1808 TParticle * part1 = (TParticle *) fParticles.At(iPart1);
1809 TParticle * part2 = (TParticle *) fParticles.At(iPart2);
1811 Int_t pdg1 = part1->GetPdgCode();
1812 Int_t pdg2 = part2->GetPdgCode();
1816 if (pdg1 == 2212 && pdg2 == 2212)
1824 if((AliPythiaRndm::GetPythiaRandom())->Uniform(0.,1.)>0.5) iPart = iPart2;
1827 else if (pdg1 == 2212)
1829 else if (pdg2 == 2212)
1838 TParticle * part = (TParticle *) fParticles.At(iPart);
1839 Double_t E= part->Energy();
1840 Double_t P= part->P();
1841 M= TMath::Sqrt((fEnergyCMS-E-P)*(fEnergyCMS-E+P));
1844 Double_t Mmin, Mmax, wSD, wDD, wND;
1845 if(!GetWeightsDiffraction(M, Mmin, Mmax, wSD, wDD, wND)) return kFALSE;
1847 if(M>-1 && M<Mmin) return kFALSE;
1850 Int_t procType=fPythia->GetMSTI(1);
1852 if(procType== 94) proc0=1;
1853 if(procType== 92 || procType== 93) proc0=0;
1857 else if(proc0==1) proc=1;
1859 if(proc==1 && (AliPythiaRndm::GetPythiaRandom())->Uniform(0.,1.) > wDD) return kFALSE;
1860 if(proc==2 && (AliPythiaRndm::GetPythiaRandom())->Uniform(0.,1.) > wND) return kFALSE;
1863 // if(proc==1 || proc==2) return kFALSE;
1866 if(proc0!=0) fProcDiff = procType;
1867 else fProcDiff = 95;
1871 if(wSD<0) AliError("wSD<0 ! \n");
1873 if((AliPythiaRndm::GetPythiaRandom())->Uniform(0.,1.)> wSD) return kFALSE;
1875 // printf("iPart = %d\n", iPart);
1877 if(iPart==iPart1) fProcDiff=93;
1878 else if(iPart==iPart2) fProcDiff=92;
1880 printf("EROOR: iPart!=iPart1 && iPart!=iPart2\n");
1889 Bool_t AliGenPythia::GetWeightsDiffraction(Double_t M, Double_t &Mmin, Double_t &Mmax,
1890 Double_t &wSD, Double_t &wDD, Double_t &wND)
1894 if(TMath::Abs(fEnergyCMS-900)<1 ){
1896 const Int_t nbin=400;
1898 1.080000, 1.577300, 2.074600, 2.571900, 3.069200, 3.566500,
1899 4.063800, 4.561100, 5.058400, 5.555700, 6.053000, 6.550300,
1900 7.047600, 7.544900, 8.042200, 8.539500, 9.036800, 9.534100,
1901 10.031400, 10.528700, 11.026000, 11.523300, 12.020600, 12.517900,
1902 13.015200, 13.512500, 14.009800, 14.507100, 15.004400, 15.501700,
1903 15.999000, 16.496300, 16.993600, 17.490900, 17.988200, 18.485500,
1904 18.982800, 19.480100, 19.977400, 20.474700, 20.972000, 21.469300,
1905 21.966600, 22.463900, 22.961200, 23.458500, 23.955800, 24.453100,
1906 24.950400, 25.447700, 25.945000, 26.442300, 26.939600, 27.436900,
1907 27.934200, 28.431500, 28.928800, 29.426100, 29.923400, 30.420700,
1908 30.918000, 31.415300, 31.912600, 32.409900, 32.907200, 33.404500,
1909 33.901800, 34.399100, 34.896400, 35.393700, 35.891000, 36.388300,
1910 36.885600, 37.382900, 37.880200, 38.377500, 38.874800, 39.372100,
1911 39.869400, 40.366700, 40.864000, 41.361300, 41.858600, 42.355900,
1912 42.853200, 43.350500, 43.847800, 44.345100, 44.842400, 45.339700,
1913 45.837000, 46.334300, 46.831600, 47.328900, 47.826200, 48.323500,
1914 48.820800, 49.318100, 49.815400, 50.312700, 50.810000, 51.307300,
1915 51.804600, 52.301900, 52.799200, 53.296500, 53.793800, 54.291100,
1916 54.788400, 55.285700, 55.783000, 56.280300, 56.777600, 57.274900,
1917 57.772200, 58.269500, 58.766800, 59.264100, 59.761400, 60.258700,
1918 60.756000, 61.253300, 61.750600, 62.247900, 62.745200, 63.242500,
1919 63.739800, 64.237100, 64.734400, 65.231700, 65.729000, 66.226300,
1920 66.723600, 67.220900, 67.718200, 68.215500, 68.712800, 69.210100,
1921 69.707400, 70.204700, 70.702000, 71.199300, 71.696600, 72.193900,
1922 72.691200, 73.188500, 73.685800, 74.183100, 74.680400, 75.177700,
1923 75.675000, 76.172300, 76.669600, 77.166900, 77.664200, 78.161500,
1924 78.658800, 79.156100, 79.653400, 80.150700, 80.648000, 81.145300,
1925 81.642600, 82.139900, 82.637200, 83.134500, 83.631800, 84.129100,
1926 84.626400, 85.123700, 85.621000, 86.118300, 86.615600, 87.112900,
1927 87.610200, 88.107500, 88.604800, 89.102100, 89.599400, 90.096700,
1928 90.594000, 91.091300, 91.588600, 92.085900, 92.583200, 93.080500,
1929 93.577800, 94.075100, 94.572400, 95.069700, 95.567000, 96.064300,
1930 96.561600, 97.058900, 97.556200, 98.053500, 98.550800, 99.048100,
1931 99.545400, 100.042700, 100.540000, 101.037300, 101.534600, 102.031900,
1932 102.529200, 103.026500, 103.523800, 104.021100, 104.518400, 105.015700,
1933 105.513000, 106.010300, 106.507600, 107.004900, 107.502200, 107.999500,
1934 108.496800, 108.994100, 109.491400, 109.988700, 110.486000, 110.983300,
1935 111.480600, 111.977900, 112.475200, 112.972500, 113.469800, 113.967100,
1936 114.464400, 114.961700, 115.459000, 115.956300, 116.453600, 116.950900,
1937 117.448200, 117.945500, 118.442800, 118.940100, 119.437400, 119.934700,
1938 120.432000, 120.929300, 121.426600, 121.923900, 122.421200, 122.918500,
1939 123.415800, 123.913100, 124.410400, 124.907700, 125.405000, 125.902300,
1940 126.399600, 126.896900, 127.394200, 127.891500, 128.388800, 128.886100,
1941 129.383400, 129.880700, 130.378000, 130.875300, 131.372600, 131.869900,
1942 132.367200, 132.864500, 133.361800, 133.859100, 134.356400, 134.853700,
1943 135.351000, 135.848300, 136.345600, 136.842900, 137.340200, 137.837500,
1944 138.334800, 138.832100, 139.329400, 139.826700, 140.324000, 140.821300,
1945 141.318600, 141.815900, 142.313200, 142.810500, 143.307800, 143.805100,
1946 144.302400, 144.799700, 145.297000, 145.794300, 146.291600, 146.788900,
1947 147.286200, 147.783500, 148.280800, 148.778100, 149.275400, 149.772700,
1948 150.270000, 150.767300, 151.264600, 151.761900, 152.259200, 152.756500,
1949 153.253800, 153.751100, 154.248400, 154.745700, 155.243000, 155.740300,
1950 156.237600, 156.734900, 157.232200, 157.729500, 158.226800, 158.724100,
1951 159.221400, 159.718700, 160.216000, 160.713300, 161.210600, 161.707900,
1952 162.205200, 162.702500, 163.199800, 163.697100, 164.194400, 164.691700,
1953 165.189000, 165.686300, 166.183600, 166.680900, 167.178200, 167.675500,
1954 168.172800, 168.670100, 169.167400, 169.664700, 170.162000, 170.659300,
1955 171.156600, 171.653900, 172.151200, 172.648500, 173.145800, 173.643100,
1956 174.140400, 174.637700, 175.135000, 175.632300, 176.129600, 176.626900,
1957 177.124200, 177.621500, 178.118800, 178.616100, 179.113400, 179.610700,
1958 180.108000, 180.605300, 181.102600, 181.599900, 182.097200, 182.594500,
1959 183.091800, 183.589100, 184.086400, 184.583700, 185.081000, 185.578300,
1960 186.075600, 186.572900, 187.070200, 187.567500, 188.064800, 188.562100,
1961 189.059400, 189.556700, 190.054000, 190.551300, 191.048600, 191.545900,
1962 192.043200, 192.540500, 193.037800, 193.535100, 194.032400, 194.529700,
1963 195.027000, 195.524300, 196.021600, 196.518900, 197.016200, 197.513500,
1964 198.010800, 198.508100, 199.005400, 199.502700, 200.000000};
1966 1.000000, 0.643457, 0.645609, 0.648347, 0.604563, 0.605002,
1967 0.602819, 0.611473, 0.576412, 0.562354, 0.550216, 0.529285,
1968 0.534558, 0.534364, 0.530358, 0.518475, 0.489253, 0.469754,
1969 0.469825, 0.450513, 0.455849, 0.435312, 0.437210, 0.456686,
1970 0.413577, 0.427093, 0.426894, 0.418834, 0.409475, 0.388483,
1971 0.412208, 0.388912, 0.389611, 0.382572, 0.389220, 0.370964,
1972 0.380463, 0.370873, 0.363701, 0.369363, 0.357361, 0.365759,
1973 0.348566, 0.337062, 0.348190, 0.332330, 0.359001, 0.335836,
1974 0.339154, 0.335599, 0.336035, 0.335204, 0.353440, 0.337836,
1975 0.333874, 0.307120, 0.294963, 0.324978, 0.313359, 0.317688,
1976 0.323758, 0.319304, 0.335317, 0.301765, 0.317257, 0.356331,
1977 0.323119, 0.297732, 0.303188, 0.326102, 0.316467, 0.294728,
1978 0.308135, 0.288857, 0.325692, 0.312493, 0.291100, 0.325921,
1979 0.313317, 0.295980, 0.308481, 0.328380, 0.313081, 0.296763,
1980 0.295431, 0.317325, 0.320462, 0.286918, 0.316035, 0.335208,
1981 0.283193, 0.333945, 0.292534, 0.294164, 0.330680, 0.296992,
1982 0.285509, 0.317260, 0.298269, 0.311299, 0.312129, 0.286822,
1983 0.287442, 0.319139, 0.283314, 0.318454, 0.297727, 0.301597,
1984 0.282483, 0.294792, 0.305569, 0.290957, 0.297817, 0.282908,
1985 0.272401, 0.305584, 0.300220, 0.297020, 0.298781, 0.278008,
1986 0.277727, 0.323777, 0.287419, 0.342074, 0.287259, 0.303658,
1987 0.302668, 0.279622, 0.280586, 0.313630, 0.276068, 0.257051,
1988 0.309996, 0.265534, 0.297138, 0.281738, 0.294610, 0.292882,
1989 0.286860, 0.312686, 0.293244, 0.293744, 0.271375, 0.278734,
1990 0.280308, 0.304739, 0.287907, 0.285261, 0.311180, 0.313476,
1991 0.289660, 0.289282, 0.319505, 0.271285, 0.272008, 0.289245,
1992 0.281038, 0.285284, 0.295836, 0.281416, 0.283501, 0.295417,
1993 0.304372, 0.297764, 0.291378, 0.321530, 0.315604, 0.329507,
1994 0.282609, 0.275576, 0.283721, 0.311714, 0.283760, 0.273188,
1995 0.312193, 0.264347, 0.281532, 0.301226, 0.281718, 0.336408,
1996 0.283157, 0.332010, 0.289974, 0.290256, 0.301569, 0.332228,
1997 0.288282, 0.326339, 0.313653, 0.300361, 0.289470, 0.264830,
1998 0.298659, 0.272359, 0.278878, 0.306001, 0.328168, 0.294991,
1999 0.327737, 0.278056, 0.302435, 0.284183, 0.279270, 0.307279,
2000 0.307917, 0.315196, 0.283803, 0.313333, 0.315730, 0.304818,
2001 0.307171, 0.295223, 0.333741, 0.346911, 0.310143, 0.336686,
2002 0.275459, 0.334781, 0.295405, 0.275816, 0.301678, 0.327242,
2003 0.320717, 0.309230, 0.292145, 0.294489, 0.305088, 0.300969,
2004 0.277438, 0.326159, 0.297065, 0.301177, 0.303843, 0.275382,
2005 0.304019, 0.284166, 0.289610, 0.331611, 0.317131, 0.310880,
2006 0.360456, 0.294052, 0.342694, 0.327166, 0.336797, 0.298040,
2007 0.295767, 0.260053, 0.325544, 0.335310, 0.320182, 0.301072,
2008 0.313117, 0.283407, 0.299206, 0.293525, 0.305067, 0.255978,
2009 0.327055, 0.316382, 0.317700, 0.278993, 0.283120, 0.314000,
2010 0.274396, 0.291208, 0.348813, 0.319603, 0.313076, 0.289155,
2011 0.343988, 0.311426, 0.322896, 0.328726, 0.337062, 0.389093,
2012 0.284122, 0.312184, 0.304008, 0.319170, 0.320778, 0.288390,
2013 0.337272, 0.356273, 0.343310, 0.312209, 0.330709, 0.297977,
2014 0.346146, 0.369162, 0.324385, 0.339831, 0.337037, 0.318739,
2015 0.343157, 0.277720, 0.368407, 0.321330, 0.338997, 0.314220,
2016 0.328861, 0.321824, 0.328013, 0.356925, 0.359144, 0.296314,
2017 0.345415, 0.396711, 0.347032, 0.294928, 0.343799, 0.322331,
2018 0.328656, 0.326098, 0.337338, 0.337038, 0.300179, 0.351391,
2019 0.324337, 0.330896, 0.302842, 0.310522, 0.337052, 0.359989,
2020 0.383250, 0.359355, 0.315382, 0.333113, 0.342598, 0.355348,
2021 0.320751, 0.320475, 0.351762, 0.351475, 0.338358, 0.326153,
2022 0.302507, 0.340048, 0.318685, 0.381646, 0.339320, 0.299453,
2023 0.426599, 0.393515, 0.353929, 0.328435, 0.413976, 0.292558,
2024 0.379340, 0.358344, 0.409259, 0.313821, 0.336675, 0.324521,
2025 0.408382, 0.346273, 0.312939, 0.362453, 0.343152, 0.330577,
2026 0.332831, 0.353299, 0.347745, 0.334818, 0.332234, 0.385585,
2027 0.395483, 0.395316, 0.326972, 0.349434, 0.317604, 0.328980,
2028 0.375056, 0.317290, 0.357083, 0.346165, 0.310444, 0.356873,
2029 0.359523, 0.365308, 0.365122, 0.383685, 0.370975, 0.396928,
2030 0.407654, 0.307755, 0.323033, 0.350580, 0.345231, 0.342462,
2031 0.400000, 0.318309, 0.403570, 0.322856, 0.383053, 0.422252,
2032 0.386112, 0.364314, 0.434375, 0.334629};
2039 if(M<Mmin || M>Mmax) return kTRUE;
2042 for(Int_t i=1; i<=nbin; i++)
2045 // printf("Mi> %f && Mi< %f\n", bin[i-1], bin[i]);
2051 else if(TMath::Abs(fEnergyCMS-2760)<1 ){
2053 const Int_t nbin=400;
2055 1.080000, 1.577300, 2.074600, 2.571900, 3.069200, 3.566500,
2056 4.063800, 4.561100, 5.058400, 5.555700, 6.053000, 6.550300,
2057 7.047600, 7.544900, 8.042200, 8.539500, 9.036800, 9.534100,
2058 10.031400, 10.528700, 11.026000, 11.523300, 12.020600, 12.517900,
2059 13.015200, 13.512500, 14.009800, 14.507100, 15.004400, 15.501700,
2060 15.999000, 16.496300, 16.993600, 17.490900, 17.988200, 18.485500,
2061 18.982800, 19.480100, 19.977400, 20.474700, 20.972000, 21.469300,
2062 21.966600, 22.463900, 22.961200, 23.458500, 23.955800, 24.453100,
2063 24.950400, 25.447700, 25.945000, 26.442300, 26.939600, 27.436900,
2064 27.934200, 28.431500, 28.928800, 29.426100, 29.923400, 30.420700,
2065 30.918000, 31.415300, 31.912600, 32.409900, 32.907200, 33.404500,
2066 33.901800, 34.399100, 34.896400, 35.393700, 35.891000, 36.388300,
2067 36.885600, 37.382900, 37.880200, 38.377500, 38.874800, 39.372100,
2068 39.869400, 40.366700, 40.864000, 41.361300, 41.858600, 42.355900,
2069 42.853200, 43.350500, 43.847800, 44.345100, 44.842400, 45.339700,
2070 45.837000, 46.334300, 46.831600, 47.328900, 47.826200, 48.323500,
2071 48.820800, 49.318100, 49.815400, 50.312700, 50.810000, 51.307300,
2072 51.804600, 52.301900, 52.799200, 53.296500, 53.793800, 54.291100,
2073 54.788400, 55.285700, 55.783000, 56.280300, 56.777600, 57.274900,
2074 57.772200, 58.269500, 58.766800, 59.264100, 59.761400, 60.258700,
2075 60.756000, 61.253300, 61.750600, 62.247900, 62.745200, 63.242500,
2076 63.739800, 64.237100, 64.734400, 65.231700, 65.729000, 66.226300,
2077 66.723600, 67.220900, 67.718200, 68.215500, 68.712800, 69.210100,
2078 69.707400, 70.204700, 70.702000, 71.199300, 71.696600, 72.193900,
2079 72.691200, 73.188500, 73.685800, 74.183100, 74.680400, 75.177700,
2080 75.675000, 76.172300, 76.669600, 77.166900, 77.664200, 78.161500,
2081 78.658800, 79.156100, 79.653400, 80.150700, 80.648000, 81.145300,
2082 81.642600, 82.139900, 82.637200, 83.134500, 83.631800, 84.129100,
2083 84.626400, 85.123700, 85.621000, 86.118300, 86.615600, 87.112900,
2084 87.610200, 88.107500, 88.604800, 89.102100, 89.599400, 90.096700,
2085 90.594000, 91.091300, 91.588600, 92.085900, 92.583200, 93.080500,
2086 93.577800, 94.075100, 94.572400, 95.069700, 95.567000, 96.064300,
2087 96.561600, 97.058900, 97.556200, 98.053500, 98.550800, 99.048100,
2088 99.545400, 100.042700, 100.540000, 101.037300, 101.534600, 102.031900,
2089 102.529200, 103.026500, 103.523800, 104.021100, 104.518400, 105.015700,
2090 105.513000, 106.010300, 106.507600, 107.004900, 107.502200, 107.999500,
2091 108.496800, 108.994100, 109.491400, 109.988700, 110.486000, 110.983300,
2092 111.480600, 111.977900, 112.475200, 112.972500, 113.469800, 113.967100,
2093 114.464400, 114.961700, 115.459000, 115.956300, 116.453600, 116.950900,
2094 117.448200, 117.945500, 118.442800, 118.940100, 119.437400, 119.934700,
2095 120.432000, 120.929300, 121.426600, 121.923900, 122.421200, 122.918500,
2096 123.415800, 123.913100, 124.410400, 124.907700, 125.405000, 125.902300,
2097 126.399600, 126.896900, 127.394200, 127.891500, 128.388800, 128.886100,
2098 129.383400, 129.880700, 130.378000, 130.875300, 131.372600, 131.869900,
2099 132.367200, 132.864500, 133.361800, 133.859100, 134.356400, 134.853700,
2100 135.351000, 135.848300, 136.345600, 136.842900, 137.340200, 137.837500,
2101 138.334800, 138.832100, 139.329400, 139.826700, 140.324000, 140.821300,
2102 141.318600, 141.815900, 142.313200, 142.810500, 143.307800, 143.805100,
2103 144.302400, 144.799700, 145.297000, 145.794300, 146.291600, 146.788900,
2104 147.286200, 147.783500, 148.280800, 148.778100, 149.275400, 149.772700,
2105 150.270000, 150.767300, 151.264600, 151.761900, 152.259200, 152.756500,
2106 153.253800, 153.751100, 154.248400, 154.745700, 155.243000, 155.740300,
2107 156.237600, 156.734900, 157.232200, 157.729500, 158.226800, 158.724100,
2108 159.221400, 159.718700, 160.216000, 160.713300, 161.210600, 161.707900,
2109 162.205200, 162.702500, 163.199800, 163.697100, 164.194400, 164.691700,
2110 165.189000, 165.686300, 166.183600, 166.680900, 167.178200, 167.675500,
2111 168.172800, 168.670100, 169.167400, 169.664700, 170.162000, 170.659300,
2112 171.156600, 171.653900, 172.151200, 172.648500, 173.145800, 173.643100,
2113 174.140400, 174.637700, 175.135000, 175.632300, 176.129600, 176.626900,
2114 177.124200, 177.621500, 178.118800, 178.616100, 179.113400, 179.610700,
2115 180.108000, 180.605300, 181.102600, 181.599900, 182.097200, 182.594500,
2116 183.091800, 183.589100, 184.086400, 184.583700, 185.081000, 185.578300,
2117 186.075600, 186.572900, 187.070200, 187.567500, 188.064800, 188.562100,
2118 189.059400, 189.556700, 190.054000, 190.551300, 191.048600, 191.545900,
2119 192.043200, 192.540500, 193.037800, 193.535100, 194.032400, 194.529700,
2120 195.027000, 195.524300, 196.021600, 196.518900, 197.016200, 197.513500,
2121 198.010800, 198.508100, 199.005400, 199.502700, 200.000000};
2123 1.000000, 0.692593, 0.673384, 0.666273, 0.657285, 0.637723,
2124 0.625881, 0.643590, 0.606100, 0.589007, 0.567824, 0.578705,
2125 0.538530, 0.517937, 0.528278, 0.515918, 0.539461, 0.466186,
2126 0.489869, 0.468402, 0.465017, 0.453336, 0.460769, 0.474638,
2127 0.456347, 0.434471, 0.427478, 0.435435, 0.410934, 0.366431,
2128 0.382240, 0.379513, 0.394249, 0.386837, 0.353103, 0.382138,
2129 0.377497, 0.389479, 0.378736, 0.347933, 0.354605, 0.352077,
2130 0.324443, 0.358792, 0.339968, 0.359052, 0.330734, 0.318291,
2131 0.333703, 0.358644, 0.335819, 0.332213, 0.309051, 0.309975,
2132 0.331626, 0.304407, 0.309819, 0.312097, 0.312462, 0.320411,
2133 0.280401, 0.302311, 0.315863, 0.281479, 0.310003, 0.296911,
2134 0.313676, 0.281071, 0.294163, 0.306500, 0.283462, 0.274867,
2135 0.307149, 0.270555, 0.282264, 0.287373, 0.307849, 0.278675,
2136 0.286990, 0.278269, 0.300105, 0.286799, 0.265674, 0.275140,
2137 0.285702, 0.257352, 0.267714, 0.248204, 0.252220, 0.255678,
2138 0.282946, 0.268257, 0.282375, 0.262675, 0.275564, 0.248345,
2139 0.236259, 0.291914, 0.259936, 0.241338, 0.267389, 0.285044,
2140 0.289419, 0.253594, 0.284568, 0.231840, 0.260008, 0.268527,
2141 0.275363, 0.224115, 0.281260, 0.257913, 0.295152, 0.264399,
2142 0.232287, 0.282533, 0.223431, 0.255756, 0.244471, 0.221695,
2143 0.272450, 0.284244, 0.253682, 0.270717, 0.275403, 0.240323,
2144 0.245081, 0.241859, 0.216340, 0.244789, 0.220291, 0.238478,
2145 0.224691, 0.244058, 0.266117, 0.271478, 0.242012, 0.267321,
2146 0.248494, 0.253343, 0.255606, 0.235458, 0.241079, 0.233223,
2147 0.226813, 0.259224, 0.234239, 0.258606, 0.210892, 0.238186,
2148 0.243271, 0.222678, 0.213437, 0.273939, 0.247966, 0.232548,
2149 0.263438, 0.222089, 0.272111, 0.248818, 0.244295, 0.238368,
2150 0.236908, 0.248776, 0.232604, 0.231194, 0.227117, 0.231152,
2151 0.282140, 0.229778, 0.232631, 0.261794, 0.216633, 0.253471,
2152 0.242157, 0.227406, 0.269335, 0.230547, 0.210618, 0.285872,
2153 0.248776, 0.229875, 0.242728, 0.227388, 0.220567, 0.222062,
2154 0.235950, 0.224087, 0.228895, 0.208287, 0.235999, 0.208696,
2155 0.230367, 0.267667, 0.220484, 0.233402, 0.233815, 0.250455,
2156 0.253120, 0.219556, 0.230980, 0.236661, 0.222395, 0.226111,
2157 0.198315, 0.210555, 0.202813, 0.208594, 0.235976, 0.221490,
2158 0.243059, 0.204901, 0.216987, 0.229039, 0.231466, 0.221975,
2159 0.231220, 0.253638, 0.250448, 0.260291, 0.328345, 0.205739,
2160 0.222014, 0.251513, 0.279427, 0.270506, 0.248409, 0.222472,
2161 0.291632, 0.227796, 0.248769, 0.276896, 0.214742, 0.200139,
2162 0.230693, 0.226031, 0.268900, 0.185160, 0.245353, 0.205843,
2163 0.231155, 0.219122, 0.214811, 0.199763, 0.274179, 0.217598,
2164 0.274988, 0.237244, 0.211820, 0.225459, 0.252799, 0.235948,
2165 0.224986, 0.245385, 0.237770, 0.213373, 0.229737, 0.215487,
2166 0.234453, 0.249684, 0.239435, 0.250422, 0.257194, 0.214762,
2167 0.212266, 0.228988, 0.253798, 0.201607, 0.239946, 0.205245,
2168 0.231670, 0.212774, 0.206768, 0.231563, 0.189388, 0.227926,
2169 0.227203, 0.237754, 0.221628, 0.211138, 0.203322, 0.200985,
2170 0.231780, 0.220294, 0.232686, 0.234243, 0.218264, 0.255870,
2171 0.223213, 0.238670, 0.213713, 0.213064, 0.246700, 0.233446,
2172 0.221503, 0.206767, 0.200722, 0.226179, 0.237425, 0.239229,
2173 0.238611, 0.240419, 0.247806, 0.215923, 0.205298, 0.232778,
2174 0.272312, 0.226773, 0.258103, 0.223287, 0.269404, 0.203398,
2175 0.223782, 0.204213, 0.229664, 0.234040, 0.228419, 0.203936,
2176 0.263686, 0.199141, 0.236127, 0.214058, 0.204611, 0.224324,
2177 0.292140, 0.190735, 0.235157, 0.213018, 0.257085, 0.190554,
2178 0.203197, 0.213044, 0.237023, 0.214243, 0.193562, 0.262403,
2179 0.206256, 0.221396, 0.233588, 0.256611, 0.249731, 0.226683,
2180 0.199330, 0.251026, 0.222596, 0.201941, 0.186374, 0.221038,
2181 0.196555, 0.222560, 0.299419, 0.231979, 0.242924, 0.198310,
2182 0.217628, 0.235458, 0.278595, 0.218624, 0.277305, 0.239109,
2183 0.205600, 0.253715, 0.221173, 0.218195, 0.277647, 0.241974,
2184 0.268748, 0.268128, 0.292258, 0.249420, 0.191034, 0.219506,
2185 0.216502, 0.250677, 0.193386, 0.201310, 0.259464, 0.255351,
2186 0.269628, 0.221063, 0.294079, 0.196726, 0.233634, 0.221870,
2187 0.216236, 0.197259, 0.247433, 0.272765, 0.294079, 0.236336,
2188 0.206396, 0.238524, 0.247846, 0.269519, 0.237141, 0.230611,
2189 0.201712, 0.242225, 0.255565, 0.258738};
2196 if(M<Mmin || M>Mmax) return kTRUE;
2199 for(Int_t i=1; i<=nbin; i++)
2202 // printf("Mi> %f && Mi< %f\n", bin[i-1], bin[i]);
2210 else if(TMath::Abs(fEnergyCMS-7000)<1 ){
2211 const Int_t nbin=400;
2213 1.080000, 1.577300, 2.074600, 2.571900, 3.069200, 3.566500,
2214 4.063800, 4.561100, 5.058400, 5.555700, 6.053000, 6.550300,
2215 7.047600, 7.544900, 8.042200, 8.539500, 9.036800, 9.534100,
2216 10.031400, 10.528700, 11.026000, 11.523300, 12.020600, 12.517900,
2217 13.015200, 13.512500, 14.009800, 14.507100, 15.004400, 15.501700,
2218 15.999000, 16.496300, 16.993600, 17.490900, 17.988200, 18.485500,
2219 18.982800, 19.480100, 19.977400, 20.474700, 20.972000, 21.469300,
2220 21.966600, 22.463900, 22.961200, 23.458500, 23.955800, 24.453100,
2221 24.950400, 25.447700, 25.945000, 26.442300, 26.939600, 27.436900,
2222 27.934200, 28.431500, 28.928800, 29.426100, 29.923400, 30.420700,
2223 30.918000, 31.415300, 31.912600, 32.409900, 32.907200, 33.404500,
2224 33.901800, 34.399100, 34.896400, 35.393700, 35.891000, 36.388300,
2225 36.885600, 37.382900, 37.880200, 38.377500, 38.874800, 39.372100,
2226 39.869400, 40.366700, 40.864000, 41.361300, 41.858600, 42.355900,
2227 42.853200, 43.350500, 43.847800, 44.345100, 44.842400, 45.339700,
2228 45.837000, 46.334300, 46.831600, 47.328900, 47.826200, 48.323500,
2229 48.820800, 49.318100, 49.815400, 50.312700, 50.810000, 51.307300,
2230 51.804600, 52.301900, 52.799200, 53.296500, 53.793800, 54.291100,
2231 54.788400, 55.285700, 55.783000, 56.280300, 56.777600, 57.274900,
2232 57.772200, 58.269500, 58.766800, 59.264100, 59.761400, 60.258700,
2233 60.756000, 61.253300, 61.750600, 62.247900, 62.745200, 63.242500,
2234 63.739800, 64.237100, 64.734400, 65.231700, 65.729000, 66.226300,
2235 66.723600, 67.220900, 67.718200, 68.215500, 68.712800, 69.210100,
2236 69.707400, 70.204700, 70.702000, 71.199300, 71.696600, 72.193900,
2237 72.691200, 73.188500, 73.685800, 74.183100, 74.680400, 75.177700,
2238 75.675000, 76.172300, 76.669600, 77.166900, 77.664200, 78.161500,
2239 78.658800, 79.156100, 79.653400, 80.150700, 80.648000, 81.145300,
2240 81.642600, 82.139900, 82.637200, 83.134500, 83.631800, 84.129100,
2241 84.626400, 85.123700, 85.621000, 86.118300, 86.615600, 87.112900,
2242 87.610200, 88.107500, 88.604800, 89.102100, 89.599400, 90.096700,
2243 90.594000, 91.091300, 91.588600, 92.085900, 92.583200, 93.080500,
2244 93.577800, 94.075100, 94.572400, 95.069700, 95.567000, 96.064300,
2245 96.561600, 97.058900, 97.556200, 98.053500, 98.550800, 99.048100,
2246 99.545400, 100.042700, 100.540000, 101.037300, 101.534600, 102.031900,
2247 102.529200, 103.026500, 103.523800, 104.021100, 104.518400, 105.015700,
2248 105.513000, 106.010300, 106.507600, 107.004900, 107.502200, 107.999500,
2249 108.496800, 108.994100, 109.491400, 109.988700, 110.486000, 110.983300,
2250 111.480600, 111.977900, 112.475200, 112.972500, 113.469800, 113.967100,
2251 114.464400, 114.961700, 115.459000, 115.956300, 116.453600, 116.950900,
2252 117.448200, 117.945500, 118.442800, 118.940100, 119.437400, 119.934700,
2253 120.432000, 120.929300, 121.426600, 121.923900, 122.421200, 122.918500,
2254 123.415800, 123.913100, 124.410400, 124.907700, 125.405000, 125.902300,
2255 126.399600, 126.896900, 127.394200, 127.891500, 128.388800, 128.886100,
2256 129.383400, 129.880700, 130.378000, 130.875300, 131.372600, 131.869900,
2257 132.367200, 132.864500, 133.361800, 133.859100, 134.356400, 134.853700,
2258 135.351000, 135.848300, 136.345600, 136.842900, 137.340200, 137.837500,
2259 138.334800, 138.832100, 139.329400, 139.826700, 140.324000, 140.821300,
2260 141.318600, 141.815900, 142.313200, 142.810500, 143.307800, 143.805100,
2261 144.302400, 144.799700, 145.297000, 145.794300, 146.291600, 146.788900,
2262 147.286200, 147.783500, 148.280800, 148.778100, 149.275400, 149.772700,
2263 150.270000, 150.767300, 151.264600, 151.761900, 152.259200, 152.756500,
2264 153.253800, 153.751100, 154.248400, 154.745700, 155.243000, 155.740300,
2265 156.237600, 156.734900, 157.232200, 157.729500, 158.226800, 158.724100,
2266 159.221400, 159.718700, 160.216000, 160.713300, 161.210600, 161.707900,
2267 162.205200, 162.702500, 163.199800, 163.697100, 164.194400, 164.691700,
2268 165.189000, 165.686300, 166.183600, 166.680900, 167.178200, 167.675500,
2269 168.172800, 168.670100, 169.167400, 169.664700, 170.162000, 170.659300,
2270 171.156600, 171.653900, 172.151200, 172.648500, 173.145800, 173.643100,
2271 174.140400, 174.637700, 175.135000, 175.632300, 176.129600, 176.626900,
2272 177.124200, 177.621500, 178.118800, 178.616100, 179.113400, 179.610700,
2273 180.108000, 180.605300, 181.102600, 181.599900, 182.097200, 182.594500,
2274 183.091800, 183.589100, 184.086400, 184.583700, 185.081000, 185.578300,
2275 186.075600, 186.572900, 187.070200, 187.567500, 188.064800, 188.562100,
2276 189.059400, 189.556700, 190.054000, 190.551300, 191.048600, 191.545900,
2277 192.043200, 192.540500, 193.037800, 193.535100, 194.032400, 194.529700,
2278 195.027000, 195.524300, 196.021600, 196.518900, 197.016200, 197.513500,
2279 198.010800, 198.508100, 199.005400, 199.502700, 200.000000};
2281 1.000000, 0.798640, 0.770197, 0.717990, 0.699434, 0.692093,
2282 0.620940, 0.655294, 0.636496, 0.633483, 0.600518, 0.594548,
2283 0.550759, 0.550413, 0.528926, 0.525877, 0.506701, 0.504117,
2284 0.486973, 0.492930, 0.461833, 0.488695, 0.438874, 0.460274,
2285 0.451551, 0.454424, 0.449270, 0.387571, 0.427556, 0.428740,
2286 0.401682, 0.402260, 0.408961, 0.395417, 0.387426, 0.378602,
2287 0.357826, 0.359125, 0.348494, 0.363710, 0.356117, 0.340363,
2288 0.337637, 0.396084, 0.341777, 0.340551, 0.348838, 0.344014,
2289 0.340468, 0.317654, 0.355584, 0.326023, 0.373416, 0.312298,
2290 0.326522, 0.290540, 0.335557, 0.318689, 0.327544, 0.319501,
2291 0.331754, 0.312728, 0.282263, 0.274937, 0.303867, 0.307820,
2292 0.289344, 0.268934, 0.288908, 0.290018, 0.291369, 0.295242,
2293 0.289067, 0.277685, 0.267957, 0.267559, 0.320229, 0.265060,
2294 0.305931, 0.305352, 0.262064, 0.281879, 0.287780, 0.270033,
2295 0.270814, 0.276667, 0.271531, 0.275283, 0.258189, 0.287969,
2296 0.251247, 0.301527, 0.267230, 0.245860, 0.293125, 0.253421,
2297 0.272396, 0.243637, 0.236206, 0.278452, 0.246544, 0.263165,
2298 0.230484, 0.231102, 0.258281, 0.244707, 0.270111, 0.248295,
2299 0.246942, 0.245592, 0.272766, 0.254243, 0.199647, 0.262590,
2300 0.226710, 0.243836, 0.214153, 0.233535, 0.235207, 0.234245,
2301 0.247698, 0.248379, 0.241463, 0.265230, 0.223242, 0.236191,
2302 0.252700, 0.231865, 0.228352, 0.229101, 0.237385, 0.246485,
2303 0.254706, 0.245565, 0.224841, 0.257301, 0.240968, 0.202379,
2304 0.236943, 0.241683, 0.220809, 0.219014, 0.213015, 0.244470,
2305 0.221042, 0.198996, 0.225295, 0.264466, 0.200600, 0.228143,
2306 0.250138, 0.216784, 0.268002, 0.275734, 0.218144, 0.229866,
2307 0.235443, 0.208909, 0.215067, 0.189743, 0.216741, 0.242686,
2308 0.200234, 0.218882, 0.245991, 0.222815, 0.219576, 0.209773,
2309 0.205247, 0.203855, 0.202534, 0.192536, 0.223447, 0.225810,
2310 0.220583, 0.248421, 0.223424, 0.206033, 0.203130, 0.221225,
2311 0.223763, 0.239216, 0.252311, 0.206893, 0.228461, 0.233591,
2312 0.201512, 0.242382, 0.215147, 0.232578, 0.207116, 0.239213,
2313 0.196215, 0.184367, 0.235135, 0.189768, 0.274084, 0.206267,
2314 0.203283, 0.248828, 0.250214, 0.217465, 0.232080, 0.215150,
2315 0.207936, 0.176789, 0.191338, 0.188655, 0.196181, 0.223691,
2316 0.254257, 0.216874, 0.213536, 0.221399, 0.192024, 0.212534,
2317 0.218169, 0.226635, 0.201191, 0.212700, 0.211634, 0.232353,
2318 0.223636, 0.208605, 0.249132, 0.183681, 0.221842, 0.187136,
2319 0.203772, 0.249575, 0.217713, 0.205193, 0.193941, 0.219707,
2320 0.264226, 0.182105, 0.207183, 0.220845, 0.236571, 0.182390,
2321 0.212914, 0.186266, 0.195361, 0.217665, 0.204527, 0.188804,
2322 0.222832, 0.191193, 0.201073, 0.185616, 0.195011, 0.200828,
2323 0.181434, 0.233707, 0.202925, 0.211992, 0.173100, 0.205258,
2324 0.182695, 0.216520, 0.218422, 0.209490, 0.211257, 0.215801,
2325 0.220652, 0.211409, 0.195731, 0.194957, 0.198888, 0.234926,
2326 0.221377, 0.229822, 0.176700, 0.172322, 0.212265, 0.206133,
2327 0.170355, 0.253305, 0.198688, 0.240043, 0.225384, 0.174729,
2328 0.195820, 0.200093, 0.196912, 0.212308, 0.200490, 0.240415,
2329 0.159744, 0.173686, 0.223853, 0.213604, 0.193779, 0.179328,
2330 0.180873, 0.237481, 0.179640, 0.235857, 0.202847, 0.167690,
2331 0.245093, 0.215504, 0.205848, 0.184408, 0.201626, 0.209651,
2332 0.236306, 0.217803, 0.188534, 0.187861, 0.161663, 0.221718,
2333 0.152044, 0.190412, 0.173505, 0.208995, 0.174575, 0.180734,
2334 0.247704, 0.203388, 0.194184, 0.169679, 0.182703, 0.213402,
2335 0.191808, 0.178604, 0.178116, 0.198452, 0.174687, 0.169809,
2336 0.222851, 0.156811, 0.229170, 0.210359, 0.178557, 0.248570,
2337 0.208536, 0.192571, 0.178912, 0.224505, 0.170822, 0.205492,
2338 0.330973, 0.160924, 0.203724, 0.213255, 0.205827, 0.187162,
2339 0.181252, 0.191723, 0.238106, 0.182398, 0.202358, 0.212066,
2340 0.201255, 0.168159, 0.185219, 0.229176, 0.158222, 0.164092,
2341 0.215405, 0.200724, 0.234811, 0.184488, 0.213112, 0.198577,
2342 0.219622, 0.160671, 0.179349, 0.206681, 0.206091, 0.165618,
2343 0.165203, 0.174442, 0.179287, 0.187318, 0.163481, 0.217440,
2344 0.160381, 0.177025, 0.204385, 0.163676, 0.210733, 0.186519,
2345 0.230701, 0.191764, 0.185119, 0.190770, 0.242987, 0.186812,
2346 0.202906, 0.161935, 0.182426, 0.197922, 0.181475, 0.155903,
2347 0.175006, 0.223482, 0.202706, 0.218108};
2355 if(M<Mmin || M>Mmax) return kTRUE;
2358 for(Int_t i=1; i<=nbin; i++)
2361 // printf("Mi> %f && Mi< %f\n", bin[i-1], bin[i]);
2371 Bool_t AliGenPythia::IsFromHeavyFlavor(Int_t ipart)
2373 // Check if this is a heavy flavor decay product
2374 TParticle * part = (TParticle *) fParticles.At(ipart);
2375 Int_t mpdg = TMath::Abs(part->GetPdgCode());
2376 Int_t mfl = Int_t (mpdg / TMath::Power(10, Int_t(TMath::Log10(mpdg))));
2379 if (mfl >= 4 && mfl < 6) return kTRUE;
2381 Int_t imo = part->GetFirstMother()-1;
2382 TParticle* pm = part;
2384 // Heavy flavor hadron produced by generator
2386 pm = (TParticle*)fParticles.At(imo);
2387 mpdg = TMath::Abs(pm->GetPdgCode());
2388 mfl = Int_t (mpdg / TMath::Power(10, Int_t(TMath::Log10(mpdg))));
2389 if ((mfl > 3) && (mfl <6) && mpdg > 400) return kTRUE;
2390 imo = pm->GetFirstMother()-1;