Added file check.

[u/mrichter/AliRoot.git] / JETAN / AliJetParticlesReaderKine.cxx

diff --git a/JETAN/AliJetParticlesReaderKine.cxx b/JETAN/AliJetParticlesReaderKine.cxx
index ff111d3e3a91b95d3e09f61a6485826dfb1f36a4..1fb6b9533d74637743b8e935fcf250fdfb0ca44d 100644 (file)
--- a/JETAN/AliJetParticlesReaderKine.cxx
+++ b/JETAN/AliJetParticlesReaderKine.cxx
@@ -12,20 +12,22 @@
 /////////////////////////////////////////////////////////////////////////
 
 #include <Riostream.h>
 /////////////////////////////////////////////////////////////////////////
 
 #include <Riostream.h>

+#include <TFile.h>

 #include <TString.h>
 #include <TParticle.h>
 #include <TString.h>
 #include <TParticle.h>

+#include <TLorentzVector.h>

 #include <AliRunLoader.h>
 #include <AliStack.h>
 #include <AliRunLoader.h>
 #include <AliStack.h>

-
+#include <AliHeader.h>
+#include <AliGenEventHeader.h>
+#include <AliGenPythiaEventHeader.h>
+#include <AliGenHijingEventHeader.h>

 #include "AliJetParticle.h"
 #include "AliJetEventParticles.h"
 #include "AliJetParticlesReaderKine.h"
 
 ClassImp(AliJetParticlesReaderKine)
 
 #include "AliJetParticle.h"
 #include "AliJetEventParticles.h"
 #include "AliJetParticlesReaderKine.h"
 
 ClassImp(AliJetParticlesReaderKine)
 

-
-/**********************************************************/
-

 AliJetParticlesReaderKine::AliJetParticlesReaderKine() :
   AliJetParticlesReader(),
   fFileName("galice.root"),
 AliJetParticlesReaderKine::AliJetParticlesReaderKine() :
   AliJetParticlesReader(),
   fFileName("galice.root"),


@@ -36,9 +38,6 @@ AliJetParticlesReaderKine::AliJetParticlesReaderKine() :
   //constructor
 }
 
   //constructor
 }
 

-/**********************************************************/
-
-

 AliJetParticlesReaderKine::AliJetParticlesReaderKine(TString& fname) :
   AliJetParticlesReader(),
   fFileName(fname),
 AliJetParticlesReaderKine::AliJetParticlesReaderKine(TString& fname) :
   AliJetParticlesReader(),
   fFileName(fname),


@@ -49,8 +48,6 @@ AliJetParticlesReaderKine::AliJetParticlesReaderKine(TString& fname) :
   //constructor
 }
 
   //constructor
 }
 

-/**********************************************************/
-

 AliJetParticlesReaderKine::AliJetParticlesReaderKine(TObjArray* dirs, const Char_t *filename):
   AliJetParticlesReader(dirs),
   fFileName(filename),
 AliJetParticlesReaderKine::AliJetParticlesReaderKine(TObjArray* dirs, const Char_t *filename):
   AliJetParticlesReader(dirs),
   fFileName(filename),


@@ -61,27 +58,21 @@ AliJetParticlesReaderKine::AliJetParticlesReaderKine(TObjArray* dirs, const Char
   //constructor
 }
 
   //constructor
 }
 

-/**********************************************************/
-

 AliJetParticlesReaderKine::~AliJetParticlesReaderKine()
 {
   //destructor
   if(fRunLoader) delete fRunLoader;
 }
 
 AliJetParticlesReaderKine::~AliJetParticlesReaderKine()
 {
   //destructor
   if(fRunLoader) delete fRunLoader;
 }
 

-/**********************************************************/
-

 void AliJetParticlesReaderKine::Rewind()
 {
   //Rewinds to the beginning
   if(fRunLoader) delete fRunLoader;
 void AliJetParticlesReaderKine::Rewind()
 {
   //Rewinds to the beginning
   if(fRunLoader) delete fRunLoader;

-  fRunLoader = 0;
-  fCurrentDir = 0;
-  fNEventsRead= 0;  
+  fRunLoader   = 0;
+  fCurrentDir  = 0;
+  fNEventsRead = 0;  

 }
 
 }
 

