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():
90 fDecayer(new AliDecayerPythia()),
98 fPhiMaxJet(2.* TMath::Pi()),
99 fJetReconstruction(kCell),
103 fPhiMaxGamma(2. * TMath::Pi()),
107 fPycellThreshold(0.),
109 fPycellMinEtJet(10.),
110 fPycellMaxRadius(1.),
111 fStackFillOpt(kFlavorSelection),
113 fFragmentation(kTRUE),
120 fTriggerMultiplicity(0),
121 fTriggerMultiplicityEta(0),
122 fTriggerMultiplicityPtMin(0),
123 fCountMode(kCountAll),
127 fFragPhotonInCalo(kFALSE),
129 fPhotonInCalo(kFALSE),
133 fCheckPHOSeta(kFALSE),
134 fFragPhotonOrPi0MinPt(0),
146 // Default Constructor
148 if (!AliPythiaRndm::GetPythiaRandom())
149 AliPythiaRndm::SetPythiaRandom(GetRandom());
152 AliGenPythia::AliGenPythia(Int_t npart)
164 fInteractionRate(0.),
178 fHadronisation(kTRUE),
179 fPatchOmegaDalitz(0),
181 fReadFromFile(kFALSE),
188 fDecayer(new AliDecayerPythia()),
189 fDebugEventFirst(-1),
196 fPhiMaxJet(2.* TMath::Pi()),
197 fJetReconstruction(kCell),
201 fPhiMaxGamma(2. * TMath::Pi()),
205 fPycellThreshold(0.),
207 fPycellMinEtJet(10.),
208 fPycellMaxRadius(1.),
209 fStackFillOpt(kFlavorSelection),
211 fFragmentation(kTRUE),
218 fTriggerMultiplicity(0),
219 fTriggerMultiplicityEta(0),
220 fTriggerMultiplicityPtMin(0),
221 fCountMode(kCountAll),
225 fFragPhotonInCalo(kFALSE),
227 fPhotonInCalo(kFALSE),
231 fCheckPHOSeta(kFALSE),
232 fFragPhotonOrPi0MinPt(0),
244 // default charm production at 5. 5 TeV
246 // structure function GRVHO
250 fTitle= "Particle Generator using PYTHIA";
252 // Set random number generator
253 if (!AliPythiaRndm::GetPythiaRandom())
254 AliPythiaRndm::SetPythiaRandom(GetRandom());
257 AliGenPythia::~AliGenPythia()
260 if(fEventsTime) delete fEventsTime;
263 void AliGenPythia::SetInteractionRate(Float_t rate,Float_t timewindow)
265 // Generate pileup using user specified rate
266 fInteractionRate = rate;
267 fTimeWindow = timewindow;
271 void AliGenPythia::GeneratePileup()
273 // Generate sub events time for pileup
275 if(fInteractionRate == 0.) {
276 Warning("GeneratePileup","Zero interaction specified. Skipping pileup generation.\n");
280 Int_t npart = NumberParticles();
282 Warning("GeneratePileup","Negative number of particles. Skipping pileup generation.\n");
286 if(fEventsTime) delete fEventsTime;
287 fEventsTime = new TArrayF(npart);
288 TArrayF &array = *fEventsTime;
289 for(Int_t ipart = 0; ipart < npart; ipart++)
292 Float_t eventtime = 0.;
295 eventtime += (AliPythiaRndm::GetPythiaRandom())->Exp(1./fInteractionRate);
296 if(eventtime > fTimeWindow) break;
297 array.Set(array.GetSize()+1);
298 array[array.GetSize()-1] = eventtime;
304 eventtime -= (AliPythiaRndm::GetPythiaRandom())->Exp(1./fInteractionRate);
305 if(TMath::Abs(eventtime) > fTimeWindow) break;
306 array.Set(array.GetSize()+1);
307 array[array.GetSize()-1] = eventtime;
310 SetNumberParticles(fEventsTime->GetSize());
313 void AliGenPythia::SetPycellParameters(Float_t etamax, Int_t neta, Int_t nphi,
314 Float_t thresh, Float_t etseed, Float_t minet, Float_t r)
316 // Set pycell parameters
317 fPycellEtaMax = etamax;
320 fPycellThreshold = thresh;
321 fPycellEtSeed = etseed;
322 fPycellMinEtJet = minet;
323 fPycellMaxRadius = r;
328 void AliGenPythia::SetEventListRange(Int_t eventFirst, Int_t eventLast)
330 // Set a range of event numbers, for which a table
331 // of generated particle will be printed
332 fDebugEventFirst = eventFirst;
333 fDebugEventLast = eventLast;
334 if (fDebugEventLast==-1) fDebugEventLast=fDebugEventFirst;
337 void AliGenPythia::Init()
341 SetMC(AliPythia::Instance());
342 fPythia=(AliPythia*) fMCEvGen;
345 fParentWeight=1./Float_t(fNpart);
349 fPythia->SetCKIN(3,fPtHardMin);
350 fPythia->SetCKIN(4,fPtHardMax);
351 fPythia->SetCKIN(7,fYHardMin);
352 fPythia->SetCKIN(8,fYHardMax);
354 if (fAProjectile > 0 && fATarget > 0) fPythia->SetNuclei(fAProjectile, fATarget, fNucPdf);
356 if (fFragmentation) {
357 fPythia->SetMSTP(111,1);
359 fPythia->SetMSTP(111,0);
363 // initial state radiation
364 fPythia->SetMSTP(61,fGinit);
365 // final state radiation
366 fPythia->SetMSTP(71,fGfinal);
369 fPythia->SetMSTP(91,1);
370 fPythia->SetPARP(91,fPtKick);
371 fPythia->SetPARP(93, 4. * fPtKick);
373 fPythia->SetMSTP(91,0);
378 fRL = AliRunLoader::Open(fkFileName, "Partons");
379 fRL->LoadKinematics();
385 fPythia->ProcInit(fProcess,fEnergyCMS,fStrucFunc, fItune);
386 // Forward Paramters to the AliPythia object
387 fDecayer->SetForceDecay(fForceDecay);
388 // Switch off Heavy Flavors on request
390 // Maximum number of quark flavours used in pdf
391 fPythia->SetMSTP(58, 3);
392 // Maximum number of flavors that can be used in showers
393 fPythia->SetMSTJ(45, 3);
394 // Switch off g->QQbar splitting in decay table
395 ((AliDecayerPythia*) fDecayer)->HeavyFlavourOff();
401 // Parent and Children Selection
404 case kPyOldUEQ2ordered:
405 case kPyOldUEQ2ordered2:
409 case kPyCharmUnforced:
410 case kPyCharmPbPbMNR:
413 case kPyCharmppMNRwmi:
414 fParentSelect[0] = 411;
415 fParentSelect[1] = 421;
416 fParentSelect[2] = 431;
417 fParentSelect[3] = 4122;
418 fParentSelect[4] = 4232;
419 fParentSelect[5] = 4132;
420 fParentSelect[6] = 4332;
426 fParentSelect[0] = 421;
429 case kPyDPlusPbPbMNR:
432 fParentSelect[0] = 411;
435 case kPyDPlusStrangePbPbMNR:
436 case kPyDPlusStrangepPbMNR:
437 case kPyDPlusStrangeppMNR:
438 fParentSelect[0] = 431;
441 case kPyLambdacppMNR:
442 fParentSelect[0] = 4122;
447 case kPyBeautyPbPbMNR:
448 case kPyBeautypPbMNR:
450 case kPyBeautyppMNRwmi:
451 fParentSelect[0]= 511;
452 fParentSelect[1]= 521;
453 fParentSelect[2]= 531;
454 fParentSelect[3]= 5122;
455 fParentSelect[4]= 5132;
456 fParentSelect[5]= 5232;
457 fParentSelect[6]= 5332;
460 case kPyBeautyUnforced:
461 fParentSelect[0] = 511;
462 fParentSelect[1] = 521;
463 fParentSelect[2] = 531;
464 fParentSelect[3] = 5122;
465 fParentSelect[4] = 5132;
466 fParentSelect[5] = 5232;
467 fParentSelect[6] = 5332;
472 fParentSelect[0] = 443;
475 case kPyMbAtlasTuneMC09:
477 case kPyMbWithDirectPhoton:
490 // JetFinder for Trigger
492 // Configure detector (EMCAL like)
494 fPythia->SetPARU(51, fPycellEtaMax);
495 fPythia->SetMSTU(51, fPycellNEta);
496 fPythia->SetMSTU(52, fPycellNPhi);
498 // Configure Jet Finder
500 fPythia->SetPARU(58, fPycellThreshold);
501 fPythia->SetPARU(52, fPycellEtSeed);
502 fPythia->SetPARU(53, fPycellMinEtJet);
503 fPythia->SetPARU(54, fPycellMaxRadius);
504 fPythia->SetMSTU(54, 2);
506 // This counts the total number of calls to Pyevnt() per run.
521 Warning("Init","SetNuclei used. Use SetProjectile + SetTarget instead. fDyBoost has been reset to 0\n");
524 fPythia->SetPARJ(200, 0.0);
525 fPythia->SetPARJ(199, 0.0);
526 fPythia->SetPARJ(198, 0.0);
527 fPythia->SetPARJ(197, 0.0);
530 fPythia->InitQuenching(0., 0.1, 0.6e6, 0);
534 // Nestor's change of the splittings
535 fPythia->SetPARJ(200, 0.8);
536 fPythia->SetMSTJ(41, 1); // QCD radiation only
537 fPythia->SetMSTJ(42, 2); // angular ordering
538 fPythia->SetMSTJ(44, 2); // option to run alpha_s
539 fPythia->SetMSTJ(47, 0); // No correction back to hard scattering element
540 fPythia->SetMSTJ(50, 0); // No coherence in first branching
541 fPythia->SetPARJ(82, 1.); // Cut off for parton showers
542 } else if (fQuench == 4) {
543 // Armesto-Cunqueiro-Salgado change of the splittings.
544 AliFastGlauber* glauber = AliFastGlauber::Instance();
546 //read and store transverse almonds corresponding to differnt
548 glauber->SetCentralityClass(0.,0.1);
549 fPythia->SetPARJ(200, 1.);
550 fPythia->SetPARJ(198, fQhat);
551 fPythia->SetPARJ(199, fLength);
552 fPythia->SetMSTJ(42, 2); // angular ordering
553 fPythia->SetMSTJ(44, 2); // option to run alpha_s
554 fPythia->SetPARJ(82, 1.); // Cut off for parton showers
558 void AliGenPythia::Generate()
560 // Generate one event
561 if (!fPythia) fPythia=(AliPythia*) fMCEvGen;
562 fDecayer->ForceDecay();
564 Double_t polar[3] = {0,0,0};
565 Double_t origin[3] = {0,0,0};
567 // converts from mm/c to s
568 const Float_t kconv=0.001/2.999792458e8;
578 // Set collision vertex position
579 if (fVertexSmear == kPerEvent) Vertex();
588 // Switch hadronisation off
590 fPythia->SetMSTJ(1, 0);
594 // Quenching comes through medium-modified splitting functions.
595 AliFastGlauber::Instance()->GetRandomBHard(bimp);
596 fPythia->SetPARJ(197, bimp);
601 // Either produce new event or read partons from file
603 if (!fReadFromFile) {
609 fNpartons = fPythia->GetN();
611 printf("Loading Event %d\n",AliRunLoader::Instance()->GetEventNumber());
612 fRL->GetEvent(AliRunLoader::Instance()->GetEventNumber());
614 LoadEvent(fRL->Stack(), 0 , 1);
619 // Run quenching routine
623 } else if (fQuench == 2){
624 fPythia->Pyquen(208., 0, 0.);
625 } else if (fQuench == 3) {
626 // Quenching is via multiplicative correction of the splittings
630 // Switch hadronisation on
632 if (fHadronisation) {
633 fPythia->SetMSTJ(1, 1);
635 // .. and perform hadronisation
636 // printf("Calling hadronisation %d\n", fPythia->GetN());
638 if (fPatchOmegaDalitz) {
639 fPythia->SetMDCY(fPythia->Pycomp(111) ,1, 0);
641 fPythia->DalitzDecays();
642 fPythia->SetMDCY(fPythia->Pycomp(111) ,1, 1);
647 fPythia->ImportParticles(&fParticles,"All");
649 if (TMath::Abs(fDyBoost) > 1.e-4) Boost();
650 if(TMath::Abs(fXingAngleY) > 1.e-10) BeamCrossAngle();
658 Int_t np = fParticles.GetEntriesFast();
660 if (np == 0) continue;
664 Int_t* pParent = new Int_t[np];
665 Int_t* pSelected = new Int_t[np];
666 Int_t* trackIt = new Int_t[np];
667 for (i = 0; i < np; i++) {
673 Int_t nc = 0; // Total n. of selected particles
674 Int_t nParents = 0; // Selected parents
675 Int_t nTkbles = 0; // Trackable particles
676 if (fProcess != kPyMbDefault &&
678 fProcess != kPyMbAtlasTuneMC09 &&
679 fProcess != kPyMbWithDirectPhoton &&
680 fProcess != kPyJets &&
681 fProcess != kPyDirectGamma &&
682 fProcess != kPyMbNonDiffr &&
683 fProcess != kPyMbMSEL1 &&
686 fProcess != kPyCharmppMNRwmi &&
687 fProcess != kPyBeautyppMNRwmi &&
688 fProcess != kPyBeautyJets) {
690 for (i = 0; i < np; i++) {
691 TParticle* iparticle = (TParticle *) fParticles.At(i);
692 Int_t ks = iparticle->GetStatusCode();
693 kf = CheckPDGCode(iparticle->GetPdgCode());
694 // No initial state partons
695 if (ks==21) continue;
697 // Heavy Flavor Selection
704 if (kfl > 100000) kfl %= 100000;
705 if (kfl > 10000) kfl %= 10000;
707 if (kfl > 10) kfl/=100;
709 if (kfl > 10) kfl/=10;
710 Int_t ipa = iparticle->GetFirstMother()-1;
713 // Establish mother daughter relation between heavy quarks and mesons
715 if (kf >= fFlavorSelect && kf <= 6) {
716 Int_t idau = iparticle->GetFirstDaughter() - 1;
718 TParticle* daughter = (TParticle *) fParticles.At(idau);
719 Int_t pdgD = daughter->GetPdgCode();
720 if (pdgD == 91 || pdgD == 92) {
721 Int_t jmin = daughter->GetFirstDaughter() - 1;
722 Int_t jmax = daughter->GetLastDaughter() - 1;
723 for (Int_t jp = jmin; jp <= jmax; jp++)
724 ((TParticle *) fParticles.At(jp))->SetFirstMother(i+1);
725 } // is string or cluster
731 TParticle * mother = (TParticle *) fParticles.At(ipa);
732 kfMo = TMath::Abs(mother->GetPdgCode());
735 // What to keep in Stack?
736 Bool_t flavorOK = kFALSE;
737 Bool_t selectOK = kFALSE;
739 if (kfl >= fFlavorSelect) flavorOK = kTRUE;
741 if (kfl > fFlavorSelect) {
745 if (kfl == fFlavorSelect) flavorOK = kTRUE;
747 switch (fStackFillOpt) {
748 case kFlavorSelection:
751 case kParentSelection:
752 if (ParentSelected(kf) || kf <= 10) selectOK = kTRUE;
755 if (flavorOK && selectOK) {
757 // Heavy flavor hadron or quark
759 // Kinematic seletion on final state heavy flavor mesons
760 if (ParentSelected(kf) && !KinematicSelection(iparticle, 0))
765 if (ParentSelected(kf)) ++nParents; // Update parent count
766 // printf("\n particle (HF) %d %d %d", i, pSelected[i], kf);
768 // Kinematic seletion on decay products
769 if (fCutOnChild && ParentSelected(kfMo) && ChildSelected(kf)
770 && !KinematicSelection(iparticle, 1))
776 // Select if mother was selected and is not tracked
778 if (pSelected[ipa] &&
779 !trackIt[ipa] && // mother will be tracked ?
780 kfMo != 5 && // mother is b-quark, don't store fragments
781 kfMo != 4 && // mother is c-quark, don't store fragments
782 kf != 92) // don't store string
785 // Semi-stable or de-selected: diselect decay products:
788 if (pSelected[i] == -1 || fDecayer->GetLifetime(kf) > fMaxLifeTime)
790 Int_t ipF = iparticle->GetFirstDaughter();
791 Int_t ipL = iparticle->GetLastDaughter();
792 if (ipF > 0) for (j = ipF-1; j < ipL; j++) pSelected[j] = -1;
794 // printf("\n particle (decay) %d %d %d", i, pSelected[i], kf);
795 pSelected[i] = (pSelected[i] == -1) ? 0 : 1;
798 if (pSelected[i] == -1) pSelected[i] = 0;
799 if (!pSelected[i]) continue;
800 // Count quarks only if you did not include fragmentation
801 if (fFragmentation && kf <= 10) continue;
804 // Decision on tracking
807 // Track final state particle
808 if (ks == 1) trackIt[i] = 1;
809 // Track semi-stable particles
810 if ((ks == 1) || (fDecayer->GetLifetime(kf) > fMaxLifeTime)) trackIt[i] = 1;
811 // Track particles selected by process if undecayed.
812 if (fForceDecay == kNoDecay) {
813 if (ParentSelected(kf)) trackIt[i] = 1;
815 if (ParentSelected(kf)) trackIt[i] = 0;
817 if (trackIt[i] == 1) ++nTkbles; // Update trackable counter
821 } // particle selection loop
823 for (i = 0; i<np; i++) {
824 if (!pSelected[i]) continue;
825 TParticle * iparticle = (TParticle *) fParticles.At(i);
826 kf = CheckPDGCode(iparticle->GetPdgCode());
827 Int_t ks = iparticle->GetStatusCode();
828 p[0] = iparticle->Px();
829 p[1] = iparticle->Py();
830 p[2] = iparticle->Pz();
831 p[3] = iparticle->Energy();
833 origin[0] = fVertex[0]+iparticle->Vx()/10; // [cm]
834 origin[1] = fVertex[1]+iparticle->Vy()/10; // [cm]
835 origin[2] = fVertex[2]+iparticle->Vz()/10; // [cm]
837 Float_t tof = fTime + kconv*iparticle->T();
838 Int_t ipa = iparticle->GetFirstMother()-1;
839 Int_t iparent = (ipa > -1) ? pParent[ipa] : -1;
841 PushTrack(fTrackIt*trackIt[i], iparent, kf,
842 p[0], p[1], p[2], p[3],
843 origin[0], origin[1], origin[2], tof,
844 polar[0], polar[1], polar[2],
845 kPPrimary, nt, 1., ks);
862 switch (fCountMode) {
864 // printf(" Count all \n");
868 // printf(" Count parents \n");
871 case kCountTrackables:
872 // printf(" Count trackable \n");
876 if (jev >= fNpart || fNpart == -1) {
877 fKineBias=Float_t(fNpart)/Float_t(fTrials);
879 fQ += fPythia->GetVINT(51);
880 fX1 += fPythia->GetVINT(41);
881 fX2 += fPythia->GetVINT(42);
882 fTrialsRun += fTrials;
889 SetHighWaterMark(nt);
890 // adjust weight due to kinematic selection
893 fXsection=fPythia->GetPARI(1);
896 Int_t AliGenPythia::GenerateMB()
899 // Min Bias selection and other global selections
901 Int_t i, kf, nt, iparent;
904 Double_t polar[3] = {0,0,0};
905 Double_t origin[3] = {0,0,0};
906 // converts from mm/c to s
907 const Float_t kconv=0.001/2.999792458e8;
911 Int_t np = (fHadronisation) ? fParticles.GetEntriesFast() : fNpartons;
915 Int_t* pParent = new Int_t[np];
916 for (i=0; i< np; i++) pParent[i] = -1;
917 if (fProcess == kPyJets || fProcess == kPyDirectGamma || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi) {
918 TParticle* jet1 = (TParticle *) fParticles.At(6);
919 TParticle* jet2 = (TParticle *) fParticles.At(7);
920 if (!CheckTrigger(jet1, jet2)) {
926 // Select jets with fragmentation photon or pi0 going to PHOS or EMCAL
927 if (fProcess == kPyJets && (fFragPhotonInCalo || fPi0InCalo) ) {
932 if (fFragPhotonInCalo) pdg = 22 ; // Photon
933 else if (fPi0InCalo) pdg = 111 ; // Pi0
935 for (i=0; i< np; i++) {
936 TParticle* iparticle = (TParticle *) fParticles.At(i);
937 if(iparticle->GetStatusCode()==1 && iparticle->GetPdgCode()==pdg &&
938 iparticle->Pt() > fFragPhotonOrPi0MinPt){
939 Int_t imother = iparticle->GetFirstMother() - 1;
940 TParticle* pmother = (TParticle *) fParticles.At(imother);
942 (pdg == 22 && pmother->GetStatusCode() != 11)) //No photon from hadron decay
944 Float_t phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees
945 Float_t eta =TMath::Abs(iparticle->Eta()); //in calos etamin=-etamax
946 if((fCheckEMCAL && IsInEMCAL(phi,eta)) ||
947 (fCheckPHOS && IsInPHOS(phi,eta)) )
958 // Select beauty jets with electron in EMCAL
959 if (fProcess == kPyBeautyJets && fEleInEMCAL) {
963 Int_t pdg = 11; //electron
968 for (i=0; i< np; i++) {
969 TParticle* iparticle = (TParticle *) fParticles.At(i);
970 if(iparticle->GetStatusCode()==1 && TMath::Abs(iparticle->GetPdgCode())==pdg &&
971 iparticle->Pt() > fElectronMinPt){
972 pt = iparticle->Pt();
973 phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees
974 eta =TMath::Abs(iparticle->Eta()); //in calos etamin=-etamax
975 if(IsInEMCAL(phi,eta))
984 AliDebug(5,Form("Found an electron jet (pt,eta,phi) = (%f,%f,%f)",pt,eta,phi));
986 // Check for diffraction
988 if(fItune==320 && ( (TMath::Abs(fEnergyCMS - 900) < 1) || (TMath::Abs(fEnergyCMS - 2760) < 1) || (TMath::Abs(fEnergyCMS - 7000) < 1)) ) {
989 if(!CheckDiffraction()) {
996 // Check for minimum multiplicity
997 if (fTriggerMultiplicity > 0) {
998 Int_t multiplicity = 0;
999 for (i = 0; i < np; i++) {
1000 TParticle * iparticle = (TParticle *) fParticles.At(i);
1002 Int_t statusCode = iparticle->GetStatusCode();
1004 // Initial state particle
1005 if (statusCode != 1)
1008 if (fTriggerMultiplicityEta > 0 && TMath::Abs(iparticle->Eta()) > fTriggerMultiplicityEta)
1011 if (iparticle->Pt() < fTriggerMultiplicityPtMin)
1014 TParticlePDG* pdgPart = iparticle->GetPDG();
1015 if (pdgPart && pdgPart->Charge() == 0)
1021 if (multiplicity < fTriggerMultiplicity) {
1025 Printf("Triggered on event with multiplicity of %d >= %d", multiplicity, fTriggerMultiplicity);
1028 // Select events with a photon pt > min pt going to PHOS eta acceptance or exactly PHOS eta phi
1029 if ((fProcess == kPyJets || fProcess == kPyDirectGamma) && fPhotonInCalo && (fCheckPHOSeta || fCheckPHOS)){
1031 Bool_t okd = kFALSE;
1035 for (i=0; i< np; i++) {
1036 TParticle* iparticle = (TParticle *) fParticles.At(i);
1037 Float_t phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees
1038 Float_t eta =TMath::Abs(iparticle->Eta());//in calos etamin=-etamax
1040 if(iparticle->GetStatusCode() == 1
1041 && iparticle->GetPdgCode() == pdg
1042 && iparticle->Pt() > fPhotonMinPt
1045 // first check if the photon is in PHOS phi
1046 if(IsInPHOS(phi,eta)){
1050 if(fCheckPHOSeta) iphcand = i; // candiate photon to rotate in phi
1055 if(!okd && iphcand != -1) // execute rotation in phi
1056 RotatePhi(iphcand,okd);
1064 if (fTriggerParticle) {
1065 Bool_t triggered = kFALSE;
1066 for (i = 0; i < np; i++) {
1067 TParticle * iparticle = (TParticle *) fParticles.At(i);
1068 kf = CheckPDGCode(iparticle->GetPdgCode());
1069 if (kf != fTriggerParticle) continue;
1070 if (iparticle->Pt() == 0.) continue;
1071 if (TMath::Abs(iparticle->Eta()) > fTriggerEta) continue;
1082 // Check if there is a ccbar or bbbar pair with at least one of the two
1083 // in fYMin < y < fYMax
1085 if (fProcess == kPyCharmppMNRwmi || fProcess == kPyBeautyppMNRwmi || fProcess == kPyBeautyJets) {
1086 TParticle *partCheck;
1088 Bool_t theQ=kFALSE,theQbar=kFALSE,inYcut=kFALSE;
1089 Bool_t theChild=kFALSE;
1090 Bool_t triggered=kFALSE;
1092 Int_t pdg,mpdg,mpdgUpperFamily;
1093 for(i=0; i<np; i++) {
1094 partCheck = (TParticle*)fParticles.At(i);
1095 pdg = partCheck->GetPdgCode();
1096 if(TMath::Abs(pdg) == fFlavorSelect) { // quark
1097 if(pdg>0) { theQ=kTRUE; } else { theQbar=kTRUE; }
1098 y = 0.5*TMath::Log((partCheck->Energy()+partCheck->Pz()+1.e-13)/
1099 (partCheck->Energy()-partCheck->Pz()+1.e-13));
1100 if(y>fYMin && y<fYMax) inYcut=kTRUE;
1102 if(fTriggerParticle) {
1103 if(TMath::Abs(pdg)==fTriggerParticle) triggered=kTRUE;
1105 if(fCutOnChild && TMath::Abs(pdg) == fPdgCodeParticleforAcceptanceCut) {
1106 Int_t mi = partCheck->GetFirstMother() - 1;
1108 mother = (TParticle*)fParticles.At(mi);
1109 mpdg=TMath::Abs(mother->GetPdgCode());
1110 mpdgUpperFamily=(mpdg>1000 ? mpdg+1000 : mpdg+100); // keep e from c from b
1111 if ( ParentSelected(mpdg) ||
1112 (fFlavorSelect==5 && ParentSelected(mpdgUpperFamily))) {
1113 if (KinematicSelection(partCheck,1)) {
1119 if (!theQ || !theQbar || !inYcut) { // one of the c/b conditions not satisfied
1123 if (fCutOnChild && !theChild) { // one of the child conditions not satisfied
1127 if(fTriggerParticle && !triggered) { // particle requested is not produced
1134 //Introducing child cuts in case kPyW, kPyZ, kPyMb, and kPyMbNonDiff
1135 if ( (fProcess == kPyW ||
1137 fProcess == kPyMbDefault ||
1138 fProcess == kPyMb ||
1139 fProcess == kPyMbAtlasTuneMC09 ||
1140 fProcess == kPyMbWithDirectPhoton ||
1141 fProcess == kPyMbNonDiffr)
1142 && (fCutOnChild == 1) ) {
1143 if ( !CheckKinematicsOnChild() ) {
1150 for (i = 0; i < np; i++) {
1152 TParticle * iparticle = (TParticle *) fParticles.At(i);
1153 kf = CheckPDGCode(iparticle->GetPdgCode());
1154 Int_t ks = iparticle->GetStatusCode();
1155 Int_t km = iparticle->GetFirstMother();
1156 if ((ks == 1 && kf!=0 && KinematicSelection(iparticle, 0)) ||
1158 ((fProcess == kPyJets || fProcess == kPyBeautyJets) && ks == 21 && km == 0 && i>1)) {
1160 if (ks == 1) trackIt = 1;
1161 Int_t ipa = iparticle->GetFirstMother()-1;
1163 iparent = (ipa > -1) ? pParent[ipa] : -1;
1166 // store track information
1167 p[0] = iparticle->Px();
1168 p[1] = iparticle->Py();
1169 p[2] = iparticle->Pz();
1170 p[3] = iparticle->Energy();
1173 origin[0] = fVertex[0]+iparticle->Vx()/10; // [cm]
1174 origin[1] = fVertex[1]+iparticle->Vy()/10; // [cm]
1175 origin[2] = fVertex[2]+iparticle->Vz()/10; // [cm]
1177 Float_t tof = fTime + fEventTime + kconv * iparticle->T();
1179 PushTrack(fTrackIt*trackIt, iparent, kf,
1180 p[0], p[1], p[2], p[3],
1181 origin[0], origin[1], origin[2], tof,
1182 polar[0], polar[1], polar[2],
1183 kPPrimary, nt, 1., ks);
1187 SetHighWaterMark(nt);
1189 } // select particle
1198 void AliGenPythia::FinishRun()
1200 // Print x-section summary
1209 printf("\nTotal number of Pyevnt() calls %d\n", fTrialsRun);
1210 printf("\nMean Q, x1, x2: %f %f %f\n", fQ, fX1, fX2);
1213 void AliGenPythia::AdjustWeights() const
1215 // Adjust the weights after generation of all events
1219 Int_t ntrack=gAlice->GetMCApp()->GetNtrack();
1220 for (Int_t i=0; i<ntrack; i++) {
1221 part= gAlice->GetMCApp()->Particle(i);
1222 part->SetWeight(part->GetWeight()*fKineBias);
1227 void AliGenPythia::SetNuclei(Int_t a1, Int_t a2, Int_t pdfset)
1229 // Treat protons as inside nuclei with mass numbers a1 and a2
1233 fNucPdf = pdfset; // 0 EKS98 1 EPS08
1238 void AliGenPythia::MakeHeader()
1241 // Make header for the simulated event
1244 if (gAlice->GetEvNumber()>=fDebugEventFirst &&
1245 gAlice->GetEvNumber()<=fDebugEventLast) fPythia->Pylist(2);
1248 // Builds the event header, to be called after each event
1249 if (fHeader) delete fHeader;
1250 fHeader = new AliGenPythiaEventHeader("Pythia");
1254 // if(fProcDiff == 92 || fProcDiff == 93) printf("\n\n\n\n\n");
1255 // printf("fPythia->GetMSTI(1) = %d fProcDiff = %d\n",fPythia->GetMSTI(1), fProcDiff);
1256 ((AliGenPythiaEventHeader*) fHeader)->SetProcessType(fProcDiff);
1259 ((AliGenPythiaEventHeader*) fHeader)->SetProcessType(fPythia->GetMSTI(1));
1262 ((AliGenPythiaEventHeader*) fHeader)->SetTrials(fTrials);
1265 fHeader->SetPrimaryVertex(fVertex);
1266 fHeader->SetInteractionTime(fTime+fEventTime);
1268 // Number of primaries
1269 fHeader->SetNProduced(fNprimaries);
1271 // Jets that have triggered
1273 //Need to store jets for b-jet studies too!
1274 if (fProcess == kPyJets || fProcess == kPyDirectGamma || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi)
1277 Float_t jets[4][10];
1278 GetJets(njet, ntrig, jets);
1281 for (Int_t i = 0; i < ntrig; i++) {
1282 ((AliGenPythiaEventHeader*) fHeader)->AddJet(jets[0][i], jets[1][i], jets[2][i],
1287 // Copy relevant information from external header, if present.
1292 AliGenPythiaEventHeader* exHeader = (AliGenPythiaEventHeader*) (fRL->GetHeader()->GenEventHeader());
1293 for (Int_t i = 0; i < exHeader->NTriggerJets(); i++)
1295 printf("Adding Jet %d %d \n", i, exHeader->NTriggerJets());
1298 exHeader->TriggerJet(i, uqJet);
1299 ((AliGenPythiaEventHeader*) fHeader)->AddUQJet(uqJet[0], uqJet[1], uqJet[2], uqJet[3]);
1303 // Store quenching parameters
1306 Double_t z[4] = {0.};
1312 fPythia->GetQuenchingParameters(xp, yp, z);
1313 } else if (fQuench == 2){
1315 Double_t r1 = PARIMP.rb1;
1316 Double_t r2 = PARIMP.rb2;
1317 Double_t b = PARIMP.b1;
1318 Double_t r = 0.5 * TMath::Sqrt(2. * (r1 * r1 + r2 * r2) - b * b);
1319 Double_t phi = PARIMP.psib1;
1320 xp = r * TMath::Cos(phi);
1321 yp = r * TMath::Sin(phi);
1323 } else if (fQuench == 4) {
1327 AliFastGlauber::Instance()->GetSavedXY(xy);
1328 AliFastGlauber::Instance()->GetSavedI0I1(i0i1);
1331 ((AliGenPythiaEventHeader*) fHeader)->SetImpactParameter(fImpact);
1334 ((AliGenPythiaEventHeader*) fHeader)->SetXYJet(xp, yp);
1335 ((AliGenPythiaEventHeader*) fHeader)->SetZQuench(z);
1339 ((AliGenPythiaEventHeader*) fHeader)->SetPtHard(fPythia->GetVINT(47));
1347 Bool_t AliGenPythia::CheckTrigger(const TParticle* jet1, const TParticle* jet2)
1349 // Check the kinematic trigger condition
1352 eta[0] = jet1->Eta();
1353 eta[1] = jet2->Eta();
1355 phi[0] = jet1->Phi();
1356 phi[1] = jet2->Phi();
1358 pdg[0] = jet1->GetPdgCode();
1359 pdg[1] = jet2->GetPdgCode();
1360 Bool_t triggered = kFALSE;
1362 if (fProcess == kPyJets || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi) {
1365 Float_t jets[4][10];
1367 // Use Pythia clustering on parton level to determine jet axis
1369 GetJets(njets, ntrig, jets);
1371 if (ntrig || fEtMinJet == 0.) triggered = kTRUE;
1376 if (pdg[0] == kGamma) {
1380 //Check eta range first...
1381 if ((eta[ij] < fEtaMaxJet && eta[ij] > fEtaMinJet) &&
1382 (eta[ig] < fEtaMaxGamma && eta[ig] > fEtaMinGamma))
1384 //Eta is okay, now check phi range
1385 if ((phi[ij] < fPhiMaxJet && phi[ij] > fPhiMinJet) &&
1386 (phi[ig] < fPhiMaxGamma && phi[ig] > fPhiMinGamma))
1397 Bool_t AliGenPythia::CheckKinematicsOnChild(){
1399 //Checking Kinematics on Child (status code 1, particle code ?, kin cuts
1401 Bool_t checking = kFALSE;
1402 Int_t j, kcode, ks, km;
1403 Int_t nPartAcc = 0; //number of particles in the acceptance range
1404 Int_t numberOfAcceptedParticles = 1;
1405 if (fNumberOfAcceptedParticles != 0) { numberOfAcceptedParticles = fNumberOfAcceptedParticles; }
1406 Int_t npart = fParticles.GetEntriesFast();
1408 for (j = 0; j<npart; j++) {
1409 TParticle * jparticle = (TParticle *) fParticles.At(j);
1410 kcode = TMath::Abs( CheckPDGCode(jparticle->GetPdgCode()) );
1411 ks = jparticle->GetStatusCode();
1412 km = jparticle->GetFirstMother();
1414 if( (ks == 1) && (kcode == fPdgCodeParticleforAcceptanceCut) && (KinematicSelection(jparticle,1)) ){
1417 if( numberOfAcceptedParticles <= nPartAcc){
1426 void AliGenPythia::LoadEvent(AliStack* stack, Int_t flag, Int_t reHadr)
1429 // Load event into Pythia Common Block
1432 Int_t npart = stack -> GetNprimary();
1436 (fPythia->GetPyjets())->N = npart;
1438 n0 = (fPythia->GetPyjets())->N;
1439 (fPythia->GetPyjets())->N = n0 + npart;
1443 for (Int_t part = 0; part < npart; part++) {
1444 TParticle *mPart = stack->Particle(part);
1446 Int_t kf = mPart->GetPdgCode();
1447 Int_t ks = mPart->GetStatusCode();
1448 Int_t idf = mPart->GetFirstDaughter();
1449 Int_t idl = mPart->GetLastDaughter();
1452 if (ks == 11 || ks == 12) {
1459 Float_t px = mPart->Px();
1460 Float_t py = mPart->Py();
1461 Float_t pz = mPart->Pz();
1462 Float_t e = mPart->Energy();
1463 Float_t m = mPart->GetCalcMass();
1466 (fPythia->GetPyjets())->P[0][part+n0] = px;
1467 (fPythia->GetPyjets())->P[1][part+n0] = py;
1468 (fPythia->GetPyjets())->P[2][part+n0] = pz;
1469 (fPythia->GetPyjets())->P[3][part+n0] = e;
1470 (fPythia->GetPyjets())->P[4][part+n0] = m;
1472 (fPythia->GetPyjets())->K[1][part+n0] = kf;
1473 (fPythia->GetPyjets())->K[0][part+n0] = ks;
1474 (fPythia->GetPyjets())->K[3][part+n0] = idf + 1;
1475 (fPythia->GetPyjets())->K[4][part+n0] = idl + 1;
1476 (fPythia->GetPyjets())->K[2][part+n0] = mPart->GetFirstMother() + 1;
1480 void AliGenPythia::LoadEvent(const TObjArray* stack, Int_t flag, Int_t reHadr)
1483 // Load event into Pythia Common Block
1486 Int_t npart = stack -> GetEntries();
1490 (fPythia->GetPyjets())->N = npart;
1492 n0 = (fPythia->GetPyjets())->N;
1493 (fPythia->GetPyjets())->N = n0 + npart;
1497 for (Int_t part = 0; part < npart; part++) {
1498 TParticle *mPart = dynamic_cast<TParticle *>(stack->At(part));
1499 if (!mPart) continue;
1501 Int_t kf = mPart->GetPdgCode();
1502 Int_t ks = mPart->GetStatusCode();
1503 Int_t idf = mPart->GetFirstDaughter();
1504 Int_t idl = mPart->GetLastDaughter();
1507 if (ks == 11 || ks == 12) {
1514 Float_t px = mPart->Px();
1515 Float_t py = mPart->Py();
1516 Float_t pz = mPart->Pz();
1517 Float_t e = mPart->Energy();
1518 Float_t m = mPart->GetCalcMass();
1521 (fPythia->GetPyjets())->P[0][part+n0] = px;
1522 (fPythia->GetPyjets())->P[1][part+n0] = py;
1523 (fPythia->GetPyjets())->P[2][part+n0] = pz;
1524 (fPythia->GetPyjets())->P[3][part+n0] = e;
1525 (fPythia->GetPyjets())->P[4][part+n0] = m;
1527 (fPythia->GetPyjets())->K[1][part+n0] = kf;
1528 (fPythia->GetPyjets())->K[0][part+n0] = ks;
1529 (fPythia->GetPyjets())->K[3][part+n0] = idf + 1;
1530 (fPythia->GetPyjets())->K[4][part+n0] = idl + 1;
1531 (fPythia->GetPyjets())->K[2][part+n0] = mPart->GetFirstMother() + 1;
1536 void AliGenPythia::RecJetsUA1(Int_t& njets, Float_t jets [4][50])
1539 // Calls the Pythia jet finding algorithm to find jets in the current event
1544 Int_t n = fPythia->GetN();
1548 fPythia->Pycell(njets);
1550 for (i = 0; i < njets; i++) {
1551 Float_t px = (fPythia->GetPyjets())->P[0][n+i];
1552 Float_t py = (fPythia->GetPyjets())->P[1][n+i];
1553 Float_t pz = (fPythia->GetPyjets())->P[2][n+i];
1554 Float_t e = (fPythia->GetPyjets())->P[3][n+i];
1565 void AliGenPythia::GetJets(Int_t& nJets, Int_t& nJetsTrig, Float_t jets[4][10])
1568 // Calls the Pythia clustering algorithm to find jets in the current event
1570 Int_t n = fPythia->GetN();
1573 if (fJetReconstruction == kCluster) {
1575 // Configure cluster algorithm
1577 fPythia->SetPARU(43, 2.);
1578 fPythia->SetMSTU(41, 1);
1580 // Call cluster algorithm
1582 fPythia->Pyclus(nJets);
1584 // Loading jets from common block
1590 fPythia->Pycell(nJets);
1594 for (i = 0; i < nJets; i++) {
1595 Float_t px = (fPythia->GetPyjets())->P[0][n+i];
1596 Float_t py = (fPythia->GetPyjets())->P[1][n+i];
1597 Float_t pz = (fPythia->GetPyjets())->P[2][n+i];
1598 Float_t e = (fPythia->GetPyjets())->P[3][n+i];
1599 Float_t pt = TMath::Sqrt(px * px + py * py);
1600 Float_t phi = TMath::Pi() + TMath::ATan2(-py, -px);
1601 Float_t theta = TMath::ATan2(pt,pz);
1602 Float_t et = e * TMath::Sin(theta);
1603 Float_t eta = -TMath::Log(TMath::Tan(theta / 2.));
1605 eta > fEtaMinJet && eta < fEtaMaxJet &&
1606 phi > fPhiMinJet && phi < fPhiMaxJet &&
1607 et > fEtMinJet && et < fEtMaxJet
1610 jets[0][nJetsTrig] = px;
1611 jets[1][nJetsTrig] = py;
1612 jets[2][nJetsTrig] = pz;
1613 jets[3][nJetsTrig] = e;
1615 // printf("\n........-Jet #%d: %10.3f %10.3f %10.3f %10.3f \n", i, pt, et, eta, phi * kRaddeg);
1617 // printf("\n........-Jet #%d: %10.3f %10.3f %10.3f %10.3f \n", i, pt, et, eta, phi * kRaddeg);
1622 void AliGenPythia::GetSubEventTime()
1624 // Calculates time of the next subevent
1627 TArrayF &array = *fEventsTime;
1628 fEventTime = array[fCurSubEvent++];
1630 // printf(" Event time: %d %f %p",fCurSubEvent,fEventTime,fEventsTime);
1634 Bool_t AliGenPythia::IsInEMCAL(Float_t phi, Float_t eta) const
1636 // Is particle in EMCAL acceptance?
1637 // phi in degrees, etamin=-etamax
1638 if(phi > fEMCALMinPhi && phi < fEMCALMaxPhi &&
1645 Bool_t AliGenPythia::IsInPHOS(Float_t phi, Float_t eta) const
1647 // Is particle in PHOS acceptance?
1648 // Acceptance slightly larger considered.
1649 // phi in degrees, etamin=-etamax
1650 if(phi > fPHOSMinPhi && phi < fPHOSMaxPhi &&
1657 void AliGenPythia::RotatePhi(Int_t iphcand, Bool_t& okdd)
1659 //calculate the new position random between fPHOSMinPhi and fPHOSMaxPhi
1660 Double_t phiPHOSmin = TMath::Pi()*fPHOSMinPhi/180;
1661 Double_t phiPHOSmax = TMath::Pi()*fPHOSMaxPhi/180;
1662 Double_t phiPHOS = (AliPythiaRndm::GetPythiaRandom())->Uniform(phiPHOSmin,phiPHOSmax);
1664 //calculate deltaphi
1665 TParticle* ph = (TParticle *) fParticles.At(iphcand);
1666 Double_t phphi = ph->Phi();
1667 Double_t deltaphi = phiPHOS - phphi;
1671 //loop for all particles and produce the phi rotation
1672 Int_t np = (fHadronisation) ? fParticles.GetEntriesFast() : fNpartons;
1673 Double_t oldphi, newphi;
1674 Double_t newVx, newVy, r, vZ, time;
1675 Double_t newPx, newPy, pt, pz, e;
1676 for(Int_t i=0; i< np; i++) {
1677 TParticle* iparticle = (TParticle *) fParticles.At(i);
1678 oldphi = iparticle->Phi();
1679 newphi = oldphi + deltaphi;
1680 if(newphi < 0) newphi = 2*TMath::Pi() + newphi; // correct angle
1681 if(newphi > 2*TMath::Pi()) newphi = newphi - 2*TMath::Pi(); // correct angle
1684 newVx = r * TMath::Cos(newphi);
1685 newVy = r * TMath::Sin(newphi);
1686 vZ = iparticle->Vz(); // don't transform
1687 time = iparticle->T(); // don't transform
1689 pt = iparticle->Pt();
1690 newPx = pt * TMath::Cos(newphi);
1691 newPy = pt * TMath::Sin(newphi);
1692 pz = iparticle->Pz(); // don't transform
1693 e = iparticle->Energy(); // don't transform
1696 iparticle->SetProductionVertex(newVx, newVy, vZ, time);
1697 iparticle->SetMomentum(newPx, newPy, pz, e);
1699 } //end particle loop
1701 // now let's check that we put correctly the candidate photon in PHOS
1702 Float_t phi = ph->Phi()*180./TMath::Pi(); //Convert to degrees
1703 Float_t eta =TMath::Abs(ph->Eta());//in calos etamin=-etamax
1704 if(IsInPHOS(phi,eta))
1711 Bool_t AliGenPythia::CheckDiffraction()
1713 // use this method only with Perugia-0 tune!
1717 Int_t np = (fHadronisation) ? fParticles.GetEntriesFast() : fNpartons;
1723 Double_t y2 = -1e10;
1725 const Int_t kNstable=20;
1726 const Int_t pdgStable[20] = {
1729 12, // Electron Neutrino
1731 14, // Muon Neutrino
1742 3112, // Sigma Minus
1749 for (Int_t i = 0; i < np; i++) {
1750 TParticle * part = (TParticle *) fParticles.At(i);
1752 Int_t statusCode = part->GetStatusCode();
1754 // Initial state particle
1755 if (statusCode != 1)
1758 Int_t pdg = TMath::Abs(part->GetPdgCode());
1759 Bool_t isStable = kFALSE;
1760 for (Int_t i1 = 0; i1 < kNstable; i1++) {
1761 if (pdg == pdgStable[i1]) {
1769 Double_t y = part->Y();
1783 if(iPart1<0 || iPart2<0) return kFALSE;
1788 TParticle * part1 = (TParticle *) fParticles.At(iPart1);
1789 TParticle * part2 = (TParticle *) fParticles.At(iPart2);
1791 Int_t pdg1 = part1->GetPdgCode();
1792 Int_t pdg2 = part2->GetPdgCode();
1796 if (pdg1 == 2212 && pdg2 == 2212)
1804 if((AliPythiaRndm::GetPythiaRandom())->Uniform(0.,1.)>0.5) iPart = iPart2;
1807 else if (pdg1 == 2212)
1809 else if (pdg2 == 2212)
1818 TParticle * part = (TParticle *) fParticles.At(iPart);
1819 Double_t E= part->Energy();
1820 Double_t P= part->P();
1821 M= TMath::Sqrt((fEnergyCMS-E-P)*(fEnergyCMS-E+P));
1824 Double_t Mmin, Mmax, wSD, wDD, wND;
1825 if(!GetWeightsDiffraction(M, Mmin, Mmax, wSD, wDD, wND)) return kFALSE;
1827 if(M>-1 && M<Mmin) return kFALSE;
1830 Int_t procType=fPythia->GetMSTI(1);
1832 if(procType== 94) proc0=1;
1833 if(procType== 92 || procType== 93) proc0=0;
1837 else if(proc0==1) proc=1;
1839 if(proc==1 && (AliPythiaRndm::GetPythiaRandom())->Uniform(0.,1.) > wDD) return kFALSE;
1840 if(proc==2 && (AliPythiaRndm::GetPythiaRandom())->Uniform(0.,1.) > wND) return kFALSE;
1843 // if(proc==1 || proc==2) return kFALSE;
1846 if(proc0!=0) fProcDiff = procType;
1847 else fProcDiff = 95;
1851 if(wSD<0) AliError("wSD<0 ! \n");
1853 if((AliPythiaRndm::GetPythiaRandom())->Uniform(0.,1.)> wSD) return kFALSE;
1855 // printf("iPart = %d\n", iPart);
1857 if(iPart==iPart1) fProcDiff=93;
1858 else if(iPart==iPart2) fProcDiff=92;
1860 printf("EROOR: iPart!=iPart1 && iPart!=iPart2\n");
1869 Bool_t AliGenPythia::GetWeightsDiffraction(Double_t M, Double_t &Mmin, Double_t &Mmax,
1870 Double_t &wSD, Double_t &wDD, Double_t &wND)
1874 if(TMath::Abs(fEnergyCMS-900)<1 ){
1876 const Int_t nbin=400;
1878 1.080000, 1.577300, 2.074600, 2.571900, 3.069200, 3.566500,
1879 4.063800, 4.561100, 5.058400, 5.555700, 6.053000, 6.550300,
1880 7.047600, 7.544900, 8.042200, 8.539500, 9.036800, 9.534100,
1881 10.031400, 10.528700, 11.026000, 11.523300, 12.020600, 12.517900,
1882 13.015200, 13.512500, 14.009800, 14.507100, 15.004400, 15.501700,
1883 15.999000, 16.496300, 16.993600, 17.490900, 17.988200, 18.485500,
1884 18.982800, 19.480100, 19.977400, 20.474700, 20.972000, 21.469300,
1885 21.966600, 22.463900, 22.961200, 23.458500, 23.955800, 24.453100,
1886 24.950400, 25.447700, 25.945000, 26.442300, 26.939600, 27.436900,
1887 27.934200, 28.431500, 28.928800, 29.426100, 29.923400, 30.420700,
1888 30.918000, 31.415300, 31.912600, 32.409900, 32.907200, 33.404500,
1889 33.901800, 34.399100, 34.896400, 35.393700, 35.891000, 36.388300,
1890 36.885600, 37.382900, 37.880200, 38.377500, 38.874800, 39.372100,
1891 39.869400, 40.366700, 40.864000, 41.361300, 41.858600, 42.355900,
1892 42.853200, 43.350500, 43.847800, 44.345100, 44.842400, 45.339700,
1893 45.837000, 46.334300, 46.831600, 47.328900, 47.826200, 48.323500,
1894 48.820800, 49.318100, 49.815400, 50.312700, 50.810000, 51.307300,
1895 51.804600, 52.301900, 52.799200, 53.296500, 53.793800, 54.291100,
1896 54.788400, 55.285700, 55.783000, 56.280300, 56.777600, 57.274900,
1897 57.772200, 58.269500, 58.766800, 59.264100, 59.761400, 60.258700,
1898 60.756000, 61.253300, 61.750600, 62.247900, 62.745200, 63.242500,
1899 63.739800, 64.237100, 64.734400, 65.231700, 65.729000, 66.226300,
1900 66.723600, 67.220900, 67.718200, 68.215500, 68.712800, 69.210100,
1901 69.707400, 70.204700, 70.702000, 71.199300, 71.696600, 72.193900,
1902 72.691200, 73.188500, 73.685800, 74.183100, 74.680400, 75.177700,
1903 75.675000, 76.172300, 76.669600, 77.166900, 77.664200, 78.161500,
1904 78.658800, 79.156100, 79.653400, 80.150700, 80.648000, 81.145300,
1905 81.642600, 82.139900, 82.637200, 83.134500, 83.631800, 84.129100,
1906 84.626400, 85.123700, 85.621000, 86.118300, 86.615600, 87.112900,
1907 87.610200, 88.107500, 88.604800, 89.102100, 89.599400, 90.096700,
1908 90.594000, 91.091300, 91.588600, 92.085900, 92.583200, 93.080500,
1909 93.577800, 94.075100, 94.572400, 95.069700, 95.567000, 96.064300,
1910 96.561600, 97.058900, 97.556200, 98.053500, 98.550800, 99.048100,
1911 99.545400, 100.042700, 100.540000, 101.037300, 101.534600, 102.031900,
1912 102.529200, 103.026500, 103.523800, 104.021100, 104.518400, 105.015700,
1913 105.513000, 106.010300, 106.507600, 107.004900, 107.502200, 107.999500,
1914 108.496800, 108.994100, 109.491400, 109.988700, 110.486000, 110.983300,
1915 111.480600, 111.977900, 112.475200, 112.972500, 113.469800, 113.967100,
1916 114.464400, 114.961700, 115.459000, 115.956300, 116.453600, 116.950900,
1917 117.448200, 117.945500, 118.442800, 118.940100, 119.437400, 119.934700,
1918 120.432000, 120.929300, 121.426600, 121.923900, 122.421200, 122.918500,
1919 123.415800, 123.913100, 124.410400, 124.907700, 125.405000, 125.902300,
1920 126.399600, 126.896900, 127.394200, 127.891500, 128.388800, 128.886100,
1921 129.383400, 129.880700, 130.378000, 130.875300, 131.372600, 131.869900,
1922 132.367200, 132.864500, 133.361800, 133.859100, 134.356400, 134.853700,
1923 135.351000, 135.848300, 136.345600, 136.842900, 137.340200, 137.837500,
1924 138.334800, 138.832100, 139.329400, 139.826700, 140.324000, 140.821300,
1925 141.318600, 141.815900, 142.313200, 142.810500, 143.307800, 143.805100,
1926 144.302400, 144.799700, 145.297000, 145.794300, 146.291600, 146.788900,
1927 147.286200, 147.783500, 148.280800, 148.778100, 149.275400, 149.772700,
1928 150.270000, 150.767300, 151.264600, 151.761900, 152.259200, 152.756500,
1929 153.253800, 153.751100, 154.248400, 154.745700, 155.243000, 155.740300,
1930 156.237600, 156.734900, 157.232200, 157.729500, 158.226800, 158.724100,
1931 159.221400, 159.718700, 160.216000, 160.713300, 161.210600, 161.707900,
1932 162.205200, 162.702500, 163.199800, 163.697100, 164.194400, 164.691700,
1933 165.189000, 165.686300, 166.183600, 166.680900, 167.178200, 167.675500,
1934 168.172800, 168.670100, 169.167400, 169.664700, 170.162000, 170.659300,
1935 171.156600, 171.653900, 172.151200, 172.648500, 173.145800, 173.643100,
1936 174.140400, 174.637700, 175.135000, 175.632300, 176.129600, 176.626900,
1937 177.124200, 177.621500, 178.118800, 178.616100, 179.113400, 179.610700,
1938 180.108000, 180.605300, 181.102600, 181.599900, 182.097200, 182.594500,
1939 183.091800, 183.589100, 184.086400, 184.583700, 185.081000, 185.578300,
1940 186.075600, 186.572900, 187.070200, 187.567500, 188.064800, 188.562100,
1941 189.059400, 189.556700, 190.054000, 190.551300, 191.048600, 191.545900,
1942 192.043200, 192.540500, 193.037800, 193.535100, 194.032400, 194.529700,
1943 195.027000, 195.524300, 196.021600, 196.518900, 197.016200, 197.513500,
1944 198.010800, 198.508100, 199.005400, 199.502700, 200.000000};
1946 1.000000, 0.643457, 0.645609, 0.648347, 0.604563, 0.605002,
1947 0.602819, 0.611473, 0.576412, 0.562354, 0.550216, 0.529285,
1948 0.534558, 0.534364, 0.530358, 0.518475, 0.489253, 0.469754,
1949 0.469825, 0.450513, 0.455849, 0.435312, 0.437210, 0.456686,
1950 0.413577, 0.427093, 0.426894, 0.418834, 0.409475, 0.388483,
1951 0.412208, 0.388912, 0.389611, 0.382572, 0.389220, 0.370964,
1952 0.380463, 0.370873, 0.363701, 0.369363, 0.357361, 0.365759,
1953 0.348566, 0.337062, 0.348190, 0.332330, 0.359001, 0.335836,
1954 0.339154, 0.335599, 0.336035, 0.335204, 0.353440, 0.337836,
1955 0.333874, 0.307120, 0.294963, 0.324978, 0.313359, 0.317688,
1956 0.323758, 0.319304, 0.335317, 0.301765, 0.317257, 0.356331,
1957 0.323119, 0.297732, 0.303188, 0.326102, 0.316467, 0.294728,
1958 0.308135, 0.288857, 0.325692, 0.312493, 0.291100, 0.325921,
1959 0.313317, 0.295980, 0.308481, 0.328380, 0.313081, 0.296763,
1960 0.295431, 0.317325, 0.320462, 0.286918, 0.316035, 0.335208,
1961 0.283193, 0.333945, 0.292534, 0.294164, 0.330680, 0.296992,
1962 0.285509, 0.317260, 0.298269, 0.311299, 0.312129, 0.286822,
1963 0.287442, 0.319139, 0.283314, 0.318454, 0.297727, 0.301597,
1964 0.282483, 0.294792, 0.305569, 0.290957, 0.297817, 0.282908,
1965 0.272401, 0.305584, 0.300220, 0.297020, 0.298781, 0.278008,
1966 0.277727, 0.323777, 0.287419, 0.342074, 0.287259, 0.303658,
1967 0.302668, 0.279622, 0.280586, 0.313630, 0.276068, 0.257051,
1968 0.309996, 0.265534, 0.297138, 0.281738, 0.294610, 0.292882,
1969 0.286860, 0.312686, 0.293244, 0.293744, 0.271375, 0.278734,
1970 0.280308, 0.304739, 0.287907, 0.285261, 0.311180, 0.313476,
1971 0.289660, 0.289282, 0.319505, 0.271285, 0.272008, 0.289245,
1972 0.281038, 0.285284, 0.295836, 0.281416, 0.283501, 0.295417,
1973 0.304372, 0.297764, 0.291378, 0.321530, 0.315604, 0.329507,
1974 0.282609, 0.275576, 0.283721, 0.311714, 0.283760, 0.273188,
1975 0.312193, 0.264347, 0.281532, 0.301226, 0.281718, 0.336408,
1976 0.283157, 0.332010, 0.289974, 0.290256, 0.301569, 0.332228,
1977 0.288282, 0.326339, 0.313653, 0.300361, 0.289470, 0.264830,
1978 0.298659, 0.272359, 0.278878, 0.306001, 0.328168, 0.294991,
1979 0.327737, 0.278056, 0.302435, 0.284183, 0.279270, 0.307279,
1980 0.307917, 0.315196, 0.283803, 0.313333, 0.315730, 0.304818,
1981 0.307171, 0.295223, 0.333741, 0.346911, 0.310143, 0.336686,
1982 0.275459, 0.334781, 0.295405, 0.275816, 0.301678, 0.327242,
1983 0.320717, 0.309230, 0.292145, 0.294489, 0.305088, 0.300969,
1984 0.277438, 0.326159, 0.297065, 0.301177, 0.303843, 0.275382,
1985 0.304019, 0.284166, 0.289610, 0.331611, 0.317131, 0.310880,
1986 0.360456, 0.294052, 0.342694, 0.327166, 0.336797, 0.298040,
1987 0.295767, 0.260053, 0.325544, 0.335310, 0.320182, 0.301072,
1988 0.313117, 0.283407, 0.299206, 0.293525, 0.305067, 0.255978,
1989 0.327055, 0.316382, 0.317700, 0.278993, 0.283120, 0.314000,
1990 0.274396, 0.291208, 0.348813, 0.319603, 0.313076, 0.289155,
1991 0.343988, 0.311426, 0.322896, 0.328726, 0.337062, 0.389093,
1992 0.284122, 0.312184, 0.304008, 0.319170, 0.320778, 0.288390,
1993 0.337272, 0.356273, 0.343310, 0.312209, 0.330709, 0.297977,
1994 0.346146, 0.369162, 0.324385, 0.339831, 0.337037, 0.318739,
1995 0.343157, 0.277720, 0.368407, 0.321330, 0.338997, 0.314220,
1996 0.328861, 0.321824, 0.328013, 0.356925, 0.359144, 0.296314,
1997 0.345415, 0.396711, 0.347032, 0.294928, 0.343799, 0.322331,
1998 0.328656, 0.326098, 0.337338, 0.337038, 0.300179, 0.351391,
1999 0.324337, 0.330896, 0.302842, 0.310522, 0.337052, 0.359989,
2000 0.383250, 0.359355, 0.315382, 0.333113, 0.342598, 0.355348,
2001 0.320751, 0.320475, 0.351762, 0.351475, 0.338358, 0.326153,
2002 0.302507, 0.340048, 0.318685, 0.381646, 0.339320, 0.299453,
2003 0.426599, 0.393515, 0.353929, 0.328435, 0.413976, 0.292558,
2004 0.379340, 0.358344, 0.409259, 0.313821, 0.336675, 0.324521,
2005 0.408382, 0.346273, 0.312939, 0.362453, 0.343152, 0.330577,
2006 0.332831, 0.353299, 0.347745, 0.334818, 0.332234, 0.385585,
2007 0.395483, 0.395316, 0.326972, 0.349434, 0.317604, 0.328980,
2008 0.375056, 0.317290, 0.357083, 0.346165, 0.310444, 0.356873,
2009 0.359523, 0.365308, 0.365122, 0.383685, 0.370975, 0.396928,
2010 0.407654, 0.307755, 0.323033, 0.350580, 0.345231, 0.342462,
2011 0.400000, 0.318309, 0.403570, 0.322856, 0.383053, 0.422252,
2012 0.386112, 0.364314, 0.434375, 0.334629};
2019 if(M<Mmin || M>Mmax) return kTRUE;
2022 for(Int_t i=1; i<=nbin; i++)
2025 // printf("Mi> %f && Mi< %f\n", bin[i-1], bin[i]);
2031 else if(TMath::Abs(fEnergyCMS-2760)<1 ){
2033 const Int_t nbin=400;
2035 1.080000, 1.577300, 2.074600, 2.571900, 3.069200, 3.566500,
2036 4.063800, 4.561100, 5.058400, 5.555700, 6.053000, 6.550300,
2037 7.047600, 7.544900, 8.042200, 8.539500, 9.036800, 9.534100,
2038 10.031400, 10.528700, 11.026000, 11.523300, 12.020600, 12.517900,
2039 13.015200, 13.512500, 14.009800, 14.507100, 15.004400, 15.501700,
2040 15.999000, 16.496300, 16.993600, 17.490900, 17.988200, 18.485500,
2041 18.982800, 19.480100, 19.977400, 20.474700, 20.972000, 21.469300,
2042 21.966600, 22.463900, 22.961200, 23.458500, 23.955800, 24.453100,
2043 24.950400, 25.447700, 25.945000, 26.442300, 26.939600, 27.436900,
2044 27.934200, 28.431500, 28.928800, 29.426100, 29.923400, 30.420700,
2045 30.918000, 31.415300, 31.912600, 32.409900, 32.907200, 33.404500,
2046 33.901800, 34.399100, 34.896400, 35.393700, 35.891000, 36.388300,
2047 36.885600, 37.382900, 37.880200, 38.377500, 38.874800, 39.372100,
2048 39.869400, 40.366700, 40.864000, 41.361300, 41.858600, 42.355900,
2049 42.853200, 43.350500, 43.847800, 44.345100, 44.842400, 45.339700,
2050 45.837000, 46.334300, 46.831600, 47.328900, 47.826200, 48.323500,
2051 48.820800, 49.318100, 49.815400, 50.312700, 50.810000, 51.307300,
2052 51.804600, 52.301900, 52.799200, 53.296500, 53.793800, 54.291100,
2053 54.788400, 55.285700, 55.783000, 56.280300, 56.777600, 57.274900,
2054 57.772200, 58.269500, 58.766800, 59.264100, 59.761400, 60.258700,
2055 60.756000, 61.253300, 61.750600, 62.247900, 62.745200, 63.242500,
2056 63.739800, 64.237100, 64.734400, 65.231700, 65.729000, 66.226300,
2057 66.723600, 67.220900, 67.718200, 68.215500, 68.712800, 69.210100,
2058 69.707400, 70.204700, 70.702000, 71.199300, 71.696600, 72.193900,
2059 72.691200, 73.188500, 73.685800, 74.183100, 74.680400, 75.177700,
2060 75.675000, 76.172300, 76.669600, 77.166900, 77.664200, 78.161500,
2061 78.658800, 79.156100, 79.653400, 80.150700, 80.648000, 81.145300,
2062 81.642600, 82.139900, 82.637200, 83.134500, 83.631800, 84.129100,
2063 84.626400, 85.123700, 85.621000, 86.118300, 86.615600, 87.112900,
2064 87.610200, 88.107500, 88.604800, 89.102100, 89.599400, 90.096700,
2065 90.594000, 91.091300, 91.588600, 92.085900, 92.583200, 93.080500,
2066 93.577800, 94.075100, 94.572400, 95.069700, 95.567000, 96.064300,
2067 96.561600, 97.058900, 97.556200, 98.053500, 98.550800, 99.048100,
2068 99.545400, 100.042700, 100.540000, 101.037300, 101.534600, 102.031900,
2069 102.529200, 103.026500, 103.523800, 104.021100, 104.518400, 105.015700,
2070 105.513000, 106.010300, 106.507600, 107.004900, 107.502200, 107.999500,
2071 108.496800, 108.994100, 109.491400, 109.988700, 110.486000, 110.983300,
2072 111.480600, 111.977900, 112.475200, 112.972500, 113.469800, 113.967100,
2073 114.464400, 114.961700, 115.459000, 115.956300, 116.453600, 116.950900,
2074 117.448200, 117.945500, 118.442800, 118.940100, 119.437400, 119.934700,
2075 120.432000, 120.929300, 121.426600, 121.923900, 122.421200, 122.918500,
2076 123.415800, 123.913100, 124.410400, 124.907700, 125.405000, 125.902300,
2077 126.399600, 126.896900, 127.394200, 127.891500, 128.388800, 128.886100,
2078 129.383400, 129.880700, 130.378000, 130.875300, 131.372600, 131.869900,
2079 132.367200, 132.864500, 133.361800, 133.859100, 134.356400, 134.853700,
2080 135.351000, 135.848300, 136.345600, 136.842900, 137.340200, 137.837500,
2081 138.334800, 138.832100, 139.329400, 139.826700, 140.324000, 140.821300,
2082 141.318600, 141.815900, 142.313200, 142.810500, 143.307800, 143.805100,
2083 144.302400, 144.799700, 145.297000, 145.794300, 146.291600, 146.788900,
2084 147.286200, 147.783500, 148.280800, 148.778100, 149.275400, 149.772700,
2085 150.270000, 150.767300, 151.264600, 151.761900, 152.259200, 152.756500,
2086 153.253800, 153.751100, 154.248400, 154.745700, 155.243000, 155.740300,
2087 156.237600, 156.734900, 157.232200, 157.729500, 158.226800, 158.724100,
2088 159.221400, 159.718700, 160.216000, 160.713300, 161.210600, 161.707900,
2089 162.205200, 162.702500, 163.199800, 163.697100, 164.194400, 164.691700,
2090 165.189000, 165.686300, 166.183600, 166.680900, 167.178200, 167.675500,
2091 168.172800, 168.670100, 169.167400, 169.664700, 170.162000, 170.659300,
2092 171.156600, 171.653900, 172.151200, 172.648500, 173.145800, 173.643100,
2093 174.140400, 174.637700, 175.135000, 175.632300, 176.129600, 176.626900,
2094 177.124200, 177.621500, 178.118800, 178.616100, 179.113400, 179.610700,
2095 180.108000, 180.605300, 181.102600, 181.599900, 182.097200, 182.594500,
2096 183.091800, 183.589100, 184.086400, 184.583700, 185.081000, 185.578300,
2097 186.075600, 186.572900, 187.070200, 187.567500, 188.064800, 188.562100,
2098 189.059400, 189.556700, 190.054000, 190.551300, 191.048600, 191.545900,
2099 192.043200, 192.540500, 193.037800, 193.535100, 194.032400, 194.529700,
2100 195.027000, 195.524300, 196.021600, 196.518900, 197.016200, 197.513500,
2101 198.010800, 198.508100, 199.005400, 199.502700, 200.000000};
2103 1.000000, 0.692593, 0.673384, 0.666273, 0.657285, 0.637723,
2104 0.625881, 0.643590, 0.606100, 0.589007, 0.567824, 0.578705,
2105 0.538530, 0.517937, 0.528278, 0.515918, 0.539461, 0.466186,
2106 0.489869, 0.468402, 0.465017, 0.453336, 0.460769, 0.474638,
2107 0.456347, 0.434471, 0.427478, 0.435435, 0.410934, 0.366431,
2108 0.382240, 0.379513, 0.394249, 0.386837, 0.353103, 0.382138,
2109 0.377497, 0.389479, 0.378736, 0.347933, 0.354605, 0.352077,
2110 0.324443, 0.358792, 0.339968, 0.359052, 0.330734, 0.318291,
2111 0.333703, 0.358644, 0.335819, 0.332213, 0.309051, 0.309975,
2112 0.331626, 0.304407, 0.309819, 0.312097, 0.312462, 0.320411,
2113 0.280401, 0.302311, 0.315863, 0.281479, 0.310003, 0.296911,
2114 0.313676, 0.281071, 0.294163, 0.306500, 0.283462, 0.274867,
2115 0.307149, 0.270555, 0.282264, 0.287373, 0.307849, 0.278675,
2116 0.286990, 0.278269, 0.300105, 0.286799, 0.265674, 0.275140,
2117 0.285702, 0.257352, 0.267714, 0.248204, 0.252220, 0.255678,
2118 0.282946, 0.268257, 0.282375, 0.262675, 0.275564, 0.248345,
2119 0.236259, 0.291914, 0.259936, 0.241338, 0.267389, 0.285044,
2120 0.289419, 0.253594, 0.284568, 0.231840, 0.260008, 0.268527,
2121 0.275363, 0.224115, 0.281260, 0.257913, 0.295152, 0.264399,
2122 0.232287, 0.282533, 0.223431, 0.255756, 0.244471, 0.221695,
2123 0.272450, 0.284244, 0.253682, 0.270717, 0.275403, 0.240323,
2124 0.245081, 0.241859, 0.216340, 0.244789, 0.220291, 0.238478,
2125 0.224691, 0.244058, 0.266117, 0.271478, 0.242012, 0.267321,
2126 0.248494, 0.253343, 0.255606, 0.235458, 0.241079, 0.233223,
2127 0.226813, 0.259224, 0.234239, 0.258606, 0.210892, 0.238186,
2128 0.243271, 0.222678, 0.213437, 0.273939, 0.247966, 0.232548,
2129 0.263438, 0.222089, 0.272111, 0.248818, 0.244295, 0.238368,
2130 0.236908, 0.248776, 0.232604, 0.231194, 0.227117, 0.231152,
2131 0.282140, 0.229778, 0.232631, 0.261794, 0.216633, 0.253471,
2132 0.242157, 0.227406, 0.269335, 0.230547, 0.210618, 0.285872,
2133 0.248776, 0.229875, 0.242728, 0.227388, 0.220567, 0.222062,
2134 0.235950, 0.224087, 0.228895, 0.208287, 0.235999, 0.208696,
2135 0.230367, 0.267667, 0.220484, 0.233402, 0.233815, 0.250455,
2136 0.253120, 0.219556, 0.230980, 0.236661, 0.222395, 0.226111,
2137 0.198315, 0.210555, 0.202813, 0.208594, 0.235976, 0.221490,
2138 0.243059, 0.204901, 0.216987, 0.229039, 0.231466, 0.221975,
2139 0.231220, 0.253638, 0.250448, 0.260291, 0.328345, 0.205739,
2140 0.222014, 0.251513, 0.279427, 0.270506, 0.248409, 0.222472,
2141 0.291632, 0.227796, 0.248769, 0.276896, 0.214742, 0.200139,
2142 0.230693, 0.226031, 0.268900, 0.185160, 0.245353, 0.205843,
2143 0.231155, 0.219122, 0.214811, 0.199763, 0.274179, 0.217598,
2144 0.274988, 0.237244, 0.211820, 0.225459, 0.252799, 0.235948,
2145 0.224986, 0.245385, 0.237770, 0.213373, 0.229737, 0.215487,
2146 0.234453, 0.249684, 0.239435, 0.250422, 0.257194, 0.214762,
2147 0.212266, 0.228988, 0.253798, 0.201607, 0.239946, 0.205245,
2148 0.231670, 0.212774, 0.206768, 0.231563, 0.189388, 0.227926,
2149 0.227203, 0.237754, 0.221628, 0.211138, 0.203322, 0.200985,
2150 0.231780, 0.220294, 0.232686, 0.234243, 0.218264, 0.255870,
2151 0.223213, 0.238670, 0.213713, 0.213064, 0.246700, 0.233446,
2152 0.221503, 0.206767, 0.200722, 0.226179, 0.237425, 0.239229,
2153 0.238611, 0.240419, 0.247806, 0.215923, 0.205298, 0.232778,
2154 0.272312, 0.226773, 0.258103, 0.223287, 0.269404, 0.203398,
2155 0.223782, 0.204213, 0.229664, 0.234040, 0.228419, 0.203936,
2156 0.263686, 0.199141, 0.236127, 0.214058, 0.204611, 0.224324,
2157 0.292140, 0.190735, 0.235157, 0.213018, 0.257085, 0.190554,
2158 0.203197, 0.213044, 0.237023, 0.214243, 0.193562, 0.262403,
2159 0.206256, 0.221396, 0.233588, 0.256611, 0.249731, 0.226683,
2160 0.199330, 0.251026, 0.222596, 0.201941, 0.186374, 0.221038,
2161 0.196555, 0.222560, 0.299419, 0.231979, 0.242924, 0.198310,
2162 0.217628, 0.235458, 0.278595, 0.218624, 0.277305, 0.239109,
2163 0.205600, 0.253715, 0.221173, 0.218195, 0.277647, 0.241974,
2164 0.268748, 0.268128, 0.292258, 0.249420, 0.191034, 0.219506,
2165 0.216502, 0.250677, 0.193386, 0.201310, 0.259464, 0.255351,
2166 0.269628, 0.221063, 0.294079, 0.196726, 0.233634, 0.221870,
2167 0.216236, 0.197259, 0.247433, 0.272765, 0.294079, 0.236336,
2168 0.206396, 0.238524, 0.247846, 0.269519, 0.237141, 0.230611,
2169 0.201712, 0.242225, 0.255565, 0.258738};
2176 if(M<Mmin || M>Mmax) return kTRUE;
2179 for(Int_t i=1; i<=nbin; i++)
2182 // printf("Mi> %f && Mi< %f\n", bin[i-1], bin[i]);
2190 else if(TMath::Abs(fEnergyCMS-7000)<1 ){
2191 const Int_t nbin=400;
2193 1.080000, 1.577300, 2.074600, 2.571900, 3.069200, 3.566500,
2194 4.063800, 4.561100, 5.058400, 5.555700, 6.053000, 6.550300,
2195 7.047600, 7.544900, 8.042200, 8.539500, 9.036800, 9.534100,
2196 10.031400, 10.528700, 11.026000, 11.523300, 12.020600, 12.517900,
2197 13.015200, 13.512500, 14.009800, 14.507100, 15.004400, 15.501700,
2198 15.999000, 16.496300, 16.993600, 17.490900, 17.988200, 18.485500,
2199 18.982800, 19.480100, 19.977400, 20.474700, 20.972000, 21.469300,
2200 21.966600, 22.463900, 22.961200, 23.458500, 23.955800, 24.453100,
2201 24.950400, 25.447700, 25.945000, 26.442300, 26.939600, 27.436900,
2202 27.934200, 28.431500, 28.928800, 29.426100, 29.923400, 30.420700,
2203 30.918000, 31.415300, 31.912600, 32.409900, 32.907200, 33.404500,
2204 33.901800, 34.399100, 34.896400, 35.393700, 35.891000, 36.388300,
2205 36.885600, 37.382900, 37.880200, 38.377500, 38.874800, 39.372100,
2206 39.869400, 40.366700, 40.864000, 41.361300, 41.858600, 42.355900,
2207 42.853200, 43.350500, 43.847800, 44.345100, 44.842400, 45.339700,
2208 45.837000, 46.334300, 46.831600, 47.328900, 47.826200, 48.323500,
2209 48.820800, 49.318100, 49.815400, 50.312700, 50.810000, 51.307300,
2210 51.804600, 52.301900, 52.799200, 53.296500, 53.793800, 54.291100,
2211 54.788400, 55.285700, 55.783000, 56.280300, 56.777600, 57.274900,
2212 57.772200, 58.269500, 58.766800, 59.264100, 59.761400, 60.258700,
2213 60.756000, 61.253300, 61.750600, 62.247900, 62.745200, 63.242500,
2214 63.739800, 64.237100, 64.734400, 65.231700, 65.729000, 66.226300,
2215 66.723600, 67.220900, 67.718200, 68.215500, 68.712800, 69.210100,
2216 69.707400, 70.204700, 70.702000, 71.199300, 71.696600, 72.193900,
2217 72.691200, 73.188500, 73.685800, 74.183100, 74.680400, 75.177700,
2218 75.675000, 76.172300, 76.669600, 77.166900, 77.664200, 78.161500,
2219 78.658800, 79.156100, 79.653400, 80.150700, 80.648000, 81.145300,
2220 81.642600, 82.139900, 82.637200, 83.134500, 83.631800, 84.129100,
2221 84.626400, 85.123700, 85.621000, 86.118300, 86.615600, 87.112900,
2222 87.610200, 88.107500, 88.604800, 89.102100, 89.599400, 90.096700,
2223 90.594000, 91.091300, 91.588600, 92.085900, 92.583200, 93.080500,
2224 93.577800, 94.075100, 94.572400, 95.069700, 95.567000, 96.064300,
2225 96.561600, 97.058900, 97.556200, 98.053500, 98.550800, 99.048100,
2226 99.545400, 100.042700, 100.540000, 101.037300, 101.534600, 102.031900,
2227 102.529200, 103.026500, 103.523800, 104.021100, 104.518400, 105.015700,
2228 105.513000, 106.010300, 106.507600, 107.004900, 107.502200, 107.999500,
2229 108.496800, 108.994100, 109.491400, 109.988700, 110.486000, 110.983300,
2230 111.480600, 111.977900, 112.475200, 112.972500, 113.469800, 113.967100,
2231 114.464400, 114.961700, 115.459000, 115.956300, 116.453600, 116.950900,
2232 117.448200, 117.945500, 118.442800, 118.940100, 119.437400, 119.934700,
2233 120.432000, 120.929300, 121.426600, 121.923900, 122.421200, 122.918500,
2234 123.415800, 123.913100, 124.410400, 124.907700, 125.405000, 125.902300,
2235 126.399600, 126.896900, 127.394200, 127.891500, 128.388800, 128.886100,
2236 129.383400, 129.880700, 130.378000, 130.875300, 131.372600, 131.869900,
2237 132.367200, 132.864500, 133.361800, 133.859100, 134.356400, 134.853700,
2238 135.351000, 135.848300, 136.345600, 136.842900, 137.340200, 137.837500,
2239 138.334800, 138.832100, 139.329400, 139.826700, 140.324000, 140.821300,
2240 141.318600, 141.815900, 142.313200, 142.810500, 143.307800, 143.805100,
2241 144.302400, 144.799700, 145.297000, 145.794300, 146.291600, 146.788900,
2242 147.286200, 147.783500, 148.280800, 148.778100, 149.275400, 149.772700,
2243 150.270000, 150.767300, 151.264600, 151.761900, 152.259200, 152.756500,
2244 153.253800, 153.751100, 154.248400, 154.745700, 155.243000, 155.740300,
2245 156.237600, 156.734900, 157.232200, 157.729500, 158.226800, 158.724100,
2246 159.221400, 159.718700, 160.216000, 160.713300, 161.210600, 161.707900,
2247 162.205200, 162.702500, 163.199800, 163.697100, 164.194400, 164.691700,
2248 165.189000, 165.686300, 166.183600, 166.680900, 167.178200, 167.675500,
2249 168.172800, 168.670100, 169.167400, 169.664700, 170.162000, 170.659300,
2250 171.156600, 171.653900, 172.151200, 172.648500, 173.145800, 173.643100,
2251 174.140400, 174.637700, 175.135000, 175.632300, 176.129600, 176.626900,
2252 177.124200, 177.621500, 178.118800, 178.616100, 179.113400, 179.610700,
2253 180.108000, 180.605300, 181.102600, 181.599900, 182.097200, 182.594500,
2254 183.091800, 183.589100, 184.086400, 184.583700, 185.081000, 185.578300,
2255 186.075600, 186.572900, 187.070200, 187.567500, 188.064800, 188.562100,
2256 189.059400, 189.556700, 190.054000, 190.551300, 191.048600, 191.545900,
2257 192.043200, 192.540500, 193.037800, 193.535100, 194.032400, 194.529700,
2258 195.027000, 195.524300, 196.021600, 196.518900, 197.016200, 197.513500,
2259 198.010800, 198.508100, 199.005400, 199.502700, 200.000000};
2261 1.000000, 0.798640, 0.770197, 0.717990, 0.699434, 0.692093,
2262 0.620940, 0.655294, 0.636496, 0.633483, 0.600518, 0.594548,
2263 0.550759, 0.550413, 0.528926, 0.525877, 0.506701, 0.504117,
2264 0.486973, 0.492930, 0.461833, 0.488695, 0.438874, 0.460274,
2265 0.451551, 0.454424, 0.449270, 0.387571, 0.427556, 0.428740,
2266 0.401682, 0.402260, 0.408961, 0.395417, 0.387426, 0.378602,
2267 0.357826, 0.359125, 0.348494, 0.363710, 0.356117, 0.340363,
2268 0.337637, 0.396084, 0.341777, 0.340551, 0.348838, 0.344014,
2269 0.340468, 0.317654, 0.355584, 0.326023, 0.373416, 0.312298,
2270 0.326522, 0.290540, 0.335557, 0.318689, 0.327544, 0.319501,
2271 0.331754, 0.312728, 0.282263, 0.274937, 0.303867, 0.307820,
2272 0.289344, 0.268934, 0.288908, 0.290018, 0.291369, 0.295242,
2273 0.289067, 0.277685, 0.267957, 0.267559, 0.320229, 0.265060,
2274 0.305931, 0.305352, 0.262064, 0.281879, 0.287780, 0.270033,
2275 0.270814, 0.276667, 0.271531, 0.275283, 0.258189, 0.287969,
2276 0.251247, 0.301527, 0.267230, 0.245860, 0.293125, 0.253421,
2277 0.272396, 0.243637, 0.236206, 0.278452, 0.246544, 0.263165,
2278 0.230484, 0.231102, 0.258281, 0.244707, 0.270111, 0.248295,
2279 0.246942, 0.245592, 0.272766, 0.254243, 0.199647, 0.262590,
2280 0.226710, 0.243836, 0.214153, 0.233535, 0.235207, 0.234245,
2281 0.247698, 0.248379, 0.241463, 0.265230, 0.223242, 0.236191,
2282 0.252700, 0.231865, 0.228352, 0.229101, 0.237385, 0.246485,
2283 0.254706, 0.245565, 0.224841, 0.257301, 0.240968, 0.202379,
2284 0.236943, 0.241683, 0.220809, 0.219014, 0.213015, 0.244470,
2285 0.221042, 0.198996, 0.225295, 0.264466, 0.200600, 0.228143,
2286 0.250138, 0.216784, 0.268002, 0.275734, 0.218144, 0.229866,
2287 0.235443, 0.208909, 0.215067, 0.189743, 0.216741, 0.242686,
2288 0.200234, 0.218882, 0.245991, 0.222815, 0.219576, 0.209773,
2289 0.205247, 0.203855, 0.202534, 0.192536, 0.223447, 0.225810,
2290 0.220583, 0.248421, 0.223424, 0.206033, 0.203130, 0.221225,
2291 0.223763, 0.239216, 0.252311, 0.206893, 0.228461, 0.233591,
2292 0.201512, 0.242382, 0.215147, 0.232578, 0.207116, 0.239213,
2293 0.196215, 0.184367, 0.235135, 0.189768, 0.274084, 0.206267,
2294 0.203283, 0.248828, 0.250214, 0.217465, 0.232080, 0.215150,
2295 0.207936, 0.176789, 0.191338, 0.188655, 0.196181, 0.223691,
2296 0.254257, 0.216874, 0.213536, 0.221399, 0.192024, 0.212534,
2297 0.218169, 0.226635, 0.201191, 0.212700, 0.211634, 0.232353,
2298 0.223636, 0.208605, 0.249132, 0.183681, 0.221842, 0.187136,
2299 0.203772, 0.249575, 0.217713, 0.205193, 0.193941, 0.219707,
2300 0.264226, 0.182105, 0.207183, 0.220845, 0.236571, 0.182390,
2301 0.212914, 0.186266, 0.195361, 0.217665, 0.204527, 0.188804,
2302 0.222832, 0.191193, 0.201073, 0.185616, 0.195011, 0.200828,
2303 0.181434, 0.233707, 0.202925, 0.211992, 0.173100, 0.205258,
2304 0.182695, 0.216520, 0.218422, 0.209490, 0.211257, 0.215801,
2305 0.220652, 0.211409, 0.195731, 0.194957, 0.198888, 0.234926,
2306 0.221377, 0.229822, 0.176700, 0.172322, 0.212265, 0.206133,
2307 0.170355, 0.253305, 0.198688, 0.240043, 0.225384, 0.174729,
2308 0.195820, 0.200093, 0.196912, 0.212308, 0.200490, 0.240415,
2309 0.159744, 0.173686, 0.223853, 0.213604, 0.193779, 0.179328,
2310 0.180873, 0.237481, 0.179640, 0.235857, 0.202847, 0.167690,
2311 0.245093, 0.215504, 0.205848, 0.184408, 0.201626, 0.209651,
2312 0.236306, 0.217803, 0.188534, 0.187861, 0.161663, 0.221718,
2313 0.152044, 0.190412, 0.173505, 0.208995, 0.174575, 0.180734,
2314 0.247704, 0.203388, 0.194184, 0.169679, 0.182703, 0.213402,
2315 0.191808, 0.178604, 0.178116, 0.198452, 0.174687, 0.169809,
2316 0.222851, 0.156811, 0.229170, 0.210359, 0.178557, 0.248570,
2317 0.208536, 0.192571, 0.178912, 0.224505, 0.170822, 0.205492,
2318 0.330973, 0.160924, 0.203724, 0.213255, 0.205827, 0.187162,
2319 0.181252, 0.191723, 0.238106, 0.182398, 0.202358, 0.212066,
2320 0.201255, 0.168159, 0.185219, 0.229176, 0.158222, 0.164092,
2321 0.215405, 0.200724, 0.234811, 0.184488, 0.213112, 0.198577,
2322 0.219622, 0.160671, 0.179349, 0.206681, 0.206091, 0.165618,
2323 0.165203, 0.174442, 0.179287, 0.187318, 0.163481, 0.217440,
2324 0.160381, 0.177025, 0.204385, 0.163676, 0.210733, 0.186519,
2325 0.230701, 0.191764, 0.185119, 0.190770, 0.242987, 0.186812,
2326 0.202906, 0.161935, 0.182426, 0.197922, 0.181475, 0.155903,
2327 0.175006, 0.223482, 0.202706, 0.218108};
2335 if(M<Mmin || M>Mmax) return kTRUE;
2338 for(Int_t i=1; i<=nbin; i++)
2341 // printf("Mi> %f && Mi< %f\n", bin[i-1], bin[i]);