Option for patching of the omega dalitz decay.

diff --git a/PYTHIA6/AliGenPythia.cxx b/PYTHIA6/AliGenPythia.cxx
index 6a1c784428fc9a93a8c02c47769939464751a514..ef26e76197b67448307e7f5c0b7ca3541707afd5 100644 (file)
--- a/PYTHIA6/AliGenPythia.cxx
+++ b/PYTHIA6/AliGenPythia.cxx
@@ -76,6 +76,7 @@ AliGenPythia::AliGenPythia():
     fGinit(1),
     fGfinal(1),
     fHadronisation(1),
+    fPatchOmegaDalitz(0), 
     fNpartons(0),
     fReadFromFile(0),
     fQuench(0),
@@ -172,6 +173,7 @@ AliGenPythia::AliGenPythia(Int_t npart)
      fGinit(kTRUE),
      fGfinal(kTRUE),
      fHadronisation(kTRUE),
+     fPatchOmegaDalitz(0), 
      fNpartons(0),
      fReadFromFile(kFALSE),
      fQuench(kFALSE),
@@ -623,10 +625,18 @@ void AliGenPythia::Generate()
 //
 // .. and perform hadronisation
 //     printf("Calling hadronisation %d\n", fPythia->GetN());
-           fPythia->Pyexec();  
+
+           if (fPatchOmegaDalitz) {
+             fPythia->SetMDCY(fPythia->Pycomp(111) ,1, 0);
+             fPythia->Pyexec();
+             fPythia->DalitzDecays();
+             fPythia->SetMDCY(fPythia->Pycomp(111) ,1, 1);
+           } 
+           fPythia->Pyexec();
        }
        fTrials++;
        fPythia->ImportParticles(&fParticles,"All");
+       
        if (TMath::Abs(fDyBoost) > 1.e-4) Boost();
 //
 //

diff --git a/PYTHIA6/AliGenPythia.h b/PYTHIA6/AliGenPythia.h
index 83eef4e4369d8f1e0dcf27756af2f6d5e9678501..662e6ed8d7f89fcdfe52e182184f6b774e30d3f9 100644 (file)
--- a/PYTHIA6/AliGenPythia.h
+++ b/PYTHIA6/AliGenPythia.h
@@ -147,6 +147,7 @@ class AliGenPythia : public AliGenMC
     void SetLength(Float_t length) {fLength = length;}
 
     virtual void SetHadronisation(Int_t flag = 1) {fHadronisation = flag;}
+    virtual void SetPatchOmegaDalitz(Int_t flag = 1) {fPatchOmegaDalitz = flag;}
     virtual void SetReadFromFile(const Text_t *filname) {fFileName = filname;  fReadFromFile = 1;}    
 
     //
@@ -222,6 +223,7 @@ class AliGenPythia : public AliGenMC
     Int_t       fGinit;             //initial state gluon radiation
     Int_t       fGfinal;            //final state gluon radiation
     Int_t       fHadronisation;     //hadronisation
+    Bool_t      fPatchOmegaDalitz;  //flag for omega dalitz decay patch
     Int_t       fNpartons;          //Number of partons before hadronisation
     Int_t       fReadFromFile;      //read partons from file
     Int_t       fQuench;            //Flag for quenching
@@ -297,12 +299,11 @@ class AliGenPythia : public AliGenMC
     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,9) // AliGenerator interface to Pythia
+    ClassDef(AliGenPythia, 9) // AliGenerator interface to Pythia
 };
 #endif
 

