If gamma in eta acceptance rotate event into PHOS phi acceptance.
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 15 Oct 2007 07:27:06 +0000 (07:27 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 15 Oct 2007 07:27:06 +0000 (07:27 +0000)
(R. Diaz Valdes)

PYTHIA6/AliGenPythia.cxx
PYTHIA6/AliGenPythia.h

index 043b696..d8238ce 100644 (file)
@@ -114,9 +114,12 @@ AliGenPythia::AliGenPythia():
     fFileName(0),
     fFragPhotonInCalo(kFALSE),
     fPi0InCalo(kFALSE) ,
+    fPhotonInCalo(kFALSE),
     fCheckEMCAL(kFALSE),
     fCheckPHOS(kFALSE),
+    fCheckPHOSeta(kFALSE),
     fFragPhotonOrPi0MinPt(0), 
+    fPhotonMinPt(0), 
     fPHOSMinPhi(219.),
     fPHOSMaxPhi(321.),
     fPHOSEta(0.13),
@@ -198,9 +201,12 @@ AliGenPythia::AliGenPythia(Int_t npart)
      fFileName(0),
      fFragPhotonInCalo(kFALSE),
      fPi0InCalo(kFALSE) ,
+     fPhotonInCalo(kFALSE),
      fCheckEMCAL(kFALSE),
      fCheckPHOS(kFALSE),
+     fCheckPHOSeta(kFALSE),
      fFragPhotonOrPi0MinPt(0),
+     fPhotonMinPt(0),
      fPHOSMinPhi(219.),
      fPHOSMaxPhi(321.),
      fPHOSEta(0.13),
@@ -870,6 +876,41 @@ Int_t  AliGenPythia::GenerateMB()
        return 0;
     }
     
+    
+     // 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)
+       return 0;
+    }
+    
     if (fTriggerParticle) {
        Bool_t triggered = kFALSE;
        for (i = 0; i < np; i++) {
@@ -1369,6 +1410,56 @@ Bool_t AliGenPythia::IsInPHOS(Float_t phi, Float_t eta)
     return kFALSE;
 }
 
+void AliGenPythia::RotatePhi(Int_t iphcand, Bool_t& okdd)
+{
+  //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 = gRandom->Uniform(phiPHOSmin,phiPHOSmax);
+  
+  //calculate deltaphi
+  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;
+  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);
+      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;
+}
 
 
 #ifdef never
index e4dfabb..460c282 100644 (file)
@@ -39,87 +39,104 @@ class AliGenPythia : public AliGenMC
     virtual ~AliGenPythia();
     virtual void    Generate();
     virtual void    Init();
-    // set a cut on the Z coord. of the primary vertex (cm)
-    //
+    // Range of events to be printed
     virtual void    SetEventListRange(Int_t eventFirst=-1, Int_t eventLast=-1);
-    // select process type
+    // Select process type
     virtual void    SetProcess(Process_t proc = kPyCharm) {fProcess = proc;}
-    // select structure function
+    // Select structure function
     virtual void    SetStrucFunc(StrucFunc_t func =  kCTEQ5L) {fStrucFunc = func;}
-    // select pt of hard scattering 
+    // Select pt of hard scattering 
     virtual void    SetPtHard(Float_t ptmin = 0, Float_t ptmax = 1.e10)
        {fPtHardMin = ptmin; fPtHardMax = ptmax; }
+    // y of hard scattering
     virtual void    SetYHard(Float_t ymin = -1.e10, Float_t ymax = 1.e10)
        {fYHardMin = ymin; fYHardMax = ymax; }
     // Set initial and final state gluon radiation
     virtual void    SetGluonRadiation(Int_t iIn, Int_t iFin)
        {fGinit = iIn; fGfinal = iFin;}
+    // Intrinsic kT
     virtual void    SetPtKick(Float_t kt = 1.)
        {fPtKick = kt;}
     // Use the Pythia 6.3 new multiple interations scenario
     virtual void    UseNewMultipleInteractionsScenario() {fNewMIS = kTRUE;}
     // Switch off heavy flavors
     virtual void    SwitchHFOff() {fHFoff = kTRUE;}
-    // set centre of mass energy
+    // Set centre of mass energy
     virtual void    SetEnergyCMS(Float_t energy = 5500) {fEnergyCMS = energy;}
