AliRoot Core  3dc7879 (3dc7879)
AliVertexGenFile.cxx
Go to the documentation of this file.
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 
19 // //
20 // Generator for vertices taken from a file //
21 // //
22 // The file name of the galice file is passed as argument to the //
23 // constructor. If a second argument is given, this determines the number //
24 // of events for which the same vertex is used. //
25 // //
27 
28 
29 #include <TArrayF.h>
30 #include <TFile.h>
31 #include <TTree.h>
32 
33 #include "AliVertexGenFile.h"
34 #include "AliLog.h"
35 #include "AliGenEventHeader.h"
36 #include "AliHeader.h"
37 
38 
39 ClassImp(AliVertexGenFile)
40 
41 
42 //_____________________________________________________________________________
44  fFile(NULL),
45  fTree(NULL),
46  fHeader(NULL),
47  fEventsPerEntry(0),
48  fEvent(0),
49  fLastTime(0)
50 {
51 // default constructor: initialize data members
52 
53 }
54 
55 //_____________________________________________________________________________
57  Int_t eventsPerEntry) :
58  fFile(NULL),
59  fTree(NULL),
60  fHeader(NULL),
61  fEventsPerEntry(eventsPerEntry),
62  fEvent(0),
63  fLastTime(0)
64 {
65 // main constructor:
66 // fileName is the name of the galice file containing the vertices
67 // eventsPerEntry is the number of events for which the same vertex is used
68 
69  TDirectory* dir = gDirectory;
70 
71  fFile = TFile::Open(fileName);
72  if (!fFile || !fFile->IsOpen()) {
73  AliError(Form("could not open file %s", fileName));
74  delete fFile;
75  fFile = NULL;
76  return;
77  }
78  fTree = (TTree*) fFile->Get("TE");
79  if (!fTree) {
80  AliError(Form("no header tree found in file %s", fileName));
81  dir->cd();
82  return;
83  }
84  fHeader = new AliHeader;
85  fTree->SetBranchAddress("Header", &fHeader);
86 
87  dir->cd();
88 }
89 
90 //_____________________________________________________________________________
92 {
93 // clean up
94 
95  if (fFile) fFile->Close();
96  delete fFile;
97  delete fHeader;
98 }
99 
100 //_____________________________________________________________________________
102 {
103  // get the timestamp of the last header used for the vertex
104  if (fHeader) return fHeader->GetTimeStamp();
105  else AliFatal("No header was loaded yet");
106  return 0;
107 }
108 
109 //_____________________________________________________________________________
111 {
112 // get the vertex from the event header tree
113 
114  Int_t entry = fEvent++ / fEventsPerEntry;
115  fLastTime = 0;
116  if (!fTree) {
117  AliError("no header tree");
118  return TVector3(0,0,0);
119  }
120 
121  if (fTree->GetEntry(entry) <= 0) {
122  AliError(Form("error loading entry %d", entry));
123  return TVector3(0,0,0);
124  }
125 
126  if (!fHeader->GenEventHeader()) {
127  AliError("no generator event header");
128  return TVector3(0,0,0);
129  }
130 
131  TArrayF vertex(3);
134  return TVector3(vertex[0], vertex[1], vertex[2]);
135 }
136 
137 
virtual ~AliVertexGenFile()
TFile * Open(const char *filename, Long64_t &nevents)
TTree * fTree
Definition: MakeTreeStat.C:55
time_t GetHeaderTimeStamp() const
virtual void PrimaryVertex(TArrayF &o) const
Float_t fLastTime
current event number
AliHeader * fHeader
tree with headers
TString fileName(const char *dir, int runNumber, const char *da, int i, const char *type)
virtual Float_t InteractionTime() const
virtual TVector3 GetVertex()
Int_t fEventsPerEntry
event header
virtual time_t GetTimeStamp() const
Definition: AliHeader.h:63
#define AliFatal(message)
Definition: AliLog.h:640
virtual AliGenEventHeader * GenEventHeader() const
Definition: AliHeader.cxx:244
#define AliError(message)
Definition: AliLog.h:591
TTree * fTree
galice file with vertices