#97827 AliGenPythia : Trigger on Pythia particles for decay photons and etas and...
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 1 Oct 2012 17:36:30 +0000 (17:36 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 1 Oct 2012 17:36:30 +0000 (17:36 +0000)
PYTHIA6/AliGenPythia.cxx
PYTHIA6/AliGenPythia.h

index 06658cd..f6511fa 100644 (file)
@@ -134,14 +134,20 @@ AliGenPythia::AliGenPythia():
     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),
@@ -240,14 +246,20 @@ AliGenPythia::AliGenPythia(Int_t npart)
      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),
@@ -925,7 +937,6 @@ Int_t  AliGenPythia::GenerateMB()
 //  converts from mm/c to s
     const Float_t kconv=0.001/2.999792458e8;
     
-
     
     Int_t np = (fHadronisation) ? fParticles.GetEntriesFast() : fNpartons;
 
@@ -942,55 +953,12 @@ Int_t  AliGenPythia::GenerateMB()
        }
     }
 
-  // Select events with fragmentation photon or pi0 or hadrons going to PHOS or EMCAL,
+  
+  // 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 || 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;
-    }
+  if ((fFragPhotonInCalo || fPi0InCalo || fEtaInCalo || fEleInCalo || fHadronInCalo || fDecayPhotonInCalo) && 
+      (fCheckPHOS || fCheckEMCAL || fCheckBarrel) ) {
+    Bool_t ok = TriggerOnSelectedParticles(np);
     
     if(!ok) {
       delete[] pParent;
@@ -998,34 +966,6 @@ Int_t  AliGenPythia::GenerateMB()
     }
   }
 
-  // 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)) ) {
@@ -1068,41 +1008,6 @@ Int_t  AliGenPythia::GenerateMB()
       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;
@@ -1677,7 +1582,17 @@ void AliGenPythia::GetSubEventTime()
   return;
 }
 
-Bool_t AliGenPythia::IsInEMCAL(Float_t phi, Float_t eta) const
+Bool_t AliGenPythia::IsInBarrel(const Float_t eta) const
+{
+  // Is particle in Central Barrel acceptance? 
+  // etamin=-etamax
+  if( eta < fTriggerEta  ) 
+    return kTRUE;
+  else 
+    return kFALSE;
+}
+
+Bool_t AliGenPythia::IsInEMCAL(const Float_t phi, const Float_t eta) const
 {
   // Is particle in EMCAL acceptance? 
   // phi in degrees, etamin=-etamax
@@ -1688,30 +1603,40 @@ Bool_t AliGenPythia::IsInEMCAL(Float_t phi, Float_t eta) const
     return kFALSE;
 }
 
-Bool_t AliGenPythia::IsInPHOS(Float_t phi, Float_t eta) const
+Bool_t AliGenPythia::IsInPHOS(const Float_t phi, const Float_t eta, const 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
@@ -1720,40 +1645,42 @@ void AliGenPythia::RotatePhi(Int_t iphcand, Bool_t& okdd)
   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!
@@ -2417,3 +2344,225 @@ Bool_t AliGenPythia::IsFromHeavyFlavor(Int_t ipart)
   }
   return kFALSE;
 }
