Correction in the photon generation (Gustavo)

[u/mrichter/AliRoot.git] / PYTHIA6 / AliGenPythia.cxx

diff --git a/PYTHIA6/AliGenPythia.cxx b/PYTHIA6/AliGenPythia.cxx
index d48b00c4c31d69dcb6394083ebdfa7aa085ab2bd..9ad2f82d53ff994957a812ba12dc12d35380b7a1 100644 (file)
--- a/PYTHIA6/AliGenPythia.cxx
+++ b/PYTHIA6/AliGenPythia.cxx
@@ -110,7 +110,19 @@ AliGenPythia::AliGenPythia():
     fCountMode(kCountAll),      
     fHeader(0),  
     fRL(0),      
-    fFileName(0)
+    fFileName(0),
+    fFragPhotonInCalo(kFALSE),
+    fPi0InCalo(kFALSE) ,
+    fCheckEMCAL(kFALSE),
+    fCheckPHOS(kFALSE),
+    fFragPhotonOrPi0MinPt(0), 
+    fPHOSMinPhi(219.),
+    fPHOSMaxPhi(321.),
+    fPHOSEta(0.13),
+    fEMCALMinPhi(79.),
+    fEMCALMaxPhi(191.),
+    fEMCALEta(0.71)
+
 {
 // Default Constructor
   SetNuclei(0,0);
@@ -182,7 +194,18 @@ AliGenPythia::AliGenPythia(Int_t npart)
      fCountMode(kCountAll),      
      fHeader(0),  
      fRL(0),      
-     fFileName(0)  
+     fFileName(0),
+     fFragPhotonInCalo(kFALSE),
+     fPi0InCalo(kFALSE) ,
+     fCheckEMCAL(kFALSE),
+     fCheckPHOS(kFALSE),
+     fFragPhotonOrPi0MinPt(0),
+     fPHOSMinPhi(219.),
+     fPHOSMaxPhi(321.),
+     fPHOSEta(0.13),
+     fEMCALMinPhi(79.),
+     fEMCALMaxPhi(191.),
+     fEMCALEta(0.71)
 {
 // default charm production at 5. 5 TeV
 // semimuonic decay
@@ -198,76 +221,6 @@ AliGenPythia::AliGenPythia(Int_t npart)
     SetNuclei(0,0);
  }
 
-AliGenPythia::AliGenPythia(const AliGenPythia & Pythia)
-    :AliGenMC(Pythia),
-     fProcess(kPyCharm),          
-     fStrucFunc(kCTEQ5L), 
-     fEnergyCMS(5500.),
-     fKineBias(0.),
-     fTrials(0),
-     fTrialsRun(0),
-     fQ(0.),
-     fX1(0.),
-     fX2(0.),
-     fEventTime(0.),
-     fInteractionRate(0.),
-     fTimeWindow(0.),
-     fCurSubEvent(0),
-     fEventsTime(0),
-     fNev(0),
-     fFlavorSelect(0),
-     fXsection(0.),
-     fPythia(0),
-     fPtHardMin(0.),
-     fPtHardMax(1.e4),
-     fYHardMin(-1.e10),
-     fYHardMax(1.e10),
-     fGinit(kTRUE),
-     fGfinal(kTRUE),
-     fHadronisation(kTRUE),
-     fNpartons(0),
-     fReadFromFile(kFALSE),
-     fQuench(kFALSE),
-     fPtKick(1.),
-     fFullEvent(kTRUE),
-     fDecayer(new AliDecayerPythia()),
-     fDebugEventFirst(-1),
-     fDebugEventLast(-1),
-     fEtMinJet(0.),      
-     fEtMaxJet(1.e4),      
-     fEtaMinJet(-20.),     
-     fEtaMaxJet(20.),     
-     fPhiMinJet(0.),     
-     fPhiMaxJet(2.* TMath::Pi()),     
-     fJetReconstruction(kCell),
-     fEtaMinGamma(-20.),      
-     fEtaMaxGamma(20.),      
-     fPhiMinGamma(0.),      
-     fPhiMaxGamma(2. * TMath::Pi()),      
-     fPycellEtaMax(2.),     
-     fPycellNEta(274),       
-     fPycellNPhi(432),       
-     fPycellThreshold(0.),  
-     fPycellEtSeed(4.),     
-     fPycellMinEtJet(10.),  
-     fPycellMaxRadius(1.), 
-     fStackFillOpt(kFlavorSelection),   
-     fFeedDownOpt(kTRUE),    
-     fFragmentation(kTRUE),
-     fSetNuclei(kFALSE),
-     fNewMIS(kFALSE),   
-     fHFoff(kFALSE),    
-     fTriggerParticle(0),
-     fTriggerEta(0.9),     
-     fCountMode(kCountAll),      
-     fHeader(0),  
-     fRL(0),      
-     fFileName(0)  
-{
-// copy constructor
-    Pythia.Copy(*this);
-}
-
 AliGenPythia::~AliGenPythia()
 {
 // Destructor
@@ -476,6 +429,7 @@ void AliGenPythia::Init()
     case kPyMbMSEL1:
     case kPyJets:
     case kPyDirectGamma:
+    case kPyLhwgMb:    
        break;
     case kPyW:
     case kPyZ:
@@ -873,6 +827,36 @@ Int_t  AliGenPythia::GenerateMB()
        }
     }
 