-    // treat protons as inside nuclei
+    // Treat protons as inside nuclei with mass numbers a1 and a2
     virtual void    SetNuclei(Int_t a1, Int_t a2);
+    //
+    // Trigger options
+    //
+    // Energy range for jet trigger
     virtual void    SetJetEtRange(Float_t etmin = 0., Float_t etmax = 1.e4)
        {fEtMinJet = etmin; fEtMaxJet = etmax;}
+    // Eta range for jet trigger
     virtual void    SetJetEtaRange(Float_t etamin = -20., Float_t etamax = 20.)
        {fEtaMinJet = etamin; fEtaMaxJet = etamax;}
-    virtual void    SetJetReconstructionMode(Int_t mode = kCell) {fJetReconstruction = mode;}
+    // Phi range for jet trigger
     virtual void    SetJetPhiRange(Float_t phimin = 0., Float_t phimax = 360.)
        {fPhiMinJet = TMath::Pi()*phimin/180.; fPhiMaxJet = TMath::Pi()*phimax/180.;}
+    // Jet reconstruction mode; default is cone algorithm
+    virtual void    SetJetReconstructionMode(Int_t mode = kCell) {fJetReconstruction = mode;}
+    // Eta range for gamma trigger 
     virtual void    SetGammaEtaRange(Float_t etamin = -20., Float_t etamax = 20.)
        {fEtaMinGamma = etamin; fEtaMaxGamma = etamax;}