+
+Bool_t AliGenPythia::CheckDetectorAcceptance(const Float_t phi, const Float_t eta, const 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(const 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;
+}
+
index b493df8..dc99b4a 100644 (file)
@@ -89,23 +89,56 @@ class AliGenPythia : public AliGenMC
     // Phi range for gamma trigger
     virtual void    SetGammaPhiRange(Float_t phimin = 0., Float_t phimax = 360.)
        {fPhiMinGamma = TMath::Pi()*phimin/180.; fPhiMaxGamma = TMath::Pi()*phimax/180.;}
-   // Select jets with fragmentation photon or pi0 going to PHOS or EMCAL
-    virtual void  SetFragPhotonInCalo(Bool_t b)  {fFragPhotonInCalo = b;}
-    virtual void  SetHadronInCalo    (Bool_t b)  {fHadronInCalo = b;}
-    virtual void  SetPi0InCalo       (Bool_t b)  {fPi0InCalo    = b;}
-    virtual void  SetPhotonInCalo(Bool_t b)      {fPhotonInCalo = b;}
-    virtual void  SetCheckPHOS (Bool_t b)        {fCheckPHOS    = b;}
-    virtual void  SetCheckEMCAL(Bool_t b)        {fCheckEMCAL   = b;}
-    virtual void  SetFragPhotonInEMCAL(Bool_t b) {fCheckEMCAL   = b; fFragPhotonInCalo = b;}
-    virtual void  SetFragPhotonInPHOS(Bool_t b)  {fCheckPHOS    = b; fFragPhotonInCalo = b;}
-    virtual void  SetHadronInEMCAL(Bool_t b)     {fCheckEMCAL   = b; fHadronInCalo     = b;}
-    virtual void  SetHadronInPHOS(Bool_t b)      {fCheckPHOS    = b; fHadronInCalo     = b;}
-    virtual void  SetPi0InEMCAL(Bool_t b)        {fCheckEMCAL   = b; fPi0InCalo        = b;}
-    virtual void  SetPi0InPHOS(Bool_t b)         {fCheckPHOS    = b; fPi0InCalo        = b;}
-    virtual void  SetPhotonInEMCAL(Bool_t b)     {fCheckEMCAL   = b; fPhotonInCalo     = b;}
-    virtual void  SetElectronInEMCAL(Bool_t b)   {fEleInEMCAL   = b;}
-    virtual void  SetPhotonInPHOS(Bool_t b)      {fCheckPHOS    = b; fPhotonInCalo     = b;}
+  
+    // Select events with fragmentation photon, decay photon, pi0 or eta going to PHOS or EMCAL and central barrel
+    virtual Bool_t TriggerOnSelectedParticles(const Int_t np);
+  
+    virtual void  SetCheckPHOS         (Bool_t b) {fCheckPHOS    = b;}
+    virtual void  SetCheckEMCAL        (Bool_t b) {fCheckEMCAL   = b;}
+    virtual void  SetCheckBarrel       (Bool_t b) {fCheckBarrel  = b;}
+  
+    //virtual void  SetElectronInEMCAL   (Bool_t b) {fEleInEMCAL   = b;}
+    //virtual void  SetPhotonInPHOS      (Bool_t b) {fCheckPHOS    = b; fPhotonInCalo     = b;} // Not in use
 
+    virtual void  SetFragPhotonInCalo  (Bool_t b) {                   fFragPhotonInCalo = b;}
+    virtual void  SetFragPhotonInBarrel(Bool_t b) {fCheckBarrel  = b; fFragPhotonInCalo = b;}
+    virtual void  SetFragPhotonInEMCAL (Bool_t b) {fCheckEMCAL   = b; fFragPhotonInCalo = b;}
+    virtual void  SetFragPhotonInPHOS  (Bool_t b) {fCheckPHOS    = b; fFragPhotonInCalo = b;}
+  
+    virtual void  SetHadronInCalo      (Bool_t b) {                   fHadronInCalo     = b;}
+    virtual void  SetHadronInBarrel    (Bool_t b) {fCheckBarrel  = b; fHadronInCalo     = b;}
+    virtual void  SetHadronInEMCAL     (Bool_t b) {fCheckEMCAL   = b; fHadronInCalo     = b;}
+    virtual void  SetHadronInPHOS      (Bool_t b) {fCheckPHOS    = b; fHadronInCalo     = b;}
+  
+    virtual void  SetElectronInCalo    (Bool_t b) {                   fEleInCalo        = b;}
+    virtual void  SetElectronInBarrel  (Bool_t b) {fCheckBarrel  = b; fEleInCalo        = b;}
+    virtual void  SetElectronInEMCAL   (Bool_t b) {fCheckEMCAL   = b; fEleInCalo        = b;}
+    virtual void  SetElectronInPHOS    (Bool_t b) {fCheckPHOS    = b; fEleInCalo        = b;}
+  
+    virtual void  SetDecayPhotonInCalo (Bool_t d)  {fDecayPhotonInCalo = d;}
+    virtual void  SetDecayPhotonInBarrel(Bool_t d) {fDecayPhotonInCalo = d; fCheckBarrel  = d;}
+    virtual void  SetDecayPhotonInEMCAL(Bool_t d)  {fDecayPhotonInCalo = d; fCheckEMCAL   = d;}
+    virtual void  SetDecayPhotonInPHOS (Bool_t d)  {fDecayPhotonInCalo = d; fCheckPHOS    = d;}
+  
+    virtual void  SetPi0InCalo         (Bool_t b, Bool_t f = kFALSE) {fPi0InCalo = b; fForceNeutralMeson2PhotonDecay = f;}
+    virtual void  SetPi0InBarrel       (Bool_t b, Bool_t f = kFALSE) {fPi0InCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckBarrel= b; }
+    virtual void  SetPi0InEMCAL        (Bool_t b, Bool_t f = kFALSE) {fPi0InCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckEMCAL = b; }
+    virtual void  SetPi0InPHOS         (Bool_t b, Bool_t f = kFALSE) {fPi0InCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckPHOS  = b; }
+    virtual void  SetEtaInCalo         (Bool_t b, Bool_t f = kFALSE) {fEtaInCalo = b; fForceNeutralMeson2PhotonDecay = f;}
+    virtual void  SetEtaInBarrel       (Bool_t b, Bool_t f = kFALSE) {fEtaInCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckBarrel= b; }
+    virtual void  SetEtaInEMCAL        (Bool_t b, Bool_t f = kFALSE) {fEtaInCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckEMCAL = b; }
+    virtual void  SetEtaInPHOS         (Bool_t b, Bool_t f = kFALSE) {fEtaInCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckPHOS  = b; }
+
+    virtual void  SetPi0PhotonDecayInBarrel(Bool_t b, Bool_t f = kFALSE) {fPi0InCalo = b; fDecayPhotonInCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckBarrel  = b; }
+    virtual void  SetPi0PhotonDecayInEMCAL (Bool_t b, Bool_t f = kFALSE) {fPi0InCalo = b; fDecayPhotonInCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckEMCAL   = b; }
+    virtual void  SetPi0PhotonDecayInPHOS  (Bool_t b, Bool_t f = kFALSE) {fPi0InCalo = b; fDecayPhotonInCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckPHOS    = b; }
+
+    virtual void  SetEtaPhotonDecayInBarrel(Bool_t b, Bool_t f = kFALSE) {fEtaInCalo = b; fDecayPhotonInCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckBarrel  = b; }
+    virtual void  SetEtaPhotonDecayInEMCAL (Bool_t b, Bool_t f = kFALSE) {fEtaInCalo = b; fDecayPhotonInCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckEMCAL   = b; }
+    virtual void  SetEtaPhotonDecayInPHOS  (Bool_t b, Bool_t f = kFALSE) {fEtaInCalo = b; fDecayPhotonInCalo = b; fForceNeutralMeson2PhotonDecay = f; fCheckPHOS    = b; }
+
+  
     // Trigger on a minimum multiplicity
     virtual void  SetTriggerChargedMultiplicity(Int_t multiplicity, Float_t etamax = 0, Float_t ptmin = -1.) 
     {fTriggerMultiplicity = multiplicity; fTriggerMultiplicityEta = etamax; 
@@ -113,16 +146,18 @@ class AliGenPythia : public AliGenMC
        
     // Calorimeters acceptance
     // Set Phi in degrees, and Eta coverage, should not be negative
-    virtual void  SetEMCALAcceptance(Float_t phimin, Float_t phimax, Float_t deta) {fEMCALMinPhi = phimin ; fEMCALMaxPhi = phimax ; fEMCALEta = deta ; }
-    virtual void  SetPHOSAcceptance (Float_t phimin, Float_t phimax, Float_t deta) {fPHOSMinPhi  = phimin ; fPHOSMaxPhi  = phimax ; fPHOSEta  = deta ; }
-
-    virtual void  SetPhotonInPHOSeta(Bool_t b)        {fCheckPHOSeta = b; fPhotonInCalo     = b;}
+    virtual void  SetBarrelAcceptance(Float_t deta) {fTriggerEta = deta ;}
+    virtual void  SetEMCALAcceptance (Float_t phimin, Float_t phimax, Float_t deta) {fEMCALMinPhi = phimin ; fEMCALMaxPhi = phimax ; fEMCALEta = deta ; }
+    virtual void  SetPHOSAcceptance  (Float_t phimin, Float_t phimax, Float_t deta) {fPHOSMinPhi  = phimin ; fPHOSMaxPhi  = phimax ; fPHOSEta  = deta ; }
+    virtual void  SetRotateParticleInPHOSeta(Bool_t b) {fCheckPHOSeta = b;}
+  
     virtual void  SetTriggerParticleMinPt(Float_t pt) {fTriggerParticleMinPt = pt;}
-    virtual void  SetPhotonMinPt(Float_t pt)          {fPhotonMinPt = pt;}
-    virtual void  SetElectronMinPt(Float_t pt)        {fElectronMinPt = pt;}
+//    virtual void  SetPhotonMinPt(Float_t pt)          {fPhotonMinPt = pt;}
+//    virtual void  SetElectronMinPt(Float_t pt)        {fElectronMinPt = pt;}
     // Trigger and rotate event 
-    void RotatePhi(Int_t iphcand, Bool_t& okdd);
-    // Trigger on a single particle
+    void RotatePhi(Bool_t& okdd);
+  
+    // Trigger on a single particle (not related to calorimeter trigger above)
     virtual void    SetTriggerParticle(Int_t particle = 0, Float_t etamax = 0.9, Float_t ptmin = -1, Float_t ptmax = 1000) 
        {fTriggerParticle = particle; fTriggerEta = etamax; fTriggerMinPt = ptmin; fTriggerMaxPt = ptmax;}
 
@@ -195,8 +230,10 @@ class AliGenPythia : public AliGenMC
     virtual void         GetGammaPhiRange(Float_t& phimin, Float_t& phimax) const
        {phimin = fPhiMinGamma*180./TMath::Pi(); phimax = fPhiMaxGamma*180./TMath::Pi();}
     //
-    Bool_t IsInEMCAL(Float_t phi, Float_t eta) const;
-    Bool_t IsInPHOS(Float_t phi, Float_t eta) const;
+    Bool_t CheckDetectorAcceptance(const Float_t phi, const Float_t eta, const Int_t iparticle);
+    Bool_t IsInEMCAL (const Float_t phi, const Float_t eta) const;
+    Bool_t IsInPHOS  (const Float_t phi, const Float_t eta, const Int_t iparticle) ;
+    Bool_t IsInBarrel(const Float_t eta) const;
     Bool_t IsFromHeavyFlavor(Int_t ipart);
     //
     virtual void FinishRun();
@@ -308,14 +345,20 @@ class AliGenPythia : public AliGenMC
     Bool_t fFragPhotonInCalo; // Option to ask for Fragmentation Photon in calorimeters acceptance
     Bool_t fHadronInCalo;     // Option to ask for hadron (not pi0) in calorimeters acceptance
     Bool_t fPi0InCalo;        // Option to ask for Pi0 in calorimeters acceptance
-    Bool_t fPhotonInCalo;     // Option to ask for Decay Photon in calorimeter acceptance
-    Bool_t fEleInEMCAL;       // Option to ask for Electron in EMCAL acceptance
-    Bool_t fCheckEMCAL;       // Option to ask for FragPhoton or Pi0 in calorimeters EMCAL acceptance
-    Bool_t fCheckPHOS;        // Option to ask for FragPhoton or Pi0 in calorimeters PHOS acceptance
-    Bool_t fCheckPHOSeta;     // Option to ask for PHOS eta acceptance
+    Bool_t fEtaInCalo;        // Option to ask for Eta in calorimeters acceptance
+    Bool_t fPhotonInCalo;     // Option to ask for Photon in calorimeter acceptance (not in use)
+    Bool_t fDecayPhotonInCalo;// Option to ask for Decay Photon in calorimeter acceptance
+    Bool_t fForceNeutralMeson2PhotonDecay; // Option to ask for Pi0/Eta in calorimeters acceptance when decay into 2 photons
+    Bool_t fEleInCalo;        // Option to ask for Electron in EMCAL acceptance
+    Bool_t fEleInEMCAL;       // Option to ask for Electron in EMCAL acceptance (not in use)
+    Bool_t fCheckBarrel;      // Option to ask for FragPhoton or Pi0 or Eta or gamma decays in central barrel acceptance
+    Bool_t fCheckEMCAL;       // Option to ask for FragPhoton or Pi0 or Eta or gamma decays in calorimeters EMCAL acceptance
+    Bool_t fCheckPHOS;        // Option to ask for FragPhoton or Pi0 or Eta or gamma decays in calorimeters PHOS acceptance
+    Bool_t fCheckPHOSeta;     // Option to ask for rotate event particles in phi to have in PHOS acceptance a requested particle that previously had the good eta
+    Int_t  fPHOSRotateCandidate; // Internal member to select the particle candidate to trigger the event phi rotation, to put it in PHOS phi acceptance
     Float_t fTriggerParticleMinPt; // Minimum momentum of Fragmentation Photon or Pi0 or other hadron
-    Float_t fPhotonMinPt;          // Minimum momentum of Photon 
-    Float_t fElectronMinPt;        // Minimum momentum of Electron 
+    Float_t fPhotonMinPt;          // Minimum momentum of Photon  (not in use)
+    Float_t fElectronMinPt;        // Minimum momentum of Electron (not in use) 
     //Calorimeters eta-phi acceptance 
     Float_t fPHOSMinPhi;           // Minimum phi PHOS, degrees
     Float_t fPHOSMaxPhi;           // Maximum phi PHOS, degrees
@@ -335,7 +378,7 @@ class AliGenPythia : public AliGenMC
     Bool_t GetWeightsDiffraction(Double_t M, Double_t &Mmin, Double_t &Mmax, 
                                               Double_t &wSD, Double_t &wDD, Double_t &wND);
 
-    ClassDef(AliGenPythia, 13) // AliGenerator interface to Pythia
+    ClassDef(AliGenPythia, 14) // AliGenerator interface to Pythia
 };
 #endif