-/**********************************************************/
-

 Int_t AliJetParticlesReaderKine::ReadNext()
 {
   //Reads Kinematics Tree
 Int_t AliJetParticlesReaderKine::ReadNext()
 {
   //Reads Kinematics Tree


@@ -91,8 +82,11 @@ Int_t AliJetParticlesReaderKine::ReadNext()
 
   while(fCurrentDir < GetNumberOfDirs())
     { 
 
   while(fCurrentDir < GetNumberOfDirs())
     { 

-      if (!OpenFile(fCurrentDir))
+
+      if (!OpenFile(fCurrentDir)) 

       { 
       { 

+       delete fRunLoader; //close current session
+       fRunLoader = 0;    //assure pointer is null

        fCurrentDir++;
        continue;
       }
        fCurrentDir++;
        continue;
       }


@@ -106,28 +100,118 @@ Int_t AliJetParticlesReaderKine::ReadNext()
          continue; 
        }
      
          continue; 
        }
      

-      //Info("ReadNext","Reading Event %d",fCurrentEvent);
-
+      Info("ReadNext","Reading Event %d",fCurrentEvent);

       fRunLoader->GetEvent(fCurrentEvent);
       AliStack* stack = fRunLoader->Stack();
       if (!stack)
        {
       fRunLoader->GetEvent(fCurrentEvent);
       AliStack* stack = fRunLoader->Stack();
       if (!stack)
        {

-         Error("ReadNext","Can not get stack for event %d",fCurrentEvent);
+         Error("ReadNext","Can't get stack for event %d",fCurrentEvent);

          continue;
        }
 
          continue;
        }
 

+      //clear old event
+      Int_t nprim = stack->GetNprimary();
+      Int_t npart = nprim;
+      if(fUseTracks)
+       npart = stack->GetNtrack();
+      fEventParticles->Reset(npart);
+
+      TString headdesc="";
+      Int_t gentype=-1;
+      AliHeader *header=fRunLoader->GetHeader();
+      if(!header) {
+       Warning("ReadNext","Header not found in event %d",fCurrentEvent);
+      } else {
+       AliGenEventHeader *genheader=header->GenEventHeader();
+       TString hname=genheader->GetName();
+       if(hname.CompareTo("HIJINGparam")==0) {
+         gentype=1;      ;
+         headdesc+="HIJINGparam";
+       } else if(hname.CompareTo("Hijing")==0) {
+         gentype=2;
+         AliGenHijingEventHeader *hheader=(AliGenHijingEventHeader*)header->GenEventHeader();
+         headdesc+="Hijing";
+         if(!hheader) {
+           Warning("ReadNext","Hijing-Header not found in event %d",fCurrentEvent);
+         } else {
+           Int_t ntrials=hheader->Trials();
+           fEventParticles->SetTrials(ntrials);
+           fEventParticles->SetImpact(hheader->ImpactParameter());
+           fEventParticles->SetNhard(hheader->HardScatters());
+           fEventParticles->SetNpart(hheader->NwNw());
+
+           TLorentzVector  jet1; 
+           TLorentzVector  jet2; 
+           TLorentzVector  jetFsr1;
+           TLorentzVector  jetFsr2;
+           hheader->GetJets(jet1,jet2,jetFsr1,jetFsr2);
+           fEventParticles->AddHard(0,jet1.X(),jet1.Y(),jet1.Z(),jet1.T(),0);
+           fEventParticles->AddHard(1,jet2.X(),jet2.Y(),jet2.Z(),jet2.T(),0);
+         }
+       } else if(hname.CompareTo("Pythia")==0){
+         gentype=3;
+         headdesc+="Pythia";
+         AliGenPythiaEventHeader *pheader=(AliGenPythiaEventHeader*)header->GenEventHeader();
+         if(!pheader) {
+           Warning("ReadNext","Pythia-Header not found in event %d",fCurrentEvent);
+         } else {
+           Int_t ntruq=0;
+#ifndef NOUQHEADERINFO
+           ntruq=pheader->NUQTriggerJets();
+           if(ntruq){
+             Double_t x0=pheader->GetXJet();
+             Double_t y0=pheader->GetYJet();
+             if(x0==y0==-1){
+               x0=y0=0.;
+             }
+             Double_t zquench[4];
+             pheader->GetZQuench(zquench);
+             fEventParticles->SetXYJet(x0,y0);
+             fEventParticles->SetZQuench(zquench);
+             for(Int_t j=0;j<ntruq;j++){
+               Float_t pjet[4];
+               pheader->UQJet(j,pjet);
+               fEventParticles->AddUQJet(pjet);
+             }
+           }
+#endif
+           //Int_t ptyp=pheader->ProcessType();
+           Int_t ntrials=pheader->Trials();
+           fEventParticles->SetTrials(ntrials);
+
+           Int_t ntr=pheader->NTriggerJets();
+           if(ntr){
+             for(Int_t j=0;j<ntr;j++){
+               Float_t pjet[4];
+               pheader->TriggerJet(j,pjet);
+               fEventParticles->AddJet(pjet);
+               if(!ntruq) fEventParticles->AddUQJet(pjet);
+             }
+           }
+           for(Int_t i=6;i<=7;i++){
+             TParticle *MP = stack->Particle(i);
+             if(!MP) break;
+             Int_t type=0;
+             if(MP->GetPdgCode()==21) type=2;
+             else type=1;
+             fEventParticles->AddHard(i-6,MP->Px(),MP->Py(),MP->Pz(),MP->Energy(),type);
+           }
+         }
+       }
+       headdesc+=" Run ";
+       headdesc+=header->GetRun();
+       headdesc+=": Ev ";
+       headdesc+=header->GetEventNrInRun();
+      }
+      fEventParticles->SetHeader(headdesc);
+

       //get vertex
       const TParticle *kv = stack->Particle(0);
       if(kv) {
        fEventParticles->SetVertex(kv->Vx(),kv->Vy(),kv->Vz());
       }
 
       //get vertex
       const TParticle *kv = stack->Particle(0);
       if(kv) {
        fEventParticles->SetVertex(kv->Vx(),kv->Vy(),kv->Vz());
       }
 

-      Int_t nprim = stack->GetNprimary();
-      Int_t npart = nprim;
-      if(fUseTracks)
-       npart = stack->GetNtrack();
-
-      fEventParticles->Reset(npart);
+      //loop over particles

       for (Int_t i = 0;i<npart; i++)
        {
          TParticle *p = stack->Particle(i);
       for (Int_t i = 0;i<npart; i++)
        {
          TParticle *p = stack->Particle(i);


@@ -135,43 +219,61 @@ Int_t AliJetParticlesReaderKine::ReadNext()
          Int_t child1 = p->GetFirstDaughter();
          //Int_t child2 = p->GetLastDaughter();        
          //Int_t mother = p->GetFirstMother();    
          Int_t child1 = p->GetFirstDaughter();
          //Int_t child2 = p->GetLastDaughter();        
          //Int_t mother = p->GetFirstMother();    

-         if((child1>=0) && (child1<nprim)) continue; 

          //cout << child1 << " " << child2 << " " << mother << endl;
          //cout << child1 << " " << child2 << " " << mother << endl;

+         if((child1>=0) && (child1<nprim)) continue; 

 
 

-         if(IsAcceptedParticle(p)) //put particle in event
+         //check status code depending on gentype
+         if(gentype==1){ // Hijing Param
+           if(p->GetStatusCode()!=0) continue;
+         } else if(gentype==2){ //Hijing
+           //if( p->GetStatusCode()!=0){
+           //  cout <<p->GetStatusCode() << " ";p->Print(); 
+           //}
+         } else if(gentype==3){ //Pythia
+           if(p->GetStatusCode()!=1) continue;
+           p->SetWeight(-123); //mark particle
+         }
+
+         //kinematic cuts
+         if(IsAcceptedParticle(p)){ //put particle in event
+           //if(p->Pt()>20){cout <<p->GetStatusCode() << " ";p->Print(); }

            fEventParticles->AddParticle(p,i); 
            fEventParticles->AddParticle(p,i); 

+         }

        }
       fCurrentEvent++;
        }
       fCurrentEvent++;

-      fCurrentDir++;

       fNEventsRead++;
       return kTRUE;
     }
       fNEventsRead++;
       return kTRUE;
     }

