X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PYTHIA6%2FAliGenPythia.cxx;h=175c7d89ac8718e76c7e1db5c1a37d46acd3a1ec;hb=e7c989e4b7fb4940e1b7175670324ca85b25f6fc;hp=af5924524d73975f3e96b1c66a860df4143004f1;hpb=4c4eac974221101f281302423f406706c4a3f3f6;p=u%2Fmrichter%2FAliRoot.git diff --git a/PYTHIA6/AliGenPythia.cxx b/PYTHIA6/AliGenPythia.cxx index af5924524d7..175c7d89ac8 100644 --- a/PYTHIA6/AliGenPythia.cxx +++ b/PYTHIA6/AliGenPythia.cxx @@ -51,7 +51,6 @@ AliGenPythia::AliGenPythia(): AliGenMC(), fProcess(kPyCharm), fStrucFunc(kCTEQ5L), - fEnergyCMS(5500.), fKineBias(0.), fTrials(0), fTrialsRun(0), @@ -108,6 +107,8 @@ AliGenPythia::AliGenPythia(): fHFoff(kFALSE), fTriggerParticle(0), fTriggerEta(0.9), + fTriggerMultiplicity(0), + fTriggerMultiplicityEta(0), fCountMode(kCountAll), fHeader(0), fRL(0), @@ -129,6 +130,7 @@ AliGenPythia::AliGenPythia(): { // Default Constructor + fEnergyCMS = 5500.; SetNuclei(0,0); if (!AliPythiaRndm::GetPythiaRandom()) AliPythiaRndm::SetPythiaRandom(GetRandom()); @@ -138,7 +140,6 @@ AliGenPythia::AliGenPythia(Int_t npart) :AliGenMC(npart), fProcess(kPyCharm), fStrucFunc(kCTEQ5L), - fEnergyCMS(5500.), fKineBias(0.), fTrials(0), fTrialsRun(0), @@ -195,6 +196,8 @@ AliGenPythia::AliGenPythia(Int_t npart) fHFoff(kFALSE), fTriggerParticle(0), fTriggerEta(0.9), + fTriggerMultiplicity(0), + fTriggerMultiplicityEta(0), fCountMode(kCountAll), fHeader(0), fRL(0), @@ -218,13 +221,13 @@ AliGenPythia::AliGenPythia(Int_t npart) // semimuonic decay // structure function GRVHO // + fEnergyCMS = 5500.; fName = "Pythia"; fTitle= "Particle Generator using PYTHIA"; SetForceDecay(); // Set random number generator if (!AliPythiaRndm::GetPythiaRandom()) AliPythiaRndm::SetPythiaRandom(GetRandom()); - fParticles = new TClonesArray("TParticle",1000); SetNuclei(0,0); } @@ -573,7 +576,7 @@ void AliGenPythia::Generate() // printf("Calling hadronisation %d\n", fPythia->GetN()); fPythia->Pyexec(); fTrials++; - fPythia->ImportParticles(fParticles,"All"); + fPythia->ImportParticles(&fParticles,"All"); Boost(); // // @@ -581,7 +584,7 @@ void AliGenPythia::Generate() Int_t i; fNprimaries = 0; - Int_t np = fParticles->GetEntriesFast(); + Int_t np = fParticles.GetEntriesFast(); if (np == 0) continue; // @@ -611,7 +614,7 @@ void AliGenPythia::Generate() fProcess != kPyBeautyppMNRwmi) { for (i = 0; i < np; i++) { - TParticle* iparticle = (TParticle *) fParticles->At(i); + TParticle* iparticle = (TParticle *) fParticles.At(i); Int_t ks = iparticle->GetStatusCode(); kf = CheckPDGCode(iparticle->GetPdgCode()); // No initial state partons @@ -638,20 +641,20 @@ void AliGenPythia::Generate() if (kf >= fFlavorSelect && kf <= 6) { Int_t idau = iparticle->GetFirstDaughter() - 1; if (idau > -1) { - TParticle* daughter = (TParticle *) fParticles->At(idau); + TParticle* daughter = (TParticle *) fParticles.At(idau); Int_t pdgD = daughter->GetPdgCode(); if (pdgD == 91 || pdgD == 92) { Int_t jmin = daughter->GetFirstDaughter() - 1; Int_t jmax = daughter->GetLastDaughter() - 1; - for (Int_t j = jmin; j <= jmax; j++) - ((TParticle *) fParticles->At(j))->SetFirstMother(i+1); + for (Int_t jp = jmin; jp <= jmax; jp++) + ((TParticle *) fParticles.At(jp))->SetFirstMother(i+1); } // is string or cluster } // has daughter } // heavy quark if (ipa > -1) { - TParticle * mother = (TParticle *) fParticles->At(ipa); + TParticle * mother = (TParticle *) fParticles.At(ipa); kfMo = TMath::Abs(mother->GetPdgCode()); } @@ -745,7 +748,7 @@ void AliGenPythia::Generate() if (nc > 0) { for (i = 0; iAt(i); + TParticle * iparticle = (TParticle *) fParticles.At(i); kf = CheckPDGCode(iparticle->GetPdgCode()); Int_t ks = iparticle->GetStatusCode(); p[0] = iparticle->Px(); @@ -831,15 +834,15 @@ Int_t AliGenPythia::GenerateMB() - Int_t np = (fHadronisation) ? fParticles->GetEntriesFast() : fNpartons; + Int_t np = (fHadronisation) ? fParticles.GetEntriesFast() : fNpartons; Int_t* pParent = new Int_t[np]; for (i=0; i< np; i++) pParent[i] = -1; if (fProcess == kPyJets || fProcess == kPyDirectGamma) { - TParticle* jet1 = (TParticle *) fParticles->At(6); - TParticle* jet2 = (TParticle *) fParticles->At(7); + TParticle* jet1 = (TParticle *) fParticles.At(6); + TParticle* jet2 = (TParticle *) fParticles.At(7); if (!CheckTrigger(jet1, jet2)) { delete [] pParent; return 0; @@ -856,11 +859,11 @@ Int_t AliGenPythia::GenerateMB() else if (fPi0InCalo) pdg = 111 ; // Pi0 for (i=0; i< np; i++) { - TParticle* iparticle = (TParticle *) fParticles->At(i); + TParticle* iparticle = (TParticle *) fParticles.At(i); if(iparticle->GetStatusCode()==1 && iparticle->GetPdgCode()==pdg && iparticle->Pt() > fFragPhotonOrPi0MinPt){ Int_t imother = iparticle->GetFirstMother() - 1; - TParticle* pmother = (TParticle *) fParticles->At(imother); + TParticle* pmother = (TParticle *) fParticles.At(imother); if(pdg == 111 || (pdg == 22 && pmother->GetStatusCode() != 11))//No photon from hadron decay { @@ -876,6 +879,40 @@ Int_t AliGenPythia::GenerateMB() return 0; } + // Check for minimum multiplicity + if (fTriggerMultiplicity > 0) { + Int_t multiplicity = 0; + for (i = 0; i < np; i++) { + TParticle * iparticle = (TParticle *) fParticles.At(i); + + Int_t statusCode = iparticle->GetStatusCode(); + + // Initial state particle + if (statusCode > 20) + continue; + + // skip quarks and gluons + Int_t pdgCode = TMath::Abs(iparticle->GetPdgCode()); + if (pdgCode <= 10 || pdgCode == 21) + continue; + + if (fTriggerMultiplicityEta > 0 && TMath::Abs(iparticle->Eta()) > fTriggerMultiplicityEta) + continue; + + TParticlePDG* pdgPart = iparticle->GetPDG(); + if (pdgPart && pdgPart->Charge() == 0) + continue; + + ++multiplicity; + } + + if (multiplicity < fTriggerMultiplicity) { + delete [] pParent; + return 0; + } + + Printf("Triggered on event with multiplicity of %d > %d", multiplicity, fTriggerMultiplicity); + } // Select events with a photon pt > min pt going to PHOS eta acceptance or exactly PHOS eta phi if ((fProcess == kPyJets || fProcess == kPyDirectGamma) && fPhotonInCalo && (fCheckPHOSeta || fCheckPHOS)){ @@ -885,14 +922,14 @@ Int_t AliGenPythia::GenerateMB() Int_t pdg = 22; Int_t iphcand = -1; for (i=0; i< np; i++) { - TParticle* iparticle = (TParticle *) fParticles->At(i); + TParticle* iparticle = (TParticle *) fParticles.At(i); Float_t phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees Float_t eta =TMath::Abs(iparticle->Eta());//in calos etamin=-etamax if(iparticle->GetStatusCode() == 1 && iparticle->GetPdgCode() == pdg && iparticle->Pt() > fPhotonMinPt - && eta > fPHOSEta){ + && eta < fPHOSEta){ // first check if the photon is in PHOS phi if(IsInPHOS(phi,eta)){ @@ -914,7 +951,7 @@ Int_t AliGenPythia::GenerateMB() if (fTriggerParticle) { Bool_t triggered = kFALSE; for (i = 0; i < np; i++) { - TParticle * iparticle = (TParticle *) fParticles->At(i); + TParticle * iparticle = (TParticle *) fParticles.At(i); kf = CheckPDGCode(iparticle->GetPdgCode()); if (kf != fTriggerParticle) continue; if (iparticle->Pt() == 0.) continue; @@ -937,7 +974,7 @@ Int_t AliGenPythia::GenerateMB() Float_t yQ; Int_t pdgQ; for(i=0; iAt(i); + hvq = (TParticle*)fParticles.At(i); pdgQ = hvq->GetPdgCode(); if(TMath::Abs(pdgQ) != fFlavorSelect) continue; if(pdgQ>0) { theQ=kTRUE; } else { theQbar=kTRUE; } @@ -967,7 +1004,7 @@ Int_t AliGenPythia::GenerateMB() for (i = 0; i < np; i++) { Int_t trackIt = 0; - TParticle * iparticle = (TParticle *) fParticles->At(i); + TParticle * iparticle = (TParticle *) fParticles.At(i); kf = CheckPDGCode(iparticle->GetPdgCode()); Int_t ks = iparticle->GetStatusCode(); Int_t km = iparticle->GetFirstMother(); @@ -1215,10 +1252,10 @@ Bool_t AliGenPythia::CheckKinematicsOnChild(){ Int_t nPartAcc = 0; //number of particles in the acceptance range Int_t numberOfAcceptedParticles = 1; if (fNumberOfAcceptedParticles != 0) { numberOfAcceptedParticles = fNumberOfAcceptedParticles; } - Int_t npart = fParticles->GetEntriesFast(); + Int_t npart = fParticles.GetEntriesFast(); for (j = 0; jAt(j); + TParticle * jparticle = (TParticle *) fParticles.At(j); kcode = TMath::Abs( CheckPDGCode(jparticle->GetPdgCode()) ); ks = jparticle->GetStatusCode(); km = jparticle->GetFirstMother(); @@ -1419,19 +1456,19 @@ void AliGenPythia::RotatePhi(Int_t iphcand, Bool_t& okdd) Double_t phiPHOS = gRandom->Uniform(phiPHOSmin,phiPHOSmax); //calculate deltaphi - TParticle* ph = (TParticle *) fParticles->At(iphcand); + TParticle* ph = (TParticle *) fParticles.At(iphcand); Double_t phphi = ph->Phi(); Double_t deltaphi = phiPHOS - phphi; //loop for all particles and produce the phi rotation - Int_t np = (fHadronisation) ? fParticles->GetEntriesFast() : fNpartons; + Int_t np = (fHadronisation) ? fParticles.GetEntriesFast() : fNpartons; Double_t oldphi, newphi; Double_t newVx, newVy, R, Vz, time; Double_t newPx, newPy, pt, Pz, e; for(Int_t i=0; i< np; i++) { - TParticle* iparticle = (TParticle *) fParticles->At(i); + TParticle* iparticle = (TParticle *) fParticles.At(i); oldphi = iparticle->Phi(); newphi = oldphi + deltaphi; if(newphi < 0) newphi = 2*TMath::Pi() + newphi; // correct angle