diff --git a/PYTHIA6/AliPythia.cxx b/PYTHIA6/AliPythia.cxx
index 35cd398f9f24a3020aef2af74f65d2feb9189e97..fe60a53d116c827d9744e1aefe9346f5bf8ee440 100644 (file)
--- a/PYTHIA6/AliPythia.cxx
+++ b/PYTHIA6/AliPythia.cxx
@@ -1,4 +1,3 @@
-
 /**************************************************************************
  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  *                                                                        *
@@ -20,7 +19,9 @@
 #include "AliPythiaRndm.h"
 #include "AliFastGlauber.h"
 #include "AliQuenchingWeights.h"
+#include "AliOmegaDalitz.h"
 #include "TVector3.h"
+#include "TLorentzVector.h"
 #include "PyquenCommon.h"
 
 ClassImp(AliPythia)
@@ -74,7 +75,8 @@ AliPythia::AliPythia():
     fZmax(0.97),
     fGlauber(0),
     fQuenchingWeights(0),
-    fItune(-1)
+    fItune(-1),
+    fOmegaDalitz()
 {
 // Default Constructor
 //
@@ -97,7 +99,8 @@ AliPythia::AliPythia(const AliPythia& pythia):
     fZmax(0.97),
     fGlauber(0),
     fQuenchingWeights(0),
-    fItune(-1)
+    fItune(-1),
+    fOmegaDalitz()
 {
     // Copy Constructor
     pythia.Copy(*this);
@@ -119,7 +122,6 @@ void AliPythia::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfun
     SetMDCY(Pycomp(310) ,1,0); // K0S
     SetMDCY(Pycomp(3122),1,0); // kLambda
     SetMDCY(Pycomp(3112),1,0); // sigma -
-    SetMDCY(Pycomp(3212),1,0); // sigma 0 
     SetMDCY(Pycomp(3222),1,0); // sigma +
     SetMDCY(Pycomp(3312),1,0); // xi - 
     SetMDCY(Pycomp(3322),1,0); // xi 0
@@ -608,7 +610,7 @@ void AliPythia::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfun
 //  
     SetMSTP(41,1);   // all resonance decays switched on
     Initialize("CMS","p","p",fEcms);
-    
+    fOmegaDalitz.Init();
 }
 
 Int_t AliPythia::CheckedLuComp(Int_t kf)
@@ -1453,3 +1455,28 @@ AliPythia& AliPythia::operator=(const  AliPythia& rhs)
     Fatal("Copy","Not implemented!\n");
 }
 
+void AliPythia::DalitzDecays()
+{
+
+  //
+  // Replace all omega dalitz decays with the correct matrix element decays
+  //
+  Int_t nt = fPyjets->N;
+  for (Int_t i = 0; i < nt; i++) {
+    if (fPyjets->K[1][i] != 223)                       continue;       
+    Int_t fd = fPyjets->K[3][i] - 1;
+    Int_t ld = fPyjets->K[4][i] - 1;
+    if (fd < 0)                                        continue;
+    if ((ld - fd) != 2)                                continue;
+    if ((fPyjets->K[1][fd] != 111) ||
+       (TMath::Abs(fPyjets->K[1][fd+1]) != 11))       continue;
+    TLorentzVector omega(fPyjets->P[0][i], fPyjets->P[1][i], fPyjets->P[2][i], fPyjets->P[3][i]);
+    fOmegaDalitz.Decay(223, &omega);
+    for (Int_t j = 0; j < 3; j++) {
+      for (Int_t k = 0; k < 4; k++) {
+       TLorentzVector vec = (fOmegaDalitz.Products())[2-j];
+       fPyjets->P[k][fd+j] = vec[k];
+      }
+    }
+  }
+}

diff --git a/PYTHIA6/AliPythia.h b/PYTHIA6/AliPythia.h
index 17e3df8939e03dee866cc3225cff93d5d1242ced..727c9cb83eaa37caa7740ccaa9ac78393d025e57 100644 (file)
--- a/PYTHIA6/AliPythia.h
+++ b/PYTHIA6/AliPythia.h
@@ -9,7 +9,7 @@
 #include <AliRndm.h>
 #include <AliStructFuncType.h>
 #include "PythiaProcesses.h"
-
+#include "AliOmegaDalitz.h"
 class AliFastGlauber;
 class AliQuenchingWeights;
 
@@ -45,6 +45,7 @@ class AliPythia : public TPythia6, public AliRndm
     // return instance of the singleton
     static  AliPythia* Instance();
     virtual void Quench();
+    void DalitzDecays();
     // Assignment Operator
     AliPythia & operator=(const AliPythia & rhs);
     void Copy(TObject&) const;
@@ -62,6 +63,7 @@ class AliPythia : public TPythia6, public AliRndm
     AliFastGlauber*       fGlauber;           //  ! The Glauber model
     AliQuenchingWeights*  fQuenchingWeights;  //  ! The Quenching Weights model
     Int_t                 fItune;             //  ! Pythia tune 
+    AliOmegaDalitz        fOmegaDalitz;       //  ! omega dalitz decayer
     static AliPythia*     fgAliPythia;        // Pointer to single instance
  private: 
     AliPythia();