- Change of coordinate system (x->-x, z->-z)

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

diff --git a/PYTHIA6/AliGenPythia.cxx b/PYTHIA6/AliGenPythia.cxx
index 3ddb173b969b7a1a2bd8bbb0aea095ac753c564f..549c6ca3066b9f1208020415b74b8f0d0c22e4fc 100644 (file)
--- a/PYTHIA6/AliGenPythia.cxx
+++ b/PYTHIA6/AliGenPythia.cxx
@@ -38,8 +38,10 @@
 #include "AliPythia.h"
 #include "AliPythiaRndm.h"
 #include "AliRun.h"
+#include "AliStack.h"
+#include "AliRunLoader.h"
 
- ClassImp(AliGenPythia)
+ClassImp(AliGenPythia)
 
 AliGenPythia::AliGenPythia()
                  :AliGenMC()
@@ -47,6 +49,7 @@ AliGenPythia::AliGenPythia()
 // Default Constructor
   fParticles = 0;
   fPythia    = 0;
+  fHeader = 0;
   fDecayer   = new AliDecayerPythia();
   SetEventListRange();
   SetJetPhiRange();
@@ -55,6 +58,8 @@ AliGenPythia::AliGenPythia()
   SetGammaPhiRange();
   SetGammaEtaRange();
   SetPtKick();
+  SetQuench();
+  
   fSetNuclei = kFALSE;
   if (!AliPythiaRndm::GetPythiaRandom()) 
     AliPythiaRndm::SetPythiaRandom(GetRandom());
@@ -83,7 +88,7 @@ AliGenPythia::AliGenPythia(Int_t npart)
     fFlavorSelect   = 0;
     // Produced particles  
     fParticles = new TClonesArray("TParticle",1000);
-    fEventVertex.Set(3);
+    fHeader = 0;
     SetEventListRange();
     SetJetPhiRange();
     SetJetEtaRange();
@@ -91,6 +96,7 @@ AliGenPythia::AliGenPythia(Int_t npart)
     SetGammaPhiRange();
     SetGammaEtaRange();
     SetJetReconstructionMode();
+    SetQuench();
     SetPtKick();
     // Options determining what to keep in the stack (Heavy flavour generation)
     fStackFillOpt = kFlavorSelection; // Keep particle with selected flavor
@@ -105,6 +111,7 @@ AliGenPythia::AliGenPythia(Int_t npart)
 }
 
 AliGenPythia::AliGenPythia(const AliGenPythia & Pythia)
+    :AliGenMC(Pythia)
 {
 // copy constructor
     Pythia.Copy(*this);
@@ -144,7 +151,7 @@ void AliGenPythia::Init()
 // Initialisation
     
     SetMC(AliPythia::Instance());
-    fPythia=(AliPythia*) fgMCEvGen;
+    fPythia=(AliPythia*) fMCEvGen;
 //
     fParentWeight=1./Float_t(fNpart);
 //
@@ -276,6 +283,7 @@ void AliGenPythia::Init()
 void AliGenPythia::Generate()
 {
 // Generate one event
+    
     fDecayer->ForceDecay();
 
     Float_t polar[3]   =   {0,0,0};
@@ -301,7 +309,19 @@ void AliGenPythia::Generate()
 //  event loop    
     while(1)
     {
+       if (fQuench) {
+           fPythia->SetMSTJ(1, 0);
+       }
+
        fPythia->Pyevnt();
+
+       if (fQuench) {
+           fPythia->Quench();
+           fPythia->SetMSTJ(1, 1);
+           fPythia->Pyexec();
+       }
+       
+       
        if (gAlice->GetEvNumber()>=fDebugEventFirst &&
            gAlice->GetEvNumber()<=fDebugEventLast) fPythia->Pylist(1);
        fTrials++;
@@ -321,9 +341,9 @@ void AliGenPythia::Generate()
        Float_t distz = iparticle->Vz()/10.;
        if(TMath::Abs(distz)>fCutVertexZ*fOsigma[2]) continue;
 //
-       fEventVertex[0] = iparticle->Vx()/10.+fOrigin.At(0);
-       fEventVertex[1] = iparticle->Vy()/10.+fOrigin.At(1);
-       fEventVertex[2] = iparticle->Vz()/10.+fOrigin.At(2);
+       fVertex[0] = iparticle->Vx()/10.+fOrigin.At(0);
+       fVertex[1] = iparticle->Vy()/10.+fOrigin.At(1);
+       fVertex[2] = iparticle->Vz()/10.+fOrigin.At(2);
 //
        Int_t* pParent   = new Int_t[np];
        Int_t* pSelected = new Int_t[np];
@@ -440,7 +460,7 @@ void AliGenPythia::Generate()
 // Track final state particle
                if (ks == 1) trackIt[i] = 1;
 // Track semi-stable particles
-               if ((ks ==1) || (fDecayer->GetLifetime(kf) > fMaxLifeTime))  trackIt[i] = 1;
+               if ((ks == 1) || (fDecayer->GetLifetime(kf) > fMaxLifeTime))  trackIt[i] = 1;
 // Track particles selected by process if undecayed. 
                if (fForceDecay == kNoDecay) {
                    if (ParentSelected(kf)) trackIt[i] = 1;
@@ -467,11 +487,11 @@ void AliGenPythia::Generate()
                    Float_t tof   = kconv*iparticle->T();
                    Int_t ipa     = iparticle->GetFirstMother()-1;
                    Int_t iparent = (ipa > -1) ? pParent[ipa] : -1;
-                   SetTrack(fTrackIt*trackIt[i] ,
+                   PushTrack(fTrackIt*trackIt[i] ,
                                     iparent, kf, p, origin, polar, tof, kPPrimary, nt, 1., ks);
                    pParent[i] = nt;
                    KeepTrack(nt); 
-               } //  SetTrack loop
+               } //  PushTrack loop
            }
        } else {
            nc = GenerateMB();
@@ -512,7 +532,7 @@ void AliGenPythia::Generate()
     } // event loop
     SetHighWaterMark(nt);
 //  adjust weight due to kinematic selection
-    AdjustWeights();
+//    AdjustWeights();
 //  get cross-section
     fXsection=fPythia->GetPARI(1);
 }
