EVGEN/AliGenReaderEcalJets.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id$ */
  17 //
  18 // Realisation of AliGenReader to be used with AliGenExtFile
  19 // It reads Pythia Jet events from a ntuple like event structure.
  20 // The event format is defined in Init()
  21 // NextEvent() is used to loop over events and NextParticle() to loop over particles.
  22 // Author: andreas.morsch@cern.ch
  23 //
  24 #include <TFile.h>
  25 #include <TParticle.h>
  26 #include <TDatabasePDG.h>
  27 #include <TTree.h>
  28
  29 #include "AliGenReaderEcalJets.h"
  30
  31 ClassImp(AliGenReaderEcalJets)
  32
  33
  34 AliGenReaderEcalJets::AliGenReaderEcalJets():
  35     fNcurrent(0),
  36     fNparticle(0),
  37     fNev(0),
  38     fNpart(0),
  39     fNjet(0),
  40     fNsjet(0),
  41     fNpjet(0),
  42     fTreeNtuple(0)
  43 {
  44 // Default constructor
  45     for (Int_t i = 0; i < 200; i++) {
  46         if (i < 2) {
  47             fX[i]     = 0.;
  48             fXtyp[i]  = 0;
  49         } else if (i < 10) {
  50             fJet[i]   = 0.;
  51             fJeta[i]  = 0.;
  52             fJphi[i]  = 0.;
  53             fJset[i]  = 0.;
  54             fJseta[i] = 0.;
  55             fJsphi[i] = 0.;
  56             fJpet[i]  = 0.;
  57             fJpeta[i] = 0.;
  58             fJpphi[i] = 0.;
  59         } else {
  60             fXpt[i]  = 0.;
  61             fXeta[i] = 0.;
  62             fXphi[i] = 0.;
  63             fXid[i]  = 0;
  64         }
  65     }
  66 }
  67
  68  AliGenReaderEcalJets::AliGenReaderEcalJets(const AliGenReaderEcalJets &reader)
  69      :AliGenReader(reader),
  70       fNcurrent(0),
  71       fNparticle(0),
  72       fNev(0),
  73       fNpart(0),
  74       fNjet(0),
  75       fNsjet(0),
  76       fNpjet(0),
  77       fTreeNtuple(0)
  78 {
  79     for (Int_t i = 0; i < 200; i++) {
  80         if (i < 2) {
  81             fX[i]     = 0.;
  82             fXtyp[i]  = 0;
  83         } else if (i < 10) {
  84             fJet[i]   = 0.;
  85             fJeta[i]  = 0.;
  86             fJphi[i]  = 0.;
  87             fJset[i]  = 0.;
  88             fJseta[i] = 0.;
  89             fJsphi[i] = 0.;
  90             fJpet[i]  = 0.;
  91             fJpeta[i] = 0.;
  92             fJpphi[i] = 0.;
  93         } else {
  94             fXpt[i]  = 0.;
  95             fXeta[i] = 0.;
  96             fXphi[i] = 0.;
  97             fXid[i]  = 0;
  98         }
  99     }
 100
 101     // Copy Constructor
 102     reader.Copy(*this);
 103 }
 104
 105 void AliGenReaderEcalJets::Init()
 106 {
 107 //
 108 // reset the existing file environment and open a new root file if
 109 // the pointer to the Fluka tree is null
 110
 111     TFile *pFile=0;
 112     if (!pFile) {
 113         pFile = new TFile(fFileName);
 114         pFile->cd();
 115         printf("\n I have opened %s file \n", fFileName);
 116     }
 117 // get the tree address in the Fluka boundary source file
 118     fTreeNtuple = (TTree*)gDirectory->Get("h1");
 119     TTree *h1=fTreeNtuple;
 120     h1->SetMakeClass(1);
 121 //Set branch addresses
 122     h1->SetBranchAddress("nev",   &fNev);
 123     h1->SetBranchAddress("x",      fX);
 124     h1->SetBranchAddress("xtyp",   fXtyp);
 125     h1->SetBranchAddress("npart", &fNpart);
 126     h1->SetBranchAddress("xpt",    fXpt);
 127     h1->SetBranchAddress("xeta",   fXeta);
 128     h1->SetBranchAddress("xphi",   fXphi);
 129     h1->SetBranchAddress("xid",    fXid);
 130     h1->SetBranchAddress("njet",  &fNjet);
 131     h1->SetBranchAddress("jet",    fJet);
 132     h1->SetBranchAddress("jeta",   fJeta);
 133     h1->SetBranchAddress("jphi",   fJphi);
 134     h1->SetBranchAddress("nsjet", &fNsjet);
 135     h1->SetBranchAddress("jset",   fJset);
 136     h1->SetBranchAddress("jseta",  fJseta);
 137     h1->SetBranchAddress("jsphi",  fJsphi);
 138     h1->SetBranchAddress("npjet", &fNpjet);
 139     h1->SetBranchAddress("jpet",   fJpet);
 140     h1->SetBranchAddress("jpeta",  fJpeta);
 141     h1->SetBranchAddress("jpphi",  fJpphi);
 142 }
 143
 144 Int_t AliGenReaderEcalJets::NextEvent()
 145 {
 146 // Read the next event
 147     Int_t nTracks=0, nread=0;
 148
 149     TFile* pFile = fTreeNtuple->GetCurrentFile();
 150     pFile->cd();
 151
 152     Int_t nentries = (Int_t) fTreeNtuple->GetEntries();
 153     if (fNcurrent < nentries) {
 154         Int_t nb = (Int_t)fTreeNtuple->GetEvent(fNcurrent);
 155         nread += nb;
 156         fNcurrent++;
 157         printf("\n Next event contains %d tracks! \n", fNpjet);
 158         nTracks    = fNpjet;
 159         fNparticle = 0;
 160         return nTracks;
 161     }
 162     return 0;
 163 }
 164
 165 TParticle* AliGenReaderEcalJets::NextParticle()
 166 {
 167 // Read the next particle
 168
 169     Float_t p[4];
 170     Int_t    ipart  = fXid[fNparticle];
 171     Float_t  pt     = fXpt[fNparticle];
 172     Float_t  eta    = fXeta[fNparticle];
 173     Float_t  phi    = fXphi[fNparticle];
 174     Float_t  theta  = 2.*TMath::ATan(TMath::Exp(-eta));
 175     Double_t amass  = TDatabasePDG::Instance()->GetParticle(ipart)->Mass();
 176
 177     p[0] = pt*TMath::Sin(phi);
 178     p[1] = pt*TMath::Cos(phi);
 179     p[2] = pt/TMath::Cos(theta);
 180     p[3] = TMath::Sqrt(pt*pt+p[2]*p[2]+amass*amass);
 181
 182
 183     TParticle* particle =
 184         new TParticle(ipart, 0, -1, -1, -1, -1, p[0], p[1], p[2], p[3],
 185                       0., 0., 0., 0.);
 186     fNparticle++;
 187     return particle;
 188 }
 189
 190
 191
 192 AliGenReaderEcalJets& AliGenReaderEcalJets::operator=(const  AliGenReaderEcalJets& rhs)
 193 {
 194 // Assignment operator
 195     rhs.Copy(*this);
 196     return (*this);
 197 }
 198
 199 void AliGenReaderEcalJets::Copy(TObject&) const
 200 {
 201     //
 202     // Copy
 203     //
 204     Fatal("Copy","Not implemented!\n");
 205 }
 206
 207
 208
 209
 210
 211
 212