-      //end of loop over directories specified in fDirs Obj Array
+
+  //end of loop over directories specified in fDirs ObjArray

   return kFALSE;
 }
 
   return kFALSE;
 }
 

-/**********************************************************/
-

 Int_t AliJetParticlesReaderKine::OpenFile(Int_t n)
 {
   //opens file with kine tree
 
 Int_t AliJetParticlesReaderKine::OpenFile(Int_t n)
 {
   //opens file with kine tree
 

+  if(fRunLoader){
+    if(fCurrentEvent < fRunLoader->GetNumberOfEvents()) return kTRUE;
+    else return kFALSE;
+  }
+

   const TString& dirname = GetDirName(n);
   if (dirname == "")
     { 
   const TString& dirname = GetDirName(n);
   if (dirname == "")
     { 

-      Error("OpenNextFile","Can not get directory name with index %d",n);
+      Error("OpenNextFile","Can't get directory name with index %d",n);

       return kFALSE;
     }
 
   TString filename = dirname +"/"+ fFileName;
       return kFALSE;
     }
 
   TString filename = dirname +"/"+ fFileName;

-  if(fRunLoader) delete fRunLoader;
-  fRunLoader = AliRunLoader::Open(filename.Data()); 
-
-  if (fRunLoader == 0)
+  TFile file(filename);
+  if ( file.IsOpen() == 0)

     {
       Error("OpenNextFile","Can't open session from file %s",filename.Data());
       return kFALSE;
     }
     {
       Error("OpenNextFile","Can't open session from file %s",filename.Data());
       return kFALSE;
     }

+  file.Close();
+  fRunLoader = AliRunLoader::Open(filename.Data()); 

   
   if (fRunLoader->GetNumberOfEvents() <= 0)
     {
   
   if (fRunLoader->GetNumberOfEvents() <= 0)
     {


@@ -180,13 +282,13 @@ Int_t AliJetParticlesReaderKine::OpenFile(Int_t n)
       fRunLoader = 0;
       return kFALSE;
     }
       fRunLoader = 0;
       return kFALSE;
     }

-  
+

   if (fRunLoader->LoadKinematics())
     {
       Error("OpenNextFile","Error occured while loading kinematics.");
       return kFALSE;
     }
   if (fRunLoader->LoadKinematics())
     {
       Error("OpenNextFile","Error occured while loading kinematics.");
       return kFALSE;
     }

-  
+

   fCurrentEvent = 0;
   return kTRUE;
 }
   fCurrentEvent = 0;
   return kTRUE;
 }