@@ -563,7 +583,7 @@ Int_t  AliGenPythia::GenerateMB()
            origin[1] = fOrigin[1]+iparticle->Vy()/10.;
            origin[2] = fOrigin[2]+iparticle->Vz()/10.;
            Float_t tof=kconv*iparticle->T();
-           SetTrack(fTrackIt*trackIt, iparent, kf, p, origin, polar,
+           PushTrack(fTrackIt*trackIt, iparent, kf, p, origin, polar,
                     tof, kPPrimary, nt, 1., ks);
            KeepTrack(nt);
            pParent[i] = nt;
@@ -615,16 +635,17 @@ void AliGenPythia::SetNuclei(Int_t a1, Int_t a2)
 void AliGenPythia::MakeHeader()
 {
 // Builds the event header, to be called after each event
-    AliGenEventHeader* header = new AliGenPythiaEventHeader("Pythia");
+    if (fHeader) delete fHeader;
+    fHeader = new AliGenPythiaEventHeader("Pythia");
 //
 // Event type  
-    ((AliGenPythiaEventHeader*) header)->SetProcessType(fPythia->GetMSTI(1));
+    ((AliGenPythiaEventHeader*) fHeader)->SetProcessType(fPythia->GetMSTI(1));
 //
 // Number of trials
-    ((AliGenPythiaEventHeader*) header)->SetTrials(fTrials);
+    ((AliGenPythiaEventHeader*) fHeader)->SetTrials(fTrials);
 //
 // Event Vertex 
-    header->SetPrimaryVertex(fEventVertex);
+    fHeader->SetPrimaryVertex(fVertex);
 //
 // Jets that have triggered
     if (fProcess == kPyJets)
@@ -634,11 +655,11 @@ void AliGenPythia::MakeHeader()
        GetJets(njet, ntrig, jets);
        
        for (Int_t i = 0; i < ntrig; i++) {
-           ((AliGenPythiaEventHeader*) header)->AddJet(jets[0][i], jets[1][i], jets[2][i], 
+           ((AliGenPythiaEventHeader*) fHeader)->AddJet(jets[0][i], jets[1][i], jets[2][i], 
                                                        jets[3][i]);
        }
     }
-    gAlice->SetGenEventHeader(header);
+    gAlice->SetGenEventHeader(fHeader);
 }
        
 
@@ -693,6 +714,7 @@ Bool_t AliGenPythia::CheckTrigger(TParticle* jet1, TParticle* jet2)
 AliGenPythia& AliGenPythia::operator=(const  AliGenPythia& rhs)
 {
 // Assignment operator
+    rhs.Copy(*this);
     return *this;
 }
 
@@ -702,12 +724,12 @@ void  AliGenPythia::LoadEvent()
 // Load event into Pythia Common Block
 //
  
-
-    Int_t npart = (Int_t) (gAlice->TreeK())->GetEntries(); 
+    AliRunLoader* rl = AliRunLoader::GetRunLoader();
+    Int_t npart = (rl->Stack())-> GetNprimary();
    (fPythia->GetPyjets())->N = npart;
 
     for (Int_t part = 0; part < npart; part++) {
-       TParticle *MPart = gAlice->Particle(part);
+       TParticle *MPart = (rl->Stack())->Particle(part);
        Int_t kf     = MPart->GetPdgCode();
        Int_t ks     = MPart->GetStatusCode();
        Float_t px = MPart->Px();
@@ -729,8 +751,7 @@ void  AliGenPythia::LoadEvent()
     }
 }
 
-void AliGenPythia::RecJetsUA1(Float_t eCellMin, Float_t eCellSeed, Float_t eMin, Float_t rMax, 
-                             Int_t& njets, Float_t jets [4][50])
+void AliGenPythia::RecJetsUA1(Int_t& njets, Float_t jets [4][50])
 {
 //
 //  Calls the Pythia jet finding algorithm to find jets in the current event
@@ -801,7 +822,7 @@ void  AliGenPythia::GetJets(Int_t& nJets, Int_t& nJetsTrig, Float_t jets[4][10])
 
        if (
            eta > fEtaMinJet && eta < fEtaMaxJet && 
-           phi > fPhiMinJet && eta < fPhiMaxJet &&
+           phi > fPhiMinJet && phi < fPhiMaxJet &&
            et  > fEtMinJet  && et  < fEtMaxJet     
            ) 
        {