+    // Select jets with fragmentation photon or pi0 going to PHOS or EMCAL
+    if (fProcess == kPyJets && (fFragPhotonInCalo || fPi0InCalo) ) {
+
+      Bool_t ok = kFALSE;
+
+      Int_t pdg  = 0; 
+      if (fFragPhotonInCalo) pdg = 22   ; // Photon
+      else if (fPi0InCalo) pdg = 111 ; // Pi0
+
+      for (i=0; i< np; i++) {
+       TParticle* iparticle = (TParticle *) fParticles->At(i);
+       if(iparticle->GetStatusCode()==1 && iparticle->GetPdgCode()==pdg && 
+          iparticle->Pt() > fFragPhotonOrPi0MinPt){
+         Int_t imother = iparticle->GetFirstMother();
+         TParticle* pmother = (TParticle *) fParticles->At(imother);
+         if(pdg == 111 || 
+            (pdg == 22 && pmother->GetStatusCode() != 11))//No photon from hadron decay
+           {
+             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;
+           }
+       }
+      }
+      if(!ok)
+       return 0;
+    }
+    
     if (fTriggerParticle) {
        Bool_t triggered = kFALSE;
        for (i = 0; i < np; i++) {
@@ -1198,14 +1182,6 @@ Bool_t AliGenPythia::CheckKinematicsOnChild(){
     return checking;
 }
 
-         
-AliGenPythia& AliGenPythia::operator=(const  AliGenPythia& rhs)
-{
-// Assignment operator
-    rhs.Copy(*this);
-    return *this;
-}
-
 void  AliGenPythia::LoadEvent(AliStack* stack, Int_t flag, Int_t reHadr)
 {
 //
@@ -1359,6 +1335,31 @@ void AliGenPythia::GetSubEventTime()
   return;
 }
 
+Bool_t AliGenPythia::IsInEMCAL(Float_t phi, Float_t eta)
+{
+  // Is particle in EMCAL acceptance? 
+  // phi in degrees, etamin=-etamax
+  if(phi > fEMCALMinPhi  && phi < fEMCALMaxPhi && 
+     eta < fEMCALEta  ) 
+    return kTRUE;
+  else 
+    return kFALSE;
+}
+
+Bool_t AliGenPythia::IsInPHOS(Float_t phi, Float_t eta)
+{
+  // Is particle in PHOS acceptance? 
+  // Acceptance slightly larger considered.
+  // phi in degrees, etamin=-etamax
+  if(phi > fPHOSMinPhi  && phi < fPHOSMaxPhi && 
+     eta < fPHOSEta  ) 
+    return kTRUE;
+  else 
+    return kFALSE;
+}
+
+
+
 #ifdef never
 void AliGenPythia::Streamer(TBuffer &R__b)
 {