3 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights
6 * See cxx source for full Copyright notice
8 //___________________________________________________________________
10 // The classes defined here, are utility classes for reading in data
11 // for the FMD. They are put in a seperate library to not polute the
12 // normal libraries. The classes are intended to be used as base
13 // classes for customized class that do some sort of analysis on the
14 // various types of data produced by the FMD.
15 /** @file AliFMDInput.h
16 @author Christian Holm Christensen <cholm@nbi.dk>
17 @date Mon Mar 27 12:42:40 2006
18 @brief FMD utility classes for reading FMD data
20 //___________________________________________________________________
21 /** @defgroup FMD_util Utility classes.
23 The classes defined here, are utility classes for reading in data
24 for the FMD. They are put in a seperate library to not polute the
25 normal libraries. The classes are intended to be used as base
26 classes for customized class that do some sort of analysis on the
27 various types of data produced by the FMD.
52 //___________________________________________________________________
53 /** @class AliFMDInput
54 @brief Base class for reading in various FMD data.
55 The class loops over all found events. For each event the
56 specified data is read in. The class then loops over all
57 elements of the read data, and process these with user defined
60 struct DigitInput : public AliFMDInput
67 fHist = new TH1F("adc", "ADC spectra", 1024, -.5, 1023.5);
70 Bool_t ProcessDigit(AliFMDDigit* d)
72 fHist->Fill(d->Counts());
75 // After processing all events, display spectrum
89 This class allows for writing small scripts, that can be compiled
90 with AcLIC, to do all sorts of tests, quick prototyping, and so
91 on. It has proven to be quiet useful. One can load more than
92 one type of data in one derived class, to for example to make
93 comparisons between hits and reconstructed points. See also the
94 various scripts in @c FMD/scripts.
97 class AliFMDInput : public TObject
100 /** The kinds of data that can be read in. */
103 kKinematics, // Kinematics (from sim)
105 kSDigits, // Summable digits
106 kHeader, // Header information
107 kRecPoints, // Reconstructed points
109 kRaw, // Read raw data
110 kGeometry, // Not really a tree
111 kTracks // Hits and tracs - for BG study
116 @param gAliceFile galice file */
117 AliFMDInput(const char* gAliceFile);
119 virtual ~AliFMDInput() {}
121 /** Add a data type to load
122 @param tree Data to load */
123 virtual void AddLoad(ETrees tree) { SETBIT(fTreeMask, tree); }
124 /** Remove a data type to load
125 @param tree Data to @e not load */
126 virtual void RemoveLoad(ETrees tree) { CLRBIT(fTreeMask, tree); }
127 /** @return # of available events */
128 virtual Int_t NEvents() const;
130 /** Initialize the class. If a user class overloads this member
131 function, then this @e must be explicitly called
132 @return @c false on error */
133 virtual Bool_t Init();
134 /** Callled at the beginning of each event. If a user class
135 overloads this member function, then this @e must be explicitly
137 @param event Event number
138 @return @c false on error */
139 virtual Bool_t Begin(Int_t event);
140 /** Process one event. This loops over all the loaded data. Users
141 can overload this member function, but then it's @e strongly
142 recommended to explicitly call this classes version.
143 @return @c false on error */
144 virtual Bool_t Event();
145 /** Called at the end of each event.
146 @return @c false on error */
147 virtual Bool_t End();
148 /** Called at the end of the run.
149 @return @c false on error */
150 virtual Bool_t Finish() { return kTRUE; }
152 @return @c false on error */
153 virtual Bool_t Run();
155 /** Loop over all hits, and call ProcessHit with that hit, and
156 optionally the corresponding kinematics track.
157 @return @c false on error */
158 virtual Bool_t ProcessHits();
159 /** Loop over all tracks, and call ProcessTrack with each hit for
161 @return @c false on error */
162 virtual Bool_t ProcessTracks();
163 /** Loop over all digits, and call ProcessDigit for each digit.
164 @return @c false on error */
165 virtual Bool_t ProcessDigits();
166 /** Loop over all summable digits, and call ProcessSDigit for each
168 @return @c false on error */
169 virtual Bool_t ProcessSDigits();
170 /** Loop over all digits read from raw data files, and call
171 ProcessRawDigit for each digit.
172 @return @c false on error */
173 virtual Bool_t ProcessRawDigits();
174 /** Loop over all reconstructed points, and call ProcessRecPoint for
175 each reconstructed point.
176 @return @c false on error */
177 virtual Bool_t ProcessRecPoints();
178 /** Loop over all ESD data, and call ProcessESD for each entry.
179 @return @c false on error */
180 virtual Bool_t ProcessESDs();
182 /** Process one hit, and optionally it's corresponding kinematics
183 track. Users should over this to process each hit.
185 @param p Associated track
186 @return @c false on error */
187 virtual Bool_t ProcessHit(AliFMDHit* h, TParticle* p);
188 /** Process one hit per track. Users should over this to process
190 @param i Track number
192 @param h Associated Hit
193 @return @c false on error */
194 virtual Bool_t ProcessTrack(Int_t i, TParticle* p, AliFMDHit* h);
195 /** Process one digit. Users should over this to process each
198 @return @c false on error */
199 virtual Bool_t ProcessDigit(AliFMDDigit* digit);
200 /** Process one summable digit. Users should over this to process
202 @param sdigit Summable digit
203 @return @c false on error */
204 virtual Bool_t ProcessSDigit(AliFMDSDigit* sdigit);
205 /** Process one digit from raw data files. Users should over this
206 to process each raw digit.
207 @param digit Raw digit
208 @return @c false on error */
209 virtual Bool_t ProcessRawDigit(AliFMDDigit* digit);
210 /** Process one reconstructed point. Users should over this to
211 process each reconstructed point.
212 @param point Reconstructed point
213 @return @c false on error */
214 virtual Bool_t ProcessRecPoint(AliFMDRecPoint* point);
215 /** Process ESD data for the FMD. Users should overload this to
217 @param d Detector number (1-3)
218 @param r Ring identifier ('I' or 'O')
219 @param s Sector number (0-19, or 0-39)
220 @param t Strip number (0-511, or 0-255)
221 @param eta Psuedo-rapidity
222 @param mult Psuedo-multiplicity
223 @return @c false on error */
224 virtual Bool_t ProcessESD(UShort_t, Char_t, UShort_t, UShort_t,
229 @param o Object to copy from */
230 AliFMDInput(const AliFMDInput& o)
258 /** Assignement operator
259 @return REference to this */
260 AliFMDInput& operator=(const AliFMDInput&) { return *this; }
262 TString fGAliceFile; // File name of gAlice file
263 AliRunLoader* fLoader; // Loader of FMD data
264 AliRun* fRun; // Run information
265 AliStack* fStack; // Stack of particles
266 AliLoader* fFMDLoader; // Loader of FMD data
267 AliRawReader* fReader; // Raw data reader
268 AliFMD* fFMD; // FMD object
269 AliESD* fMainESD; // ESD Object
270 AliESDFMD* fESD; // FMD ESD data
271 TTree* fTreeE; // Header tree
272 TTree* fTreeH; // Hits tree
273 TTree* fTreeD; // Digit tree
274 TTree* fTreeS; // SDigit tree
275 TTree* fTreeR; // RecPoint tree
276 TTree* fTreeA; // Raw data tree
277 TChain* fChainE; // Chain of ESD's
278 TClonesArray* fArrayE; // Event info array
279 TClonesArray* fArrayH; // Hit info array
280 TClonesArray* fArrayD; // Digit info array
281 TClonesArray* fArrayS; // SDigit info array
282 TClonesArray* fArrayR; // Rec points info array
283 TClonesArray* fArrayA; // Raw data (digits) info array
284 TGeoManager* fGeoManager; // Geometry manager
285 Int_t fTreeMask; // Which tree's to load
286 Bool_t fIsInit; // Have we been initialized
287 ClassDef(AliFMDInput,0) //Hits for detector FMD
290 inline Bool_t AliFMDInput::ProcessHit(AliFMDHit*,TParticle*) { return kTRUE; }
291 inline Bool_t AliFMDInput::ProcessTrack(Int_t,TParticle*,
292 AliFMDHit*) { return kTRUE; }
293 inline Bool_t AliFMDInput::ProcessDigit(AliFMDDigit*) { return kTRUE; }
294 inline Bool_t AliFMDInput::ProcessSDigit(AliFMDSDigit*) { return kTRUE; }
295 inline Bool_t AliFMDInput::ProcessRawDigit(AliFMDDigit*) { return kTRUE; }
296 inline Bool_t AliFMDInput::ProcessRecPoint(AliFMDRecPoint*) { return kTRUE; }
297 inline Bool_t AliFMDInput::ProcessESD(UShort_t,Char_t,UShort_t,UShort_t,
298 Float_t,Float_t) { return kTRUE; }
302 //____________________________________________________________________