+    // 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  SetPi0InCalo       (Bool_t b)      {fPi0InCalo = 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  SetPi0InEMCAL(Bool_t b) 
-    {fCheckEMCAL = b ;  fPi0InCalo = b ;}
-    virtual void  SetPi0InPHOS(Bool_t b) 
-    {fCheckPHOS = b ;   fPi0InCalo = b ;}
-    
-    virtual void  SetFragPhotonOrPi0MinPt(Float_t pt)      {fFragPhotonOrPi0MinPt=pt;}
-  
-    // Trigger on a particle
+    virtual void  SetFragPhotonInCalo(Bool_t b)  {fFragPhotonInCalo = 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  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  SetPhotonInPHOS(Bool_t b)      {fCheckPHOS    = b; fPhotonInCalo     = b;}
+    virtual void  SetPhotonInPHOSeta(Bool_t b)   {fCheckPHOSeta = b; fPhotonInCalo     = b;}
+    virtual void  SetFragPhotonOrPi0MinPt(Float_t pt)      {fFragPhotonOrPi0MinPt = pt;}
+    virtual void  SetPhotonMinPt(Float_t pt)     {fPhotonMinPt = pt;}
+    // Trigger and rotate event 
+    void RotatePhi(Int_t iphcand, Bool_t& okdd);
+    // Trigger on a single particle
     virtual void    SetTriggerParticle(Int_t particle = 0, Float_t etamax = 0.9) 
        {fTriggerParticle = particle; fTriggerEta = etamax;}
+    //
+    // Heavy flavor options
+    //
     // Set option for feed down from higher family
     virtual void SetFeedDownHigherFamily(Bool_t opt) {
-      fFeedDownOpt = opt;
+       fFeedDownOpt = opt;
     }
     // Set option for selecting particles kept in stack according to flavor
     // or to parent selection
     virtual void SetStackFillOpt(StackFillOpt_t opt) {
-      fStackFillOpt = opt;
+       fStackFillOpt = opt;
     }
     // Set fragmentation option
     virtual void SetFragmentation(Bool_t opt) {
-      fFragmentation = opt;
+       fFragmentation = opt;
     }
     // Set counting mode
     virtual void SetCountMode(CountMode_t mode) {
-      fCountMode = mode;
+       fCountMode = mode;
     }
+    //
+    // Quenching
+    //
     // Set quenching mode 0 = no, 1 = AM, 2 = IL
     virtual void SetQuench(Int_t flag = 0) {fQuench = flag;}
     virtual void SetHadronisation(Int_t flag = 1) {fHadronisation = flag;}
     virtual void SetReadFromFile(const Text_t *filname) {fFileName = filname;  fReadFromFile = 1;}    
-  // Get interaction rate for pileup studies
+
+    //
+    // Pile-up
+    //
+    // Get interaction rate for pileup studies
     virtual void    SetInteractionRate(Float_t rate,Float_t timewindow = 90.e-6);
     virtual Float_t GetInteractionRate() const {return fInteractionRate;}
     // get cross section of process
@@ -154,11 +171,8 @@ class AliGenPythia : public AliGenMC
     //
     virtual void FinishRun();
     Bool_t CheckTrigger(TParticle* jet1, TParticle* jet2);
-
     //Used in some processes to selected child properties
     Bool_t CheckKinematicsOnChild();
-    
-
     void     GetSubEventTime();
 
  protected:
@@ -223,43 +237,47 @@ class AliGenPythia : public AliGenMC
                                     // parents and their decays
     Bool_t fFeedDownOpt;            // Option to set feed down from higher
                                     // quark families (e.g. b->c)
-    Bool_t  fFragmentation;          // Option to activate fragmentation by Pythia
-    Bool_t  fSetNuclei;              // Flag indicating that SetNuclei has been called
-    Bool_t  fNewMIS;                 // Flag for the new multipple interactions scenario
-    Bool_t  fHFoff;                  // Flag for switching heafy flavor production off
-    Int_t   fTriggerParticle;        // Trigger on this particle ...
-    Float_t fTriggerEta;             // .. within |eta| < fTriggerEta
-    //
-    CountMode_t fCountMode;            // Options for counting when the event will be finished.
-    AliGenPythiaEventHeader* fHeader;  //! Event header
-    AliRunLoader*            fRL;      //! Run Loader
-    const Text_t* fFileName;           //! Name of file to read from
-     
+    Bool_t  fFragmentation;         // Option to activate fragmentation by Pythia
+    Bool_t  fSetNuclei;             // Flag indicating that SetNuclei has been called
+    Bool_t  fNewMIS;                // Flag for the new multipple interactions scenario
+    Bool_t  fHFoff;                 // Flag for switching heafy flavor production off
+    Int_t   fTriggerParticle;       // Trigger on this particle ...
+    Float_t fTriggerEta;            // .. within |eta| < fTriggerEta
+    CountMode_t fCountMode;         // Options for counting when the event will be finished.     
     // fCountMode = kCountAll         --> All particles that end up in the
     //                                    stack are counted
     // fCountMode = kCountParents     --> Only selected parents are counted
     // fCountMode = kCountTrackabless --> Only particles flagged for tracking
     //                                     are counted
     //
+    //
+
+    AliGenPythiaEventHeader* fHeader;  //! Event header
+    AliRunLoader*            fRL;      //! Run Loader
+    const Text_t* fFileName;           //! Name of file to read from
+
 
-    Bool_t fFragPhotonInCalo ; // Option ask for Fragmentation Photon in calorimeters acceptance
-    Bool_t fPi0InCalo ; // Option ask for Pi0 in calorimeters acceptance
-    Bool_t fCheckEMCAL ; //  Option ask for FragPhoton or Pi0 in calorimeters EMCAL acceptance
-    Bool_t fCheckPHOS ; //  Option ask for FragPhoton or Pi0 in calorimeters PHOS acceptance
-    Float_t fFragPhotonOrPi0MinPt ; // Minimum momentum of Fragmentation Photon or Pi0 
+    Bool_t fFragPhotonInCalo; // Option to ask for Fragmentation Photon 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 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
+    Float_t fFragPhotonOrPi0MinPt; // Minimum momentum of Fragmentation Photon or Pi0
+    Float_t fPhotonMinPt;          // Minimum momentum of Photon 
     //Calorimeters eta-phi acceptance 
-    Float_t fPHOSMinPhi ;
-    Float_t fPHOSMaxPhi ;
-    Float_t fPHOSEta ;
-    Float_t fEMCALMinPhi ;
-    Float_t fEMCALMaxPhi ;
-    Float_t fEMCALEta ;
+    Float_t fPHOSMinPhi;           // Minimum phi PHOS
+    Float_t fPHOSMaxPhi;           // Maximum phi PHOS
+    Float_t fPHOSEta;              // Minimum eta PHOS
+    Float_t fEMCALMinPhi;          // Minimum phi EMCAL
+    Float_t fEMCALMaxPhi;          // Maximum phi EMCAL
+    Float_t fEMCALEta;             // Maximum eta EMCAL
 
  private:
     AliGenPythia(const AliGenPythia &Pythia);
     AliGenPythia & operator=(const AliGenPythia & rhs);
 
-    ClassDef(AliGenPythia,5) // AliGenerator interface to Pythia
+    ClassDef(AliGenPythia,6) // AliGenerator interface to Pythia
 };
 #endif