fFragPhotonInCalo(kFALSE),
fHadronInCalo(kFALSE) ,
fPi0InCalo(kFALSE) ,
- fPhotonInCalo(kFALSE),
- fEleInEMCAL(kFALSE),
+ fEtaInCalo(kFALSE) ,
+ fPhotonInCalo(kFALSE), // not in use
+ fDecayPhotonInCalo(kFALSE),
+ fForceNeutralMeson2PhotonDecay(kFALSE),
+ fEleInCalo(kFALSE),
+ fEleInEMCAL(kFALSE), // not in use
+ fCheckBarrel(kFALSE),
fCheckEMCAL(kFALSE),
fCheckPHOS(kFALSE),
fCheckPHOSeta(kFALSE),
+ fPHOSRotateCandidate(-1),
fTriggerParticleMinPt(0),
- fPhotonMinPt(0),
- fElectronMinPt(0),
+ fPhotonMinPt(0), // not in use
+ fElectronMinPt(0), // not in use
fPHOSMinPhi(219.),
fPHOSMaxPhi(321.),
fPHOSEta(0.13),
fFragPhotonInCalo(kFALSE),
fHadronInCalo(kFALSE) ,
fPi0InCalo(kFALSE) ,
- fPhotonInCalo(kFALSE),
- fEleInEMCAL(kFALSE),
+ fEtaInCalo(kFALSE) ,
+ fPhotonInCalo(kFALSE), // not in use
+ fDecayPhotonInCalo(kFALSE),
+ fForceNeutralMeson2PhotonDecay(kFALSE),
+ fEleInCalo(kFALSE),
+ fEleInEMCAL(kFALSE), // not in use
+ fCheckBarrel(kFALSE),
fCheckEMCAL(kFALSE),
fCheckPHOS(kFALSE),
fCheckPHOSeta(kFALSE),
+ fPHOSRotateCandidate(-1),
fTriggerParticleMinPt(0),
- fPhotonMinPt(0),
- fElectronMinPt(0),
+ fPhotonMinPt(0), // not in use
+ fElectronMinPt(0), // not in use
fPHOSMinPhi(219.),
fPHOSMaxPhi(321.),
fPHOSEta(0.13),
break;
case kPyW:
case kPyZ:
+ case kPyMBRSingleDiffraction:
+ case kPyMBRDoubleDiffraction:
+ case kPyMBRCentralDiffraction:
break;
}
//
// converts from mm/c to s
const Float_t kconv=0.001/2.999792458e8;
-
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 || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi) {
- TParticle* jet1 = (TParticle *) fParticles.At(6);
+ TParticle* jet1 = (TParticle *) fParticles.At(6);
TParticle* jet2 = (TParticle *) fParticles.At(7);
- if (!CheckTrigger(jet1, jet2)) {
+
+ if (!jet1 || ! jet2 || !CheckTrigger(jet1, jet2)) {
delete [] pParent;
return 0;
}
}
- // Select jets with fragmentation photon or pi0 or hadrons going to PHOS or EMCAL
- if ( fProcess == kPyJets &&
- (fFragPhotonInCalo || fPi0InCalo || fHadronInCalo) &&
- (fCheckPHOS || fCheckEMCAL) ) {
-
- Bool_t ok = kFALSE;
-
- for (i=0; i< np; i++) {
-
- TParticle* iparticle = (TParticle *) fParticles.At(i);
-
- Int_t pdg = iparticle->GetPdgCode();
- Int_t status = iparticle->GetStatusCode();
- ok = kFALSE;
-
- if (fFragPhotonInCalo && pdg == 22 && status == 1)
- {
- Int_t imother = iparticle->GetFirstMother() - 1;
- TParticle* pmother = (TParticle *) fParticles.At(imother);
-
- if(pmother->GetStatusCode() != 11) ok = kTRUE ; // No photon from hadron decay
- }
- else if (fPi0InCalo && pdg == 111) // pi0 status can be 1 or 11 depending on decay settings
- {
- //printf("Pi0! pdg %d, status %d, pt %2.2f\n",pdg,status,iparticle->Pt());
- ok = kTRUE;
- }
- else if (fHadronInCalo && status == 1)
- {
- if(TMath::Abs(pdg) > 23 && pdg !=221 && pdg != 111) // avoid photons, electrons, muons, neutrinos and eta or pi0
- // (in case neutral mesons were declared stable)
- ok = kTRUE;
- }
-
- if(ok && iparticle->Pt() > fTriggerParticleMinPt)
- {
- Float_t phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees
- Float_t eta =TMath::Abs(iparticle->Eta()); //in calos etamin=-etamax
-
- if((fCheckEMCAL && IsInEMCAL(phi,eta)) ||
- (fCheckPHOS && IsInPHOS (phi,eta)) )
- {
- ok =kTRUE;
- AliDebug(1,Form("Selected trigger pdg %d, status %d, pt %2.2f, eta %2.2f, phi %2.2f\n",pdg,status,iparticle->Pt(), eta, phi));
- break;
- }
- }
- else ok = kFALSE;
- }
+
+ // Select events with fragmentation photon, decay photon, pi0, eta or other hadrons going to PHOS or EMCAL or central barrel,
+ // implemented primaryly for kPyJets, but extended to any kind of process.
+ if ((fFragPhotonInCalo || fPi0InCalo || fEtaInCalo || fEleInCalo || fHadronInCalo || fDecayPhotonInCalo) &&
+ (fCheckPHOS || fCheckEMCAL || fCheckBarrel) ) {
+ Bool_t ok = TriggerOnSelectedParticles(np);
if(!ok) {
delete[] pParent;
}
}
- // Select beauty jets with electron in EMCAL
- if (fProcess == kPyBeautyJets && fEleInEMCAL) {
-
- Bool_t ok = kFALSE;
-
- Int_t pdg = 11; //electron
-
- Float_t pt = 0.;
- Float_t eta = 0.;
- Float_t phi = 0.;
- for (i=0; i< np; i++) {
- TParticle* iparticle = (TParticle *) fParticles.At(i);
- if(iparticle->GetStatusCode()==1 && TMath::Abs(iparticle->GetPdgCode())==pdg &&
- iparticle->Pt() > fElectronMinPt){
- pt = iparticle->Pt();
- phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees
- eta =TMath::Abs(iparticle->Eta()); //in calos etamin=-etamax
- if(IsInEMCAL(phi,eta))
- ok =kTRUE;
- }
- }
- if(!ok) {
- delete[] pParent;
- return 0;
- }
-
- AliDebug(5,Form("Found an electron jet (pt,eta,phi) = (%f,%f,%f)",pt,eta,phi));
- }
// Check for diffraction
if(fkTuneForDiff) {
if(fItune==320 && ( (TMath::Abs(fEnergyCMS - 900) < 1) || (TMath::Abs(fEnergyCMS - 2760) < 1) || (TMath::Abs(fEnergyCMS - 7000) < 1)) ) {
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)){
-
- Bool_t okd = kFALSE;
-
- Int_t pdg = 22;
- Int_t iphcand = -1;
- for (i=0; i< np; 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){
-
- // first check if the photon is in PHOS phi
- if(IsInPHOS(phi,eta)){
- okd = kTRUE;
- break;
- }
- if(fCheckPHOSeta) iphcand = i; // candiate photon to rotate in phi
-
- }
- }
-
- if(!okd && iphcand != -1) // execute rotation in phi
- RotatePhi(iphcand,okd);
-
- if(!okd) {
- delete [] pParent;
- return 0;
- }
- }
if (fTriggerParticle) {
Bool_t triggered = kFALSE;
fAProjectile = a1;
fATarget = a2;
- fNucPdf = pdfset; // 0 EKS98 1 EPS08
+ fNucPdf = pdfset; // 0 EKS98 9 EPS09LO 19 EPS09NLO
fSetNuclei = kTRUE;
}
return;
}
+Bool_t AliGenPythia::IsInBarrel(Float_t eta) const
+{
+ // Is particle in Central Barrel acceptance?
+ // etamin=-etamax
+ if( eta < fTriggerEta )
+ return kTRUE;
+ else
+ return kFALSE;
+}
+
Bool_t AliGenPythia::IsInEMCAL(Float_t phi, Float_t eta) const
{
// Is particle in EMCAL acceptance?
return kFALSE;
}
-Bool_t AliGenPythia::IsInPHOS(Float_t phi, Float_t eta) const
+Bool_t AliGenPythia::IsInPHOS(Float_t phi, Float_t eta, Int_t iparticle)
{
// Is particle in PHOS acceptance?
// Acceptance slightly larger considered.
// phi in degrees, etamin=-etamax
+ // iparticle is the index of the particle to be checked, for PHOS rotation case
+
if(phi > fPHOSMinPhi && phi < fPHOSMaxPhi &&
eta < fPHOSEta )
return kTRUE;
else
+ {
+ if( fCheckPHOSeta && eta < fPHOSEta) fPHOSRotateCandidate = iparticle;
+
return kFALSE;
+ }
}
-void AliGenPythia::RotatePhi(Int_t iphcand, Bool_t& okdd)
+void AliGenPythia::RotatePhi(Bool_t& okdd)
{
+ //Rotate event in phi to enhance events in PHOS acceptance
+
+ if(fPHOSRotateCandidate < 0) return ;
+
//calculate the new position random between fPHOSMinPhi and fPHOSMaxPhi
Double_t phiPHOSmin = TMath::Pi()*fPHOSMinPhi/180;
Double_t phiPHOSmax = TMath::Pi()*fPHOSMaxPhi/180;
Double_t phiPHOS = (AliPythiaRndm::GetPythiaRandom())->Uniform(phiPHOSmin,phiPHOSmax);
//calculate deltaphi
- TParticle* ph = (TParticle *) fParticles.At(iphcand);
+ TParticle* ph = (TParticle *) fParticles.At(fPHOSRotateCandidate);
Double_t phphi = ph->Phi();
Double_t deltaphi = phiPHOS - phphi;
-
+
//loop for all particles and produce the phi rotation
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);
- oldphi = iparticle->Phi();
- newphi = oldphi + deltaphi;
- if(newphi < 0) newphi = 2*TMath::Pi() + newphi; // correct angle
- if(newphi > 2*TMath::Pi()) newphi = newphi - 2*TMath::Pi(); // correct angle
-
- r = iparticle->R();
- newVx = r * TMath::Cos(newphi);
- newVy = r * TMath::Sin(newphi);
- vZ = iparticle->Vz(); // don't transform
- time = iparticle->T(); // don't transform
-
- pt = iparticle->Pt();
- newPx = pt * TMath::Cos(newphi);
- newPy = pt * TMath::Sin(newphi);
- pz = iparticle->Pz(); // don't transform
- e = iparticle->Energy(); // don't transform
-
- // apply rotation
- iparticle->SetProductionVertex(newVx, newVy, vZ, time);
- iparticle->SetMomentum(newPx, newPy, pz, e);
-
+ TParticle* iparticle = (TParticle *) fParticles.At(i);
+ oldphi = iparticle->Phi();
+ newphi = oldphi + deltaphi;
+ if(newphi < 0) newphi = 2*TMath::Pi() + newphi; // correct angle
+ if(newphi > 2*TMath::Pi()) newphi = newphi - 2*TMath::Pi(); // correct angle
+
+ r = iparticle->R();
+ newVx = r * TMath::Cos(newphi);
+ newVy = r * TMath::Sin(newphi);
+ vZ = iparticle->Vz(); // don't transform
+ time = iparticle->T(); // don't transform
+
+ pt = iparticle->Pt();
+ newPx = pt * TMath::Cos(newphi);
+ newPy = pt * TMath::Sin(newphi);
+ pz = iparticle->Pz(); // don't transform
+ e = iparticle->Energy(); // don't transform
+
+ // apply rotation
+ iparticle->SetProductionVertex(newVx, newVy, vZ, time);
+ iparticle->SetMomentum(newPx, newPy, pz, e);
+
} //end particle loop
- // now let's check that we put correctly the candidate photon in PHOS
- Float_t phi = ph->Phi()*180./TMath::Pi(); //Convert to degrees
- Float_t eta =TMath::Abs(ph->Eta());//in calos etamin=-etamax
- if(IsInPHOS(phi,eta))
- okdd = kTRUE;
+ // now let's check that we put correctly the candidate photon in PHOS
+ Float_t phi = ph->Phi()*180./TMath::Pi(); //Convert to degrees
+ Float_t eta =TMath::Abs(ph->Eta());//in calos etamin=-etamax
+ if(IsInPHOS(phi,eta,-1))
+ okdd = kTRUE;
+
+ // reset the value for next event
+ fPHOSRotateCandidate = -1;
+
}
-
-
Bool_t AliGenPythia::CheckDiffraction()
{
// use this method only with Perugia-0 tune!
}
return kFALSE;
}
+
+Bool_t AliGenPythia::CheckDetectorAcceptance(Float_t phi, Float_t eta, Int_t iparticle)
+{
+ // check the eta/phi correspond to the detectors acceptance
+ // iparticle is the index of the particle to be checked, for PHOS rotation case
+ if (fCheckPHOS && IsInPHOS (phi,eta,iparticle)) return kTRUE;
+ else if(fCheckEMCAL && IsInEMCAL (phi,eta)) return kTRUE;
+ else if(fCheckBarrel && IsInBarrel( eta)) return kTRUE;
+ else return kFALSE;
+}
+
+Bool_t AliGenPythia::TriggerOnSelectedParticles(Int_t np)
+{
+ // Select events with fragmentation photon, decay photon, pi0, eta or other hadrons going to PHOS or EMCAL or central barrel,
+ // implemented primaryly for kPyJets, but extended to any kind of process.
+
+ //printf("Check: frag photon %d, pi0 %d, eta %d, electron %d, hadron %d, decay %d, PHOS %d, EMCAL %d, Barrel %d \n",
+ // fFragPhotonInCalo,fPi0InCalo, fEtaInCalo,fEleInCalo,fHadronInCalo,fDecayPhotonInCalo,fCheckPHOS,fCheckEMCAL, fCheckBarrel);
+
+ Bool_t ok = kFALSE;
+ for (Int_t i=0; i< np; i++) {
+
+ TParticle* iparticle = (TParticle *) fParticles.At(i);
+
+ Int_t pdg = iparticle->GetPdgCode();
+ Int_t status = iparticle->GetStatusCode();
+ Int_t imother = iparticle->GetFirstMother() - 1;
+
+ TParticle* pmother = 0x0;
+ Int_t momStatus = -1;
+ Int_t momPdg = -1;
+ if(imother > 0 ){
+ pmother = (TParticle *) fParticles.At(imother);
+ momStatus = pmother->GetStatusCode();
+ momPdg = pmother->GetPdgCode();
+ }
+
+ ok = kFALSE;
+
+ //
+ // Check the particle type: hadron (not pi0 or eta), electron, decay photon (from pi0 or eta or any), pi0 or eta
+ //
+ // Hadron
+ if (fHadronInCalo && status == 1)
+ {
+ if(TMath::Abs(pdg) > 23 && pdg !=221 && pdg != 111) // avoid photons, electrons, muons, neutrinos and eta or pi0
+ // (in case neutral mesons were declared stable)
+ ok = kTRUE;
+ }
+ //Electron
+ else if (fEleInCalo && status == 1 && TMath::Abs(pdg) == 11)
+ {
+ ok = kTRUE;
+ }
+ //Fragmentation photon
+ else if (fFragPhotonInCalo && pdg == 22 && status == 1)
+ {
+ if(momStatus != 11) ok = kTRUE ; // No photon from hadron decay
+ }
+ // Decay photon
+ else if (fDecayPhotonInCalo && !fForceNeutralMeson2PhotonDecay && pdg == 22) // pi0 status can be 1 or 11 depending on decay settings, work only for 11
+ {
+ if( momStatus == 11)
+ {
+ //if(iparticle->Pt() > fTriggerParticleMinPt) printf("Decay photon! pdg %d, status %d, pt %2.2f, mom: pdg %d, pt %2.2f\n",
+ // pdg,status,iparticle->Pt(),momPdg,pmother->Pt());
+ ok = kTRUE ; // photon from hadron decay
+
+ //In case only decays from pi0 or eta requested
+ if(fPi0InCalo && momPdg!=111) ok = kFALSE;
+ if(fEtaInCalo && momPdg!=221) ok = kFALSE;
+ }
+
+ }
+ // Pi0 or Eta particle
+ else if ((fPi0InCalo || fEtaInCalo))
+ {
+ if(fDecayPhotonInCalo && !fForceNeutralMeson2PhotonDecay ) continue ;
+
+ if (fPi0InCalo && pdg == 111) // pi0 status can be 1 or 11 depending on decay settings
+ {
+ //if(iparticle->Pt() > fTriggerParticleMinPt) printf("Pi0! pdg %d, status %d, pt %2.2f\n",pdg,status,iparticle->Pt());
+ ok = kTRUE;
+ }
+ else if (fEtaInCalo && pdg == 221)
+ {
+ //if(iparticle->Pt() > fTriggerParticleMinPt) printf("Eta! pdg %d, status %d, pt %2.2f\n",pdg,status,iparticle->Pt());
+ ok = kTRUE;
+ }
+
+ }// pi0 or eta
+
+ //
+ // Check that the selected particle is in the calorimeter acceptance
+ //
+ if(ok && iparticle->Pt() > fTriggerParticleMinPt)
+ {
+ //Just check if the selected particle falls in the acceptance
+ if(!fForceNeutralMeson2PhotonDecay )
+ {
+ //printf("\t Check acceptance! \n");
+ Float_t phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees
+ Float_t eta =TMath::Abs(iparticle->Eta()); //in calos etamin=-etamax
+
+ if(CheckDetectorAcceptance(phi,eta,i))
+ {
+ ok =kTRUE;
+ AliDebug(1,Form("Selected trigger pdg %d, status %d, pt %2.2f, eta %2.2f, phi %2.2f\n",pdg,status,iparticle->Pt(), eta, phi));
+ //printf("\t Accept \n");
+ break;
+ }
+ else ok = kFALSE;
+ }
+ //Mesons have several decay modes, select only those decaying into 2 photons
+ else if(fForceNeutralMeson2PhotonDecay && (fPi0InCalo || fEtaInCalo))
+ {
+ // In case we want the pi0/eta trigger,
+ // check the decay mode (2 photons)
+
+ //printf("\t Force decay 2 gamma\n");
+
+ Int_t ndaughters = iparticle->GetNDaughters();
+ if(ndaughters != 2){
+ ok=kFALSE;
+ continue;
+ }
+
+ TParticle*d1 = (TParticle *) fParticles.At(iparticle->GetDaughter(0)-1);
+ TParticle*d2 = (TParticle *) fParticles.At(iparticle->GetDaughter(1)-1);
+ if(!d1 || !d2) {
+ ok=kFALSE;
+ continue;
+ }
+
+ //iparticle->Print();
+ //d1->Print();
+ //d2->Print();
+
+ Int_t pdgD1 = d1->GetPdgCode();
+ Int_t pdgD2 = d2->GetPdgCode();
+ //printf("\t \t 2 daughters pdg = %d - %d\n",pdgD1,pdgD2);
+ //printf("mother %d - %d\n",d1->GetFirstMother(),d2->GetFirstMother());
+
+ if(pdgD1 != 22 || pdgD2 != 22){
+ ok = kFALSE;
+ continue;
+ }
+
+ //printf("\t accept decay\n");
+
+ //Trigger on the meson, not on the daughter
+ if(!fDecayPhotonInCalo){
+
+ Float_t phi = iparticle->Phi()*180./TMath::Pi(); //Convert to degrees
+ Float_t eta =TMath::Abs(iparticle->Eta()); //in calos etamin=-etamax
+
+ if(CheckDetectorAcceptance(phi,eta,i))
+ {
+ //printf("\t Accept meson pdg %d\n",pdg);
+ ok =kTRUE;
+ AliDebug(1,Form("Selected trigger pdg %d (decay), status %d, pt %2.2f, eta %2.2f, phi %2.2f\n",pdg,status,iparticle->Pt(), eta, phi));
+ break;
+ } else {
+ ok=kFALSE;
+ continue;
+ }
+ }
+
+ //printf("Check daughters acceptance\n");
+
+ //Trigger on the meson daughters
+ //Photon 1
+ Float_t phi = d1->Phi()*180./TMath::Pi(); //Convert to degrees
+ Float_t eta =TMath::Abs(d1->Eta()); //in calos etamin=-etamax
+ if(d1->Pt() > fTriggerParticleMinPt)
+ {
+ //printf("\t Check acceptance photon 1! \n");
+ if(CheckDetectorAcceptance(phi,eta,i))
+ {
+ //printf("\t Accept Photon 1\n");
+ ok =kTRUE;
+ AliDebug(1,Form("Selected trigger pdg %d (decay), status %d, pt %2.2f, eta %2.2f, phi %2.2f\n",pdg,status,iparticle->Pt(), eta, phi));
+ break;
+ }
+ else ok = kFALSE;
+ } // pt cut
+ else ok = kFALSE;
+
+ //Photon 2
+ phi = d2->Phi()*180./TMath::Pi(); //Convert to degrees
+ eta =TMath::Abs(d2->Eta()); //in calos etamin=-etamax
+
+ if(d2->Pt() > fTriggerParticleMinPt)
+ {
+ //printf("\t Check acceptance photon 2! \n");
+ if(CheckDetectorAcceptance(phi,eta,i))
+ {
+ //printf("\t Accept Photon 2\n");
+ ok =kTRUE;
+ AliDebug(1,Form("Selected trigger pdg %d (decay), status %d, pt %2.2f, eta %2.2f, phi %2.2f\n",pdg,status,iparticle->Pt(), eta, phi));
+ break;
+ }
+ else ok = kFALSE;
+ } // pt cut
+ else ok = kFALSE;
+ } // force 2 photon daughters in pi0/eta decays
+ else ok = kFALSE;
+ } else ok = kFALSE; // check acceptance
+ } // primary loop
+
+ //
+ // If requested, rotate the particles event in phi to enhance/speed PHOS selection
+ // A particle passing all trigger conditions except phi position in PHOS, is used as reference
+ //
+ if(fCheckPHOSeta)
+ {
+ RotatePhi(ok);
+ }
+
+ return ok;